\documentclass[english, a4paper, 10pt, headings=small, DIV11]{scrartcl}
\usepackage[utf8]{inputenc} % Sweave doesn't find the encoding definition in vignettes.defs
\SweaveInput{vignettes.defs}
\SweaveOpts{prefix.string=laser-fig-}
\hypersetup{pdftitle={laser},
 pdfauthor={C. Beleites},
 pdfsubject={Vignette on hyperSpec's laser data set},
 pdfkeywords={hyperSpec, laser, time series}}

% \VignetteIndexEntry{laser: Vignette on hyperSpec's laser data set.}
% \VignetteKeywords{hyperSpec, laser, time series}
% \VignettePackage{hyperSpec}
\begin{document}
\SweaveOpts{concordance=TRUE}
\title{Unstable Laser Emission\\
Vignette for the Data Set \Robject{laser} of the R package \Rpackage{hyperSpec} }
\maketitle

\warnbox{Reproducing the Examples in this Vignette}{
All spectra used in this manual are installed automatically with \phy. 

The source data files can be found with the command:\\
\texttt{\small\textsl{> system.file ("doc/src/rawdata/laser.txt.gz", package = "hyperSpec")}}.

In order to reproduce the examples by typing in the commands, have a look at the definitions in \texttt{vignettes.defs}.}

\tableofcontents
\warnbox[blue!50!black]{Suggested Packages}{
\begin{labeling}{latticeExtra}
<<check-required, echo = FALSE, results = tex>>=
check.req.pkg ("rgl", donothing = "")
@
\end{labeling}
}

\section{Introduction}
This data set consists of a time series of \Sexpr{nrow (laser)} spectra of an
unstable laser emission at 405 nm recorded during ca. 1.5 h. 

The spectra were recorded  during the installation of the 405~nm laser at a Raman spectrometer.
There is no Raman scattering involved in this data, but the Raman software recorded the
abscissa of the spectra as Raman shift in wavenumbers.

This document shows
\begin{itemize}
\item How to convert the wavelength axis (spectral abscissa)
\item How to display time series data as intensity over time diagram
\item How to display time series data as 3d and false colour image of intensity as function of time and wavelength.
\end{itemize}

\section{Loading the Data and Preprocessing}

<<rawspc, fig = TRUE, width = 5, height = 3, include = FALSE>>=
laser <- scan.txt.Renishaw ("rawdata/laser.txt.gz", data = "ts")
plot (laser, "spcprctl5") 
@

As the laser emission was recorded with a Raman spectrometer, the wavelength axis initially is the Raman shift in wavenumbers (cm\textsuperscript{-1}).

As most of the spectra do not show any signal (fig. \ref{fig:rawspc}), so the spectral range
can be cut to -75~--~0~cm\textsuperscript{-1}. 
Note that negative numbers in the spectral range specification with the tilde do not exclude the spectral range but rather mean negative values of the wavelength axis.
The results are shown in figure \ref{fig:cutspc}.
<<cut, fig = TRUE, width = 5, height = 3, include = FALSE>>=
laser <- laser [,,-75~0]
plot (laser, "spcprctl5") 
@
\begin{figure}[tbh]
  \centering
  \subfloat[\label{fig:rawspc} The raw spectra.]{\includegraphics[width=.33\textwidth]{laser-fig--rawspc}}
  \subfloat[\label{fig:cutspc} The cut spectra.]{\includegraphics[width=.33\textwidth]{laser-fig--cut}}  
  \subfloat[\label{fig:nmspc} The spectra with wavelength axis.]{\includegraphics[width=.33\textwidth]{laser-fig--wlspc}}
  \caption{The laser emission spectra.}
  \label{fig:raw}
\end{figure}

The wavelength axis was recorded as Raman shift from 405 nm. However,
the spectra were taken before calibrating the wavelength axis. 
The band at -50 cm\textsuperscript{-1} is known to be at 405 nm. 
<<wlspc1>>=
wl (laser) <- wl (laser) + 50
@ 
Furthermore, as the spectra are not Raman shift but emission, the wavelength axis should be converted to proper wavelengths in nm.

The Raman shift is calculated from the wavelength as follows:
\[ \Delta\tilde\nu = \frac{1}{\lambda_0} -  \frac{1}{\lambda}\]
with $\Delta\tilde\nu$ being the Raman shift, and $\lambda_0$
the excitation wavelength for a Raman process, here 405~nm.

The wavelengths corresponding to the wavenumbers are thus:
\[ \lambda = \frac{1}{ \frac{1}{\lambda_0} - \Delta\tilde\nu}\]

Taking into account that 1~cm = 10$^7$ nm, we arrive at the new
wavelength axis:
<<wlcalc>>=
wl (laser) <- list (
   wl = 1e7 / (1/405e-7 - wl (laser)),
   label = expression (lambda / nm)
   )
@ 
<<wlspc, fig = TRUE, width = 5, height = 3, include = FALSE >>=
plot (laser, "spcprctl5") 
@ 
Note that the new wavelength axis label is immediately assigned as well.

Now, save \Robject{laser} as the \Robject{laser} data set shipped with \Rpackage{hyperSpec}.
<<save>>=
save (laser, file = "laser.rda") 
laser
@ 

\section{Inspecting the time dependency of the laser emission}
The maxima of the different emission lines encountered during this
measurement are at 405.0, 405.1, 405.3, and 405.4\,nm (fig. \ref{fig:markedspc}).

Alternatively they can be extracted from the graph using \Rfunction{locator} which reads out the coordinates of the points the user clicks with the mouse (use middle or right click to end the input):
<<locator, eval=FALSE>>=
wls <- locator()$x
@ 
<<echo = FALSE, results = hide>>=
wls <-  c(405.0063, 405.1121, 405.2885, 405.3591)
@
<<markspc, fig = TRUE, width = 5, height = 3, include = FALSE>>=
plot (laser, "spcmeansd")
cols <- c("black", "blue", "red", "darkgreen")
abline (v = wls, col = cols )  
@ 

\Rfunction{plotc} can also be used to plot time-series. In that case,
the abscissa needs to be specified in parameter \Rcommand{use.c}. The collection
time is stored in column \Rcode{\$t} in seconds from start of the
measurement, and can be handed over as the column name. The resulting
time series are shown in figure \ref{fig:ts}
variable, 
<<ts, fig = TRUE, width = 10, height = 3, include = FALSE>>=
plotc (laser [,, wls], spc ~ t, groups = .wavelength, type = "b",  cex = 0.3, col = cols)
@ 
\begin{figure}[t]
  \centering
  \subfloat[\label{fig:markedspc}]{\includegraphics[width=.33\textwidth]{laser-fig--markspc}}
  \subfloat[\label{fig:ts}]{\includegraphics[width=.66\textwidth]{laser-fig--ts}}  
  \caption[]{The laser emission time series. \ref{fig:markedspc} shows
    the spectral position of the bands. The time series are plotted in
  corresponding colors in \ref{fig:ts}.}
  \label{fig:tsf}
\end{figure}

Another option is to condition the plot on $\lambda$:
<<tsextra, fig = TRUE, width = 12, include = FALSE>>=
plotc (laser [,, wls], spc ~ t | .wavelength, type = "b", cex = 0.3, col = "black")
@ 
\begin{figure}[t]
  \centering
  \includegraphics[width=\textwidth]{laser-fig--tsextra}
  \caption[]{The time series plots can also be conditioned on \Rcode{\$.wavelength}.}
  \label{fig:tsextra}
\end{figure}

\section{False-colour plot of the spectral intensity over wavelength and time}
\label{sec:plotmat}
\phy supplies functions to draw the spectral matrix using \Rpackage{lattice}'s \Rfunction{levelplot}.
<<plotmatr, fig = TRUE, include = FALSE>>=
plot (laser, "mat", contour = TRUE, col = "#00000060")
@

\phy's \Rfunction{levelplot} method can be used to display the spectra matrix over a data column (instead of the row number): figure~\ref{fig:tsmat}. Note that the \chy object is the \emph{second} argument to the function (according to the notation in \Rfunction{levelplot}).
<<plotmatt, fig = TRUE, include = FALSE>>=
levelplot (spc ~ .wavelength * t, laser, contour = TRUE, col ="#00000080")
@
\begin{figure}[tbh]
  \centering
  \subfloat[\label{fig:row}]{\includegraphics[width=.5\textwidth]{laser-fig--plotmatr}}
  \subfloat[\label{fig:t}]{\includegraphics[width=.5\textwidth]{laser-fig--plotmatt}}  
  \caption[]{The spectra matrix of the \Robject{laser} data set. The ordinate of the plot may be the number of the spectrum accessed by \Rcode{\$.row} \subref{fig:row} or any other extra data column, here \Rcode{\$t} \subref{fig:t}.}
  \label{fig:tsmat}
\end{figure}


\section{3d plot of the spectral intensity over wavelength and time}
\label{sec:3d-plot-spectral}
\chy objects  can easily be drawn with \Rpackage{rgl}\Sexpr{cite.pkg ("rgl")}:

<<libraryrgl>>=
require (rgl)
@
\Rpackage{rgl}'s \Rfunction{persp3d} plots a surface in 3d defined by points in x, y, and z. Handing over the appropriate data columns of the \chy object is easy (fig.~\ref{fig:laser3d}):

<<rgl-plot, eval = FALSE>>=
message ("plot chunk")
laser <- laser [,,404.8 ~ 405.6] / 10000
laser$t <- laser$t / 3600
cols <- rep (matlab.palette (nrow (laser)), nwl (laser))

surface3d(y = wl(laser), x = laser$t, z = laser$spc, col =  cols)
surface3d(y = wl(laser), x = laser$t, z = laser$spc + .1 * min (laser), 
          col =  "black", alpha = .2,front = "lines", line_antialias = TRUE)

aspect3d (c(1, 1, 0.25))

axes3d (c ('x+-', 'y--', 'z--'))
axes3d ('y--', nticks = 25, labels= FALSE)
mtext3d ("t / h",    'x+-', line = 2.5)
mtext3d ("lambda / nm",   'y--', line = 2.5)
mtext3d ("I / a.u.", 'z--', line = 2.5)
@ 
<<rgl-do, echo = FALSE, results = hide>>=
  if (require (rgl)){
    open3d (windowRect=c(20,20,600, 350))  # this is needed only for automatically 
                                        # producing the snapshot
<<rgl-plot>>
     par3d (userMatrix = matrix (c (-0.52,  0.4, -0.75, 0, 
                                    -0.85, -0.28, 0.44, 0, 
                                    -0.04,  0.87, 0.49, 0, 
                                    -0.75,  0.75,     0, 1), ncol = 4L),
            scale = c (2.75, 5, 0.175),
            windowRect = c(20L, 50L, 520L, 330L),
            zoom = 0.75)
    rgl.snapshot ("fig-3D.png", fmt="png", top=TRUE )
    rgl.quit ()
  } else {
    png ("fig-3D.png")
    ploterrormsg ("", "rgl")
    dev.off ()
  }
@
\begin{figure}[t]
  \centering
  \includegraphics[width = .5\textwidth]{fig-3D}  
  \caption{The 3d plot of the laser data}
  \label{fig:laser3d}
\end{figure}

<<makebib, echo=FALSE, results = hide>>=
make.bib (file = "laser-pkg.bib")
@
\bibliographystyle{unsrtnat}
\bibliography{laser-pkg}

\section*{Session Info}
\begin{footnotesize}
<<echo = FALSE>>= 
make.bib (file = "laser-pkg.bib") 
<<cleanup>> 
@
\end{footnotesize}
\end{document}
