source: Sophya/trunk/SophyaLib/Manual/pmixer.tex@ 2750

\documentclass[11pt]{article}
%\documentstyle[11pt,epsfig,url]{article}
\usepackage[latin1]{inputenc}
\usepackage[T1]{fontenc}
\usepackage[english]{babel}
\usepackage{url}
\usepackage{graphicx}

% package a mettre pour faire du pdf
\usepackage{palatino}

%  Definition pour Docs Sophya
\usepackage{defsophya}



\setlength{\parskip}{\smallskipamount} %\documentclass[a4paper,10pt]{article}
\newenvironment{maliste}%
{ \begin{list}%
        {$\ast$}%
        {\setlength{\labelwidth}{20pt}%
         \setlength{\leftmargin}{25pt}%
         \setlength{\itemsep}{\parsep}}}%
{ \end{list} }
\newcommand {\ddate}    {April~12,~2000}

\begin{document}
%  The title page - top of the page with the title of the paper
\begin{titlepage}
\titrehp{The SkyMixer \\ (SkyT and PMixer modules)}
%  Authors list
\auteurs{
R. Ansari            &  ansari@lal.in2p3.fr       \\
S. Henrot-Versille   &  versille@in2p3.fr     
}
%  The title page - bottom of the page with the paper number
\titrebp{2}
\end{titlepage}

%--------------------------------------- Abstract ---------------------
\begin{abstract}
This note gives a description of the SkyMixer which makes
use of the PMixer and the SkyT
modules of the Sophya library.
In a first part, we provide a user-guide on how to use the
skymixer and the extractRS program, and how to
fill the corresponding datacards.
In a second part, we describe more technically
the different classes
of the SkyT module. 
\end{abstract}
\newpage
\tableofcontents
\newpage


\section{Skymixer and its datacard}
The skymixer program adds several sky component maps
simulated or measured at several frequencies, such as,
for example, foregrounds (synchrotron radiation, dust, radio-sources...)
and the microwave background.
It takes into account the corresponding radiation spectra
$F_\nu$ (the one of a black-body for microwave background,
a power law spectrum for synchrotron...),
and integrates them over a given detector filter frequency-re\-sponse
$W_\nu$ according to:
\begin{eqnarray}
S &=& \sum_{sky\  components} \int_0^\infty F_\nu W_\nu d\nu \ .
\end{eqnarray}
The output maps\footnote{$F_\nu$ is in $\hbox{W}/\hbox{m}^2/\hbox{sr}/\hbox{Hz}$
and $d\nu$ is in $\hbox{GHz}$, while $W_\nu$ is dimensionless.} are in $10^{-9}\hbox{W}/\hbox{m}^2/\hbox{sr}$.

The interaction between the user and the skymixer is done
using ASCII cards. The keywords one may use are described below.

The way to call the skymixer program is the following one:
\begin{verbatim}
 Usage: skymixer paramFile outputfits [outppfname]
\end{verbatim}
where paramFile
is the name of the datacard (described below) and       
outputfits is the name of the Fits file in which the
final map will be saved\footnote{PPFName stands for Planck Persistent File (which can be read using piapp or spiapp).}.
When asking ``skymixer -h'', the previous usage line will appear on the
screen.

In the datacards, the lines are activated with the $@$ symbol
at the start of line. All other lines are treated as comments.
However for clarity, we use  
$\#$ to mark commented lines.

\subsection{Input maps}
        One has to define the number of components one wants to add
        and the characteristics 
        of the input maps (namely the NSide of HealPix maps):
        \begin{verbatim}
        #  SKYMIX  NbComponents  NSide_HealPix
        @SKYMIX  1 2048
        \end{verbatim}
        Note that, within this scheme, only maps with the same NSide can
        be added.

        In the case the user wants to perform a check on the PDMTYPE
        keyword for the different inputs files, he can use the 
        above CHECKPDMTYPE command (note that by default, no check
        is performed).
        \begin{verbatim}
        @CHECKPDMTYPE
        \end{verbatim}


        In the case the files containing the 
        maps  are not in the
        working directory, it is possible to precise the path to
        their location:
        \begin{verbatim}
         #  MAPPATH  PathForFITS   
         #  path for input Sky Component FITS files , 
         #   and EmissionSpectra files 
         MAPPATH /exp/planck/data/
        \end{verbatim}

        To have access to different printout, one can play with the
        DEBUGLEVEL and PRINTLEVEL keywords:
        \begin{verbatim}
         #  Define the Debug level
         @DEBUGLEVEL  0
         #  Define the Print level
         @PRINTLEVEL  5
        \end{verbatim}
        when set to 0, almost no printout will be given.

\subsection{Detection filter}
        Then one has to define the detection frequency 
        filter response (cf. section
        \ref{specresp}). 
        For the example, for a gaussian filter (cf. section \ref{gaussFilt}), 
        one should use:
        \begin{verbatim}
         #  Defining the Detection filter pass - band
         #  GAUSSFILTER -> Gaussian filter
         #  GAUSSFILTER  Nu0 Sigma_Nu A NuMin NuMax
         @GAUSSFILTER 143 12 0.84 50. 300.
        \end{verbatim}
        Note that the frequencies have to be given 
        in $\hbox{GHz}$.
        In the case of a detection filter defined in a Fits file
        with 2$\times$n arrays (one for the frequencies and the other for
        the transmission values), the syntax would be:
        \begin{verbatim}
         #  Defining the Detection filter pass - band
         #  FILTERFITSFILE -> Filter (nu, T(nu)) from FITS 
         #  FILTERFITSFILE FileName NuMin NuMax
         @FILTERFITSFILE filterfile.fits 50. 230.
        \end{verbatim}

\subsection{Sky components}
        Once the detector filter has been defined, one has to 
        define the different sky components one wants to add
        and
        their corresponding radiation spectra. The maps are
        given the MAPFITSFILEi keyword where i corresponds
        to a given component (synchrotron radiation,
        dust...): i has to match between the map
        and the radiation spectrum. One also has to specify
        a normalization factor corresponding to each component.
        \begin{verbatim}
        #  ----        Sky components  ---------
        #  MAPFITSFILEi FITSfilename Normalisation
        @MAPFITSFILE4  dust.fits  1.0  
        \end{verbatim}   
        Here the number 4 (i=4) means that we are refering to
        the fourth sky components. And the normalization factor
        of the map is set to 1.
        
\subsection{Radiation Spectrum}
The radiation spectra corresponding to the sky component
maps can either be defined in a fits files, or by the instanciation
of the RadSpectra class (cf. \ref{radspec}), or they can depend
on the position in the map. The different cases are treated differently
as explained below. There is also the possibility of adding
a dipole component without the use of a map, as described below.

\subsubsection{A Fits file}
The user can provide the skymixer program with his own 
Fits file containing 2 $\times$ n arrays (one for the frequencies and
one for the corresponding flux values), and a frequency range.
It will be done for example using:
 \begin{verbatim}
 #      ---   Emission Spectra definition -----
 #   SPECTRAFITSFILEi  FITSfilename  NuMin NuMax 
 @SPECTRAFITSFILE4  pl.fits 100 600  
\end{verbatim}
The flux function will returns 0 outside of the 
$[ \hbox{Numin,Numax} ]$ range.

\subsubsection{The RadSpectra class}
For well-known spectra, the user does not have to provide
a Fits file, nor the function, but rather use already coded
functions, filling the desired parameter values.
For example, for a black body emission spectrum:
\begin{verbatim}
 #   Other possible definition of emission spectra
 #   BLACKBODYi  temperature
 #  The blackbody object returns flux in units of 
 #    watt/m2/sr 
 #   We thus scale the map (x 1.e17) 
 @MAPFITSFILE2  cossin.fits  1.e17 
 @BLACKBODY2 2.726
\end{verbatim}
The user has to give the temperature (here 2.726K).

Another example is given below where one wants to use a
power law spectrum for the radiation:
\begin{verbatim}
 #   For power-law spectrum  f = a ((nu-nu0)/dnu)^b)
 #   POWERLAWSPECTRAi  a nu0 dnu b Fmin Fmax
 @MAPFITSFILE3  noise.fits    0.1 
 @POWERLAWSPECTRA3 1. 150. 50. -0.5 100. 500.
\end{verbatim}

\subsubsection{Map position dependence}
In the case the emission spectrum depends on the position
in the map (for example for the case of synchrotron radiation),
one should use the BETAFITSFILE key word as follows:
\begin{verbatim}
 #  For an emission spectrum which depends  
 #  on the position in the map -> synchrotron 
 #  for instance (nu/nu0)^(-beta(theta,phi))
 #  408e-3 is the normalisation frequency of the map
 @MAPFITSFILE1 Sync/Template408MHz.fits 0.2
 #  BETAFITSFILE BetaMap NormFreq NsideRePixelization
 @BETAFITSFILE1 Sync/SpectralIndex.fits  408e-3 16 
\end{verbatim}
In this case, a second argument is needed for BETAFITSFILEi
which corresponds to the  normalization frequency value.
It is possible to perform the computation using a $\beta$ map
which has a different number of pixels than the sky map. 
It is also possible to change the NSideHealpix of the
$\beta$ map for example to low down the resolution for faster
computation. This is done using the third argument
(NsideRePixelization, here taken to be 16). For negative values
of NsideRePixelization, the original $\beta$ map resolution
is used.

One can refer to Eq.~\ref{pls} for the formula of the
radiation spectrum built from the $\beta$ map. In this equation, 
$\nu_0$ is set to 0 and $\delta \nu$ is given the value of
the normalisation
frequency of the map. 

\subsubsection{The dipole}      
The dipole can be treated in a different manner than the other
sky components,
using the DIPOLE card, as follows:
\begin{verbatim}        
# For a dipole sky map and emission spectrum
# use: DIPOLE theta phi amplitude delta_temperature
# theta and phi are in radian, delta_temp in K
@DIPOLE1 1.2  0.3  1. 3E-3
\end{verbatim}
In that case no MAPFITSFILEi is needed: the DIPOLE
card will create the sky map according to the
$(\theta,\phi)$ user requirement and integrate
the dipole radiation spectrum (DerivBlackBody - cf.~Eq.\ref{derivBB})
over the filter response.
Note that if the user prefers to give his own map, he
can do it as well.

\section{Treating point sources}
The extractRS code builds a full sky map from two data sets:
\begin{itemize}
\item{} a list of point sources with emission spectra described at a
set of discrete frequencies
\item{} maps of a small patch of the sky with emission spectrum described
at a set of discrete frequencies.
\end{itemize}

It treats the  case when one has made a simulation of a patch of the sky
producing a set of maps corresponding to the sky emission and the 
point sources emission at 
discrete frequencies,
and wants to build a random map on the basis of these
two informations.
The map produced by 
extractRS can then be given to the skymixer program to be
added to other sky components. 

As for skymixer, the user interacts with extractRS through
the use of a datacard (described below).  The way to call
it is the following one:
\begin{verbatim}
   Usage: extractRS parameterFile outputfits  [outppf]
\end{verbatim}
where parameterFile is the datacard, and outputfits
corresponds to the Fits file in which the final map will be
saved.

An example of the datacard is shown below:
\begin{verbatim}
#  EXT_RS  NbComponents  NSide_HealPix
@EXT_RS  1 256

#  RD_EXT_PATH  PathForFITS   
#  path for input Sky Component FITS files , 
#  from which the RadSpec files will be read
@RS_EXT_PATH /exp/planck/CartesRT/

#  Sky ``map'' for different frequencies 
#  From which the radspectra will be read
#  RADSPECMAPi FITSfilename 
@RADSPECMAP1    eran102_441_00015000mhz.fits

# Source point ``map'' for different freq
# From which the radspectra will be read
#  SOURCEPTMAPi FITSfilename freq(MHz)
@SOURCEPTMAP1   eran102_441_00015000mhz.s_rare

# Filter definition now !!
# define the detector (only one)
@GAUSSFILTER 143 12 0.84 50. 300.
\end{verbatim}
We first indicate in the datacard (after the EXT$\_$RS keyword):
\begin{maliste}
\item the number of frequency components we are going to
process: in the example the code will only process
one frequency component. 
\item and the NSideHealpix corresponding to the healpix output maps
\end{maliste}

Then comes the path (if any) to the maps, followed by the
name of the fits files containing the ``maps''
(these maps are not Healpix spheres but matrices). The keyword for
the frequency
components is RADSPECMAPi (with i being the sky component number)
and for the sources: SOURCEPTMAPi. 

The names of the fits
file should contain the frequency of the map in units of  $\hbox{MHz}$.
For example: 
\begin{verbatim}
@RADSPECMAP1    eran102_441_00015000mhz.fits
\end{verbatim}
corresponds to a frequency of 15000 $\hbox{Mhz}$.
The sky maps definition is followed by the point source
fits file, with the keyword SOURCEPTMAPi. 

Then comes the definition of the filter (here it is a
gaussian filter at 143 $\hbox{GHz}$).

\section{Simulating components and adding them}
For quick simulations, one can use two codes to generate maps,
detector filter responses and radiation spectra:
\begin{itemize}
\item{} {{\bf tgsky:}}\\
        this code generates a sky map in a Fits file.   
\item{} {{\bf tgrsr:}}\\
        this code writes in a Fits file the
        definition of a RadSpectraVec or of a SpecRespVec       
        (i.e. one frequency vector and the corresponding
        transmission or flux vector.)
\end{itemize}
\noindent The way to use these two functions is the following one:
\begin{itemize}
\item{} {{\bf tgsky:}}
\begin{verbatim}
  Usage: tgsky NSide TypSky Params FitsFile [PPF] 
\end{verbatim}
where PPF stands for Planck Persistent File (which can be read using piapp or spiapp),
NSide is the Healpix map parameter and TypSky and Params are described here:
\begin{maliste}
\item  if TypSky = 0, the map will be filled with a gaussian noise,
        the parameters (Params) are therefore the mean value and the sigma
        of the gaussian.
\item  if TypSky = 1, the map will be filled with a $1/\hbox{f}$ noise using
        the OOFNoise generator. The parameters to give are the
        sigma of the noise and a constant offset.
\item  if TypSky = 2, the map is filled with an analytical function of the
        form:   
        \begin{eqnarray}
        f &=& K \cos(a\theta)\sin(b\phi)+m
        \end{eqnarray}
        and the corresponding parameters are (in this order): K,a,b,m.
\item  if TypSky = 3, the map is filled with a dipole,
        the corresponding parameters are (in this order): theta,phi
        corresponding to the dipole orientation
\end{maliste}
\item{}{{\bf tgrsr:}}
\begin{verbatim}
 Usage: tgrsr choice Params NPoints FitsFileName 
                      [ImageR4FITSName] [ppfname]
\end{verbatim}
where choice, and Params are described below, and NPoints corresponds
to the number of elements in the generated vectors.
\begin{maliste}
\item if choice = 0, the code will generate a gaussian spectrum in a Fits file.
        Params are therefore the parameters of the gaussian plus the
        minimum and the maximum frequency.
        In the order:  A,Nu0,DNu,MinFreq,MaxFreq.
\item if choice = 1, the code will generate a power law spectrum in a 
        Fits file. Params are therefore, in the order:\\
        A, Nu0, dnu, b, MinFreq, MaxFreq
\end{maliste}
\end{itemize}



\section{The SkyT module}
The SkyT module is composed of two types of classes:
\begin{itemize}
\item{} one which corresponds to an emission spectrum of
radiation, which is called RadSpectra
\item{} one which corresponds to the spectral response
of a given detector (i.e. corresponding to a detector
filter in a given frequency domain), which is called
SpectralResponse.
\end{itemize}
In the coming sub-sections, we summarize the different functions
of these two classes and the corresponding subclasses.
For more informations about the classes, and the functions
the user may want to use (to have access or reset some
parameter values for example), we refer the reader to 
\cite{DOX-Sophya}.

All frequencies are expressed in $\hbox{GHz}$, which enables
us to represent wavelengths from the radio up to the optical 
range at a reasonnable accuracy using double values.

\subsection{ConvTools}
The ConvTools class is a very simple class where are
defined some needed constants for the computation of the
spectral emission function, such as:
\begin{verbatim}
// Boltzman constant  J/K
double ConvTools::kb =  1.380662e-23;     
// Mean value of the cmb temperature (K)
double ConvTools::tcmb=2.728;
// Planck constant in Js
double ConvTools::hpl=6.626176e-34;   
// Light speed  PDG 1992 m/s
double ConvTools::c=299792458.; 
\end{verbatim}

\subsection{SpectralResponse}
\label{specresp}
SpectralResponse 
is an abstract class representing 
a detector frequency filter
response, $W_\nu$ as a function of the frequency.
Several subclasses
have been written within the SkyT module which are listed
in the following subsections.
The constructor takes as arguments the
minimum and maximum values of the frequency range on
which the detector response is defined

The main public member function is the  
transmission function which has to be implemented
in the subclasses and is purely virtual at this stage:
\begin{verbatim}
virtual double transmission(double nu) const =0 ; 
\end{verbatim}
The () operator has also been redefined so that
one can access the transmission value of the detector's
filter directly writting for example myFilter(myfrequency).

Other member functions  are implemented at this level, but still
are virtual functions that can be re-implemented in the
sub-classes: they respectively returns the minimum, the
maximum, the mean and the peak values of the frequency
range (where the peak is defined as the maximum of the transmission
spectrum), but also the value of the transmission at the so-called
peak. One can also change the range in frequency of the filter
using the {{\it setMinMaxFreq}} function.
%\begin{verbatim}
% virtual double minFreq() const; 
% virtual double maxFreq() const; 
% virtual double meanFreq() const;
% virtual double peakFreq() const;
% virtual double peakTransmission() const;
% virtual void  setMinMaxFreq(double numin, double numax);
%\end{verbatim}
For all the sub-classes, $\nu$ is given in units of $GHz$ and
$W_\nu$ is dimensionless.

Since the user may want to integrate the detector response
within a given frequency range:
\begin{eqnarray}
\int_{\nu_{min}}^{\nu_{max}} F_\nu \delta \nu\ \ \ ,
\end{eqnarray}
the following functions are
available:
\begin{verbatim}
double IntegratedSpect(double numin,double numax) 
double logIntegratedSpect(double numin,double numax)
double IntegratedSpect() 
double logIntegratedSpect() 
\end{verbatim}
When no frequency range is defined the default values are used.
As, for some spectra, the calculation is more accurate when
one performs the numerical integration using a logarithmic
frequency scale.

The following subsections are devoted to the description of the
SpectralResponse subclasses which are already implemented.
\subsubsection{SquareFilter}
\begin{center}
\dclsbb{SpectralResponse}{SquareFilter}
\end{center}
The square filter is the simplest filter one can think of:
it is defined by a minimum and a maximum frequencies. 
\begin{verbatim}
SquareFilter (double numin, double numax)
\end{verbatim}
The transmission function returns 1 when the tested frequency
is within this range and 0 outside.
For the peakFrequency and peakTransmission functions, the square
filter returns respectively the frequency and the transmission
corresponding to the middle of the frequency range.
To use the SquareFilter class, one should include the squarefilt.h file.

\subsubsection{TriangleFilter}
\begin{center}\dclsbb{SpectralResponse}{TriangleFilter}\end{center}
For the triangle filter, as it was the case for the
square filter, its name itself is sufficient to define it:
the corresponding constructor is:
\begin{verbatim}
 TriangleFilter (double numin, double numax,
                double nupeak, double valueatpeak); 
\end{verbatim}
To use the TriangleFilter class, one should include the trianglefilt.h file.

\subsubsection{GaussianFilter}
\label{gaussFilt}
\begin{center}\dclsbb{SpectralResponse}{GaussianFilter}\end{center}
The constructor of the gaussian filter defines the parameters
of the corresponding gaussian:
\begin{verbatim}
 GaussianFilter (double nu0, double s, double A=1., 
           double numin=0., double numax= 9.e49); 
\end{verbatim}
where the transmission function returns:
\begin{eqnarray}
\label{gaussian}
W_\nu&=& A  e^{-{({\nu-\nu_0/ s})^2}};
\end{eqnarray}
%There are functions to reset or
%get back the value of the parameters (there are
%described in \cite{}. 
%using the 
%following functions:
%\begin{verbatim}
%  double giveNorm()
%  double giveNu0()
%  double giveDNu()
%\end{verbatim}
%or reset them with:
%\begin{verbatim}
%  virtual void setParams(double nu0, double s, double a=1.);
%\end{verbatim}
To use the GaussianFilter class, one should include the gaussfilt.h file.

\subsubsection{SpecRespVec}
\begin{center}\dclsbb{SpectralResponse}{SpecRespVec}\end{center}
One may also want to give the filter of a detector as a 
set of frequencies and the corresponding filter transmission values,
when there is no analytic expression to describe it.
In that case he/she should use the SpecRespVec class:
\begin{verbatim}
  SpecRespVec(Vector const & nu, Vector const & fdenu, 
              double numin, double numax)
\end{verbatim}
The first vector (called nu) is the frequency vector and the
second one (called fdenu) is the corresponding transmission:
note that the two vectors should have the same length, otherwise
an exception will be launched ({\it SzMismatchError}).

%One can access the different vectors with the use
%of appropriate functions (cf. \cite{}).
%Several functions are defined to access these vectors, either
%as a whole or element per element:
%\begin{verbatim}
%  Vector&  getNuVec()  { return _vecOfNu; } ;
%  Vector& getTNuVec()  { return _vecOfFDeNu; }
%  double getNuVec(int i) const { return _vecOfNu(i); }
%  double getTNuVec(int i) const { return _vecOfFDeNu(i); }
%  int NbElts() const { return _size;}
%\end{verbatim}
To use the SpecRespVec class, one should include the specrespvector.h file.

\subsection{RadSpectra}
\label{radspec}
The name RadSpectra stands for Radiation Spectrum.
The aims of RadSpectra is to define a base class
for the radiation emission spectra. It is therefore
abstract and contains one virtual member function 
which has to be implemented in the subclasses:
\begin{verbatim}
    virtual double flux(double nu) const =0;
\end{verbatim}
This function returns the value of the
flux (the spectral energy distribution), $F_\nu$,
as a function of the frequency $\nu$,
As in the SpectralResponse class,
the () operator has been redefined at this level, so that
the user can access the flux value, either by calling the
function or directly by using this operator.

The constructor takes as an argument the minimum
and the maximum frequency of the spectrum, if any.
In the case the user does not want to specify these
values, there are set respectively to 0. and $9.10^{49}$
by default.
For all the sub-classes, $\nu$ is given in units of $\hbox{GHz}$ and
the flux is returned in units of $\hbox{W}/\hbox{m}^2/\hbox{sr}/\hbox{Hz}$.

%We also defined functions to access the minimum, maximum,
%peak frequency values, or to change them 
%(see \cite{}). 
%And the user may change the
%values of the frequency range using the function:
%\begin{verbatim}
%virtual void    setMinMaxFreq(double numin, double numax);
%\end{verbatim}
One may also want to  integrate the emission
spectrum: this can be done in a given frequency
range, or in the full range of definition of the spectrum,
or  taking into account the
SpectralResponse (in frequency) of a detector (see section \ref{specresp}):
\begin{verbatim}
  virtual double integratedFlux
            (double numin, double numax)
  virtual double integratedFlux()
  virtual double filteredIntegratedFlux 
            (SpectralResponse const & filter, 
             double numin, double numax)
  virtual double filteredIntegratedFlux
            (SpectralResponse const & filter)
  virtual double logIntegratedFlux
            (double numin, double numax) 
  virtual double logIntegratedFlux() 
  virtual double filteredLogIntFlux
            (SpectralResponse const & filter, 
             double numin, double numax)
  virtual double filteredLogIntFlux
            (SpectralResponse const & filter);
\end{verbatim}

Several sub-classes of RadSpectra have been defined, which are listed
below with their corresponding characteristics.
\subsubsection{BlackBody}
\begin{center}\dclsbb{RadSpectra}{BlackBody}\end{center}
This class corresponds to the emission spectrum of a 
blackbody radiation: by default the temperature of the 
blackbody is 2.728K (as defined in convtools.cc), but when
the object is instanciated the temperature 
can be given any other value.
\begin{verbatim}
  BlackBody (double temp = ConvTools::tcmb);
\end{verbatim}
The flux function is the well known formula: 
\begin{eqnarray}
F_\nu &=& {2 h_{pl} (1.10^9*\nu)^3 / c^2 (e^{{h_{pl}(1.10^9*\nu) / kT}} -1)}
\end{eqnarray}

%In addition to the flux member function, we can get or set the
%temperature of the blackbody spectrum (see \cite{}).
%\begin{verbatim}
%  double getTemperature() 
%  void setTemperature(double newtemp)
%\end{verbatim}
To use the BlackBody class, one should include the blackbody.h file.


\subsubsection{DerivBlackBody}
\begin{center}\dclsbb{RadSpectra}{DerivBlackBody}\end{center}
This class corresponds to the emission spectrum of a 
the dipole (since its emission spectrum is the derivation
of a blackbody spectrum wrt the temperature): 
by default the temperature of the 
DerivBlackBody is set 2.728K (as defined in convtools.cc), but when
the object is instanciated the temperature 
can be given any other value.
\begin{verbatim}
  DerivBlackBody (double temp = ConvTools::tcmb);
\end{verbatim}
The flux function is simply given by: 
\begin{eqnarray}
\label{derivBB}
F_\nu &=& {2 h_{pl} (1.10^9*\nu)^3  {h_{pl}1.10^9*\nu \over k T^2}
        {e^{{h_{pl}(1.10^9*\nu) \over kT}}\over c^2 (e^{{h_{pl}(1.10^9*\nu) \over kT}} -1)^2}}
\end{eqnarray}

%In addition to the flux member function, we can get or set the
%temperature as for the spectrum (see \cite{}).:
%\begin{verbatim}
%  double getTemperature() 
%  void setTemperature(double newtemp)
%\end{verbatim}
To use the DerivBlackBody class, one should include the derivblackbody.h file.

\subsubsection{PowerLawSpectra}
\begin{center}\dclsbb{RadSpectra}{DerivBlackBody}\end{center}
This class corresponds to a power law spectrum, the
flux function being defined as follows:
\begin{eqnarray}
\label{pls}
I_\nu &=& A ({\nu-\nu_0\over \delta\nu})^\beta 
\end{eqnarray}
The constructor takes as arguments the parameters of this
function and therefore reads:
\begin{verbatim}
PowerLawSpectra(double a, double beta, double nu0, 
        double dnu, double numin=0., double numax=9.e49);
\end{verbatim}
%The intrinsic member functions of this class 
%helps getting back these parameters, once an object of this class
%has been
%instanciated:
%\begin{verbatim}
%  double giveNorm()
%  double giveNu0()
%  double giveDNu()
%  double giveExp() 
%\end{verbatim}
%One can also set them using:
%\begin{verbatim}
%  void setNorm(double myvalue)
%  void setNu0(double myvalue)
%  void setDNu(double myvalue)
%  void setExp(double myvalue) 
%\end{verbatim}
%
To use the PowerLawSpectra class, one should include the nupower.h file.


\subsubsection{RadSpectraVec}
\begin{center}\dclsbb{RadSpectra}{RadSpectraVec}\end{center}
As in the case of the SpectralResponse class, one may define the radiation
spectrum with two vectors: one for the frequencies and the second for the
values of the flux function. In that case, the class to use is RadSpectraVec
and the corresponding constructor is:
\begin{verbatim}
 RadSpectraVec (Vector const nu, Vector const fdenu, 
                double numin=0., double numax=9.e49);
\end{verbatim}
%The following functions helps getting back the input vectors, or their
%number of elements:
%\begin{verbatim}
%  Vector&  getNuVec()
%  Vector& getTNuVec()
%  double getNuVec(int i)
%  double getTNuVec(int i)
%  int NbElts() 
%\end{verbatim}
To use the RadSpectraVec class, one should include the radspecvector.h file.

\begin{thebibliography}{99}
\bibitem{DOX-Sophya}
\url{http://www.sophya.org/}
\end{thebibliography}
\end{document}