Changeset 1218 in Sophya


Ignore:
Timestamp:
Oct 3, 2000, 2:14:14 PM (25 years ago)
Author:
ansari
Message:

doc dans .cc

Location:
trunk
Files:
6 edited

Legend:

Unmodified
Added
Removed
  • trunk/SophyaExt/FitsIOServer/fitsfile.cc

    r1209 r1218  
    88
    99
     10void BnTblLine::setFormat(int dc, int fc, int ic, int cc, vector<string> names)
     11   {
     12    int nbcols = dc + fc + ic + cc;
     13    int maxName = names.size();
     14    if (nbcols != maxName)
     15      {
     16        cout << " WARNING: BnTblLine:: length of vector of column names not equal to total number of columns" << endl;
     17        maxName = nbcols < maxName ? nbcols : maxName;
     18      }
     19    ColName_ = vector<string>(nbcols);
     20     for (int k=0; k < maxName; k++) ColName_[k] = names[k];
     21    if (dc >0) ddata_ = vector<double>(dc);
     22    if (fc >0) fdata_ = vector<float>(fc);
     23    if (ic >0) idata_ = vector<int>(fc);
     24    if (cc >0) cdata_ = vector<string>(fc);
     25   }
     26
     27bool BnTblLine::sameFormat(const BnTblLine& btl) const
     28   {
     29     if (btl.ddata_.size() == ddata_.size() && btl.fdata_.size() == fdata_.size() && btl.idata_.size() == idata_.size() && btl.cdata_.size() == cdata_.size()) return true;
     30     else return false;
     31   }
     32
     33void BnTblLine::Print()
     34   {
     35     int k;
     36     cout << " ********* ligne ************* " << endl;
     37     cout << " *** noms de variables  " << endl;
     38     for (k=0; k < ColName_.size(); k++) cout << ColName_[k] << " ";
     39     cout << endl;
     40     cout << " *** variables doubles  " << endl;
     41     for (k=0; k < ddata_.size(); k++) cout << ddata_[k] << " ";
     42     cout << endl;
     43     cout << " *** variables float  " << endl;
     44     for (k=0; k < fdata_.size(); k++) cout << fdata_[k] << " ";
     45     cout << endl;
     46     cout << " *** variables int  " << endl;
     47     for (k=0; k < idata_.size(); k++) cout << idata_[k] << " ";
     48     cout << endl;
     49     cout << " *** variables string  " << endl;
     50     for (k=0; k < cdata_.size(); k++) cout << cdata_[k] << " ";
     51     cout << endl;
     52     cout << " ***************************** " << endl;
     53   }
     54
     55
     56
     57/*!
     58  \class SOPHYA::FitsIOHandler
     59The class structure is analogous to Sophya-PPersist system :
     60Each SOPHYA object XXX is associated with a object of class FITS_XXX
     61 (inheriting from FitsFileHandler), to which input/output operations with FITS
     62 files are delegated (through a class Hierarchy : FitsFile (virtual),
     63 FitsInFile, FitsOutFile) . A typical example of use is the following :
     64
     65\verbatim
     66  int m=... ;
     67  SphereHEALPix<r_8> sphere1(m);           // definition of the SOPHYA object
     68  .... fill the sphere ....
     69
     70  FITS_SphereHEALPix<r_8> fits_sph1(sphere1);
     71                                           // delegated object
     72  fits_sph.Write("myfile.fits");           // writing on FITS file
     73
     74   FITS_SphereHEALPix<r_8> fits_sph2("myfile.fits");
     75                                           // load a delegated object
     76                                           // from FITS file
     77   SphereHEALPix<r_8> sphere2=(SphereHEALPix<r_8>)fits_sph2;
     78                                           // casting the delegated object
     79                                           // into a SOPHYA object
     80\endverbatim
     81 
     82
     83*/
     84
     85/*! \fn void  SOPHYA::FitsIOHandler::Read(char flnm[],int hdunum)
     86
     87this method is called from inherited objects :
     88
     89opens a file 'flnm'
     90
     91gets parameters in extension-header (hdunum)
     92
     93calls the method 'ReadFromFits' from the inherited  object
     94*/
    1095void   FitsIOHandler::Read(char flnm[],int hdunum)
    1196{
     
    1398  Read(ifts, hdunum);
    1499}
     100
     101  /*! \fn void SOPHYA::FitsIOHandler::Read(FitsInFile& is, int hdunum)
     102Read the data on extension hdunum (or primary header, if hdunum=1) from FitsInFIle. If hdunum is not addressed, , one reads the next extension, with respect to the current position.
     103   */
    15104void FitsIOHandler::Read(FitsInFile& is, int hdunum)
    16105{
     
    21110
    22111
     112/*! \fn void SOPHYA::FitsIOHandler::Write(char flnm[])
     113this method is called from inherited objects.
     114
     115for writing a new object in a new fits-extension :
     116
     117\warning By convention, primary header does not contain fits-image data : i.e.
     118all data are fits-extensions. The first relevant header will have hdunum=2.
     119For switching off this convention use the method :
     120
     121firstImageOnPrimaryHeader() (see below)
     122
     123In that case do not forget to precise hdunum=1 when reading data on primary header.
     124
     125calls the method 'WriteToFits' from the inherited  object
     126
     127*/
    23128void FitsIOHandler::Write(char flnm[])
    24129
     
    33138}
    34139
     140
     141/*!
     142  \class SOPHYA::FitsIOHandler
     143Class (virtual) for managing FITS format files
     144*/
    35145
    36146
     
    87197}
    88198
     199/*!
     200  \class SOPHYA::FitsInFile
     201
     202class for saving  SOPHYA objects on FITS Format Files (uses cfitsio lib)
     203*/
     204
    89205FitsInFile::FitsInFile()
    90206{
    91207  InitNull();
    92208}
    93 //FitsInFile::FitsInFile(char flnm[], int hdunum)
     209
    94210FitsInFile::FitsInFile(char flnm[])
    95211{
     
    118234}
    119235
     236//////////////////////////////////////////////////////////
     237//     methods with general purpose
     238/////////////////////////////////////////////////////////
    120239
    121240int FitsInFile::NbBlocks(char flnm[])
     
    200319}
    201320
     321
     322void FitsInFile::ReadFInit(int hdunum)
     323{
     324   InitNull();
     325  // int status = 0;
     326 
     327  //  fits_open_file(&fptr_,flnm,READONLY,&status);
     328  //  if( status ) printerror( status );
     329
     330  int status = 0;
     331
     332  if (hdunum <= 1)
     333    {
     334      hdunum_ = 1;
     335      // presence of image ?
     336      int naxis= 0;
     337      fits_read_key(fptr_,TINT,"NAXIS",&naxis,NULL,&status);
     338      if( status ) printerror( status );
     339      if (naxis > 0 )       // there is an image
     340        {
     341          hdutype_ = IMAGE_HDU;
     342          GetImageParameters (fptr_, bitpix_, naxis_, naxisn_);
     343          nbData_ =  1;
     344          int k;
     345          for (k=0; k<naxis_; k++) if (naxisn_[k] > 0) nbData_ *= naxisn_[k];
     346          KeywordsIntoDVList(fptr_, dvl_,hdunum_);
     347        }
     348      else
     349        {
     350          throw PException(" first header : no image, probably error in hdunum");
     351        }
     352  //
     353    }
     354  else
     355    {
     356      hdunum_ = hdunum-1;
     357      int hdutype;
     358      fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
     359      if( status ) printerror( status,":FitsFile::ReadF : erreur movabs");
     360      moveToFollowingHeader();
     361    }
     362}
     363
     364
    202365void FitsInFile::GetImageParameters (fitsfile* fileptr,int& bitpix,int& naxis,vector<int>& naxisn)
    203366{
     
    230393
    231394
    232 void FitsInFile::ReadFInit(int hdunum)
    233 {
    234    InitNull();
    235   // int status = 0;
    236  
    237   //  fits_open_file(&fptr_,flnm,READONLY,&status);
    238   //  if( status ) printerror( status );
    239 
    240   int status = 0;
    241 
    242   if (hdunum <= 1)
    243     {
    244       hdunum_ = 1;
    245       // presence of image ?
    246       int naxis= 0;
    247       fits_read_key(fptr_,TINT,"NAXIS",&naxis,NULL,&status);
    248       if( status ) printerror( status );
    249       if (naxis > 0 )       // there is an image
    250         {
    251           hdutype_ = IMAGE_HDU;
    252           GetImageParameters (fptr_, bitpix_, naxis_, naxisn_);
    253           nbData_ =  1;
    254           int k;
    255           for (k=0; k<naxis_; k++) if (naxisn_[k] > 0) nbData_ *= naxisn_[k];
    256           KeywordsIntoDVList(fptr_, dvl_,hdunum_);
    257         }
    258       else
    259         {
    260           throw PException(" first header : no image, probably error in hdunum");
    261         }
    262   //
    263     }
    264   else
    265     {
    266       hdunum_ = hdunum-1;
    267       int hdutype;
    268       fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
    269       if( status ) printerror( status,":FitsFile::ReadF : erreur movabs");
    270       moveToFollowingHeader();
    271     }
    272 }
    273 
     395
     396
     397  /*! \fn DVList SOPHYA::FitsInFile::DVListFromPrimaryHeader() const
     398
     399   \return the keywords of primary header in a DVList
     400
     401*/
    274402DVList  FitsInFile::DVListFromPrimaryHeader() const
    275403   {
     
    308436
    309437
     438
     439
     440
     441/*! \fn int  SOPHYA::FitsInFile::NbColsFromFits() const
     442\return number of columns (return 1 if IMAGE)
     443*/
    310444int  FitsInFile::NbColsFromFits() const
    311445{
     
    320454}
    321455
     456/*! \fn  int SOPHYA::FitsInFile::NentriesFromFits(int nocol) const
     457\return number of data in the current IMAGE extension on FITS file, or number
     458 of data of column number 'nocol' of the current BINTABLE extension
     459*/
    322460int FitsInFile::NentriesFromFits(int nocol) const
    323461{
     
    334472}
    335473
     474/*! \fn char SOPHYA::FitsInFile::ColTypeFromFits(int nocol) const
     475
     476return a character denoting data type of column number 'nocol' in a BINTABLE :
     477
     478D : double
     479
     480E : float
     481
     482I : integer
     483   
     484S : character string
     485
     486  */
     487
    336488char FitsInFile::ColTypeFromFits(int nocol) const
    337489{
     
    342494  return types_[nocol];
    343495}
     496
     497
     498/*! \fn string SOPHYA::FitsInFile::ColNameFromFits(int nocol) const
     499
     500\return name of the column number 'nocol' of the current BINTABLE extension
     501   */
     502
    344503string FitsInFile::ColNameFromFits(int nocol) const
    345504{
     
    350509  return noms_[nocol];
    351510}
     511
     512/*! \fn int DSOPHYA::FitsInFile::ColStringLengthFromFits(int nocol) const
     513
     514 \return number of characters of each data  for the column number 'nocol' (if char* typed) of the current BINTABLE extension
     515*/
    352516
    353517int FitsInFile::ColStringLengthFromFits(int nocol) const
     
    365529  return  taille_des_chaines_[index];
    366530}
     531
     532
     533
     534/*! \fn void  SOPHYA::FitsInFile::GetBinTabLine(int NoLine, double* ddata, float* fdata, int* idata, char ** cdata)
     535
     536Get the NoLine-th 'line'  from the current BINTABLE extension on FITS file,
     537  */
     538
    367539void  FitsInFile::GetBinTabLine(int NoLine, double* ddata, float* fdata, int* idata, char ** cdata)
    368540{
     
    405577}
    406578
     579/*! \fn void   SOPHYA::FitsInFile::GetBinTabLine(long NoLine,  BnTblLine& ligne)
     580Get the NoLine-th 'line'  from the current BINTABLE extension on FITS file,
     581*/
    407582void   FitsInFile::GetBinTabLine(long NoLine,  BnTblLine& ligne)
    408583{
     
    447622}
    448623
     624/*! \fn void SOPHYA::FitsInFile::GetBinTabLine(int NoLine, float* fdata)
     625
     626Get the NoLine-th float 'line'  from the current BINTABLE extension on FITS file,
     627*/
    449628void FitsInFile::GetBinTabLine(int NoLine, float* fdata)
    450629{
     
    466645
    467646
     647/*! \fn void SPOPHYA::FitsInFile::GetBinTabFCol(double* valeurs,int nentries, int NoCol) const
     648
     649fill the array 'valeurs' with double data from the current BINTABLE extension on FITS file, from column number 'NoCol'
     650
     651\param <nentries>  number of data to be read
     652*/
    468653void FitsInFile::GetBinTabFCol(double* valeurs,int nentries, int NoCol) const
    469654    {
     
    499684    }
    500685
     686/*! \fn  void SOPHYA::FitsInFile::GetBinTabFCol(float* valeurs,int nentries, int NoCol) const
     687
     688 same as previous method with float data
     689*/
    501690void FitsInFile::GetBinTabFCol(float* valeurs,int nentries, int NoCol) const
    502691    {
     
    531720    }
    532721
     722/*! \fn void SOPHYA::FitsInFile::GetBinTabFCol(int* valeurs,int nentries, int NoCol) const
     723
     724 same as previous method with int data
     725*/
     726
    533727void FitsInFile::GetBinTabFCol(int* valeurs,int nentries, int NoCol) const
    534728    {
     
    561755    }
    562756
     757/*! \fn void SOPHYA::FitsInFile::GetBinTabFCol(char** valeurs, int nentries, int NoCol) const
     758
     759 same as previous method with char* data
     760*/
     761
    563762void FitsInFile::GetBinTabFCol(char** valeurs, int nentries, int NoCol) const
    564763    {
     
    587786    }
    588787
     788/*! \fn void SOPHYA::FitsInFile::GetSingleColumn(double* map, int nentries) const
     789fill the array 'map' with double data from the current extension on FITS file.
     790If the extension is BINTABLE, the first column is provided.
     791
     792\param <nentries>  number of data to be read
     793*/
    589794void FitsInFile::GetSingleColumn(double* map, int nentries) const
    590795{
     
    619824}
    620825
     826/*! \fn void SOPHYA::FitsInFile::GetSingleColumn(float* map, int nentries) const
     827same as above with float data
     828*/
    621829void FitsInFile::GetSingleColumn(float* map, int nentries) const
    622830{
     
    649857}
    650858
     859/*! \fn void SOPHYA::FitsInFile::GetSingleColumn( int* map, int nentries) const
     860 same as above with int data
     861*/
    651862void FitsInFile::GetSingleColumn( int* map, int nentries) const
    652863{
     
    8661077}
    8671078
     1079
     1080/*!
     1081  \class SOPHYA::FitsOutFile
     1082 Class for loading  SOPHYA objects from FITS Format Files (uses cfitsio lib)
     1083*/
     1084
    8681085FitsOutFile::FitsOutFile()
    8691086{
     
    8711088}
    8721089
     1090   /*! \fn SOPHYA::FitsOutFile::FitsOutFile(char flnm[], WriteMode wrm)
     1091
     1092\param <WriteMode>  enum , WriteMode = clear -> if alreadyy exists, the file will be overwritten (else created) ; WriteMode = append -> further objects will be appended to the file if it exists (else : file created). WriteMode = unknown -> file created if does not exist, else : exception. (the last situation is the default)
     1093
     1094   */
    8731095FitsOutFile::FitsOutFile(char flnm[], WriteMode wrm)
    8741096{
     
    9291151
    9301152
     1153/*! \fn void SOPHYA::FitsOutFile::makeHeaderImageOnFits(char type, int nbdim, int* naxisn,  DVList &dvl)
     1154
     1155create an IMAGE header on FITS file.
     1156\param <type> type of data (see method ColTypeFromFits)
     1157\param <nbdim>  number of dimensions : 1D, 2D, 3D etc. = NAXIS
     1158\param <naxisn>  array containind sizes of the different dimensions
     1159*/
    9311160void FitsOutFile::makeHeaderImageOnFits(char type, int nbdim, int* naxisn,  DVList &dvl)
    9321161{
     
    9681197
    9691198}
     1199
     1200
     1201/*! \fn void SOPHYA::FitsOutFile::PutImageToFits(int nbData, double* map) const
     1202
     1203write double data from array 'map'on an IMAGE extension
     1204\param <nbData>  number of data to be written
     1205*/
    9701206void FitsOutFile::PutImageToFits(int nbData, double* map) const
    9711207{
     
    9771213}
    9781214
     1215/*! \fn void SOPHYA::FitsOutFile::PutImageToFits(int nbData, float* map) const
     1216
     1217same as previous method with float data
     1218*/
    9791219void FitsOutFile::PutImageToFits(int nbData, float* map) const
    9801220{
     
    9861226
    9871227}
     1228
     1229  /*! \fn void SOPHYA::FitsOutFile::PutImageToFits( int nbData, int* map) const
     1230
     1231 same as previous method with int data */
    9881232void FitsOutFile::PutImageToFits( int nbData, int* map) const
    9891233{
     
    9981242
    9991243
     1244/*! \fn void SOPHYA::FitsOutFile::makeHeaderBntblOnFits( string fieldType, vector<string> Noms, int nentries, int tfields, DVList &dvl, string extname, vector<int> taille_des_chaines)
     1245
     1246create an BINTABLE header on FITS file.
     1247\param <fieldType> array conta
     1248ining characters denoting types of the different column (see method ColTypeFromFits)
     1249\param <Noms>  array of the names of columns
     1250\param <nentries>  number of data of each column
     1251\param <tfields> number of columns
     1252\param <dvl> a SOPHYA DVList containing keywords to be appended
     1253\param <extname> keyword EXTNAME for FITS file
     1254\param <taille_des_chaines> vector containing the number of characters of  data  for each char* typed column, with order of appearance in 'fieldType'
     1255*/
    10001256void FitsOutFile::makeHeaderBntblOnFits( string fieldType, vector<string> Noms, int nentries, int tfields, DVList &dvl, string extname, vector<int> taille_des_chaines)
    10011257{
     
    10921348}
    10931349
     1350/*! \fn void SOPHYA::FitsOutFile::PutColToFits(int nocol, int nentries, double* donnees) const
     1351
     1352write double data from array 'donnees ' on column number 'nocol' of a BINTABLE  extension.
     1353\param <nentries>  number of data to be written
     1354*/
    10941355void FitsOutFile::PutColToFits(int nocol, int nentries, double* donnees) const
    10951356{
     
    11141375  if( status )  printerror( status,"erreur ecriture du fichier fits" );
    11151376}
     1377
     1378
     1379
     1380  /*! \fn void SOPHYA::FitsOutFile::PutColToFits(int nocol, int nentries, float* donnees) const
     1381
     1382same as previous method with float data
     1383*/
    11161384void FitsOutFile::PutColToFits(int nocol, int nentries, float* donnees) const
    11171385{
     
    11401408  if( status )  printerror( status,"erreur ecriture du fichier fits" );
    11411409}
     1410
     1411
     1412/*! \fn void FitsOutFile::PutColToFits(int nocol, int nentries, int* donnees) const
     1413
     1414same as previous method with int data
     1415*/
    11421416void FitsOutFile::PutColToFits(int nocol, int nentries, int* donnees) const
    11431417{
     
    11661440  if( status )  printerror( status," ecriture du fichier fits" );
    11671441}
     1442
     1443
     1444/*! \fn void SOPHYA::FitsOutFile::PutColToFits(int nocol, int nentries, char** donnees) const
     1445same as previous method with char* data
     1446*/
    11681447void FitsOutFile::PutColToFits(int nocol, int nentries, char** donnees) const
    11691448{
     
    12281507
    12291508
     1509/* \fn void  SOPHYA::FitsOutFile::DVListIntoPrimaryHeader(DVList& dvl) const
     1510
     1511Put keywords from a DVList into the primary header of the fits-file
     1512*/
    12301513void  FitsOutFile::DVListIntoPrimaryHeader(DVList& dvl) const
    12311514{
  • trunk/SophyaExt/FitsIOServer/fitsfile.h

    r1209 r1218  
    1010#define LEN_KEYWORD 9
    1111
     12// classes for saving/loading SOPHYA objects to/from FITS files...
     13// Guy le Meur (september 2000)
     14
     15
    1216namespace SOPHYA {
    1317
    14 struct BnTblLine
    15  {
    16 BnTblLine() {}
    17 void setFormat(int dc, int fc, int ic, int cc, vector<string> names)
    18    {
    19     int nbcols = dc + fc + ic + cc;
    20     int maxName = names.size();
    21     if (nbcols != maxName)
    22       {
    23         cout << " WARNING: BnTblLine:: length of vector of column names not equal to total number of columns" << endl;
    24         maxName = nbcols < maxName ? nbcols : maxName;
    25       }
    26     ColName_ = vector<string>(nbcols);
    27      for (int k=0; k < maxName; k++) ColName_[k] = names[k];
    28     if (dc >0) ddata_ = vector<double>(dc);
    29     if (fc >0) fdata_ = vector<float>(fc);
    30     if (ic >0) idata_ = vector<int>(fc);
    31     if (cc >0) cdata_ = vector<string>(fc);
    32    }
    33 
    34 bool sameFormat(const BnTblLine& btl) const
    35    {
    36      if (btl.ddata_.size() == ddata_.size() && btl.fdata_.size() == fdata_.size() && btl.idata_.size() == idata_.size() && btl.cdata_.size() == cdata_.size()) return true;
    37      else return false;
    38    }
    39 
    40 void Print()
    41    {
    42      int k;
    43      cout << " ********* ligne ************* " << endl;
    44      cout << " *** noms de variables  " << endl;
    45      for (k=0; k < ColName_.size(); k++) cout << ColName_[k] << " ";
    46      cout << endl;
    47      cout << " *** variables doubles  " << endl;
    48      for (k=0; k < ddata_.size(); k++) cout << ddata_[k] << " ";
    49      cout << endl;
    50      cout << " *** variables float  " << endl;
    51      for (k=0; k < fdata_.size(); k++) cout << fdata_[k] << " ";
    52      cout << endl;
    53      cout << " *** variables int  " << endl;
    54      for (k=0; k < idata_.size(); k++) cout << idata_[k] << " ";
    55      cout << endl;
    56      cout << " *** variables string  " << endl;
    57      for (k=0; k < cdata_.size(); k++) cout << cdata_[k] << " ";
    58      cout << endl;
    59      cout << " ***************************** " << endl;
    60    }
    61   vector<double> ddata_;
    62   vector<float>  fdata_;
    63   vector<int>    idata_;
    64   vector<string>  cdata_;
    65   vector<string> ColName_;
    66 
    67  };
    68 
    69  class FitsFile;
    70  class FitsInFile;
    71  class FitsOutFile;
    72  enum WriteMode {append, clear, unknown};
    73 
     18  struct BnTblLine;
     19  class FitsFile;
     20  class FitsInFile;
     21  class FitsOutFile;
     22  enum WriteMode {append, clear, unknown};
     23 
    7424
    7525//
    7626//! Class for managing Interface for SOPHYA objects to FITS Format Files (uses cfitsio lib)
    7727
    78 /*!
    79 The class structure is analogous to Sophya-PPersist system :
    80 Each SOPHYA object XXX is associated with a object of class FITS_XXX
    81  (inheriting from FitsFileHandler), to which input/output operations with FITS
    82  files are delegated (through a class Hierarchy : FitsFile (virtual),
    83  FitsInFile, FitsOutFile) . A typical example of use is the following :
    84 
    85 \verbatim
    86   int m=... ;
    87   SphereHEALPix<r_8> sphere1(m);           // definition of the SOPHYA object
    88   .... fill the sphere ....
    89 
    90   FITS_SphereHEALPix<r_8> fits_sph1(sphere1);
    91                                            // delegated object
    92   fits_sph.Write("myfile.fits");           // writing on FITS file
    93 
    94    FITS_SphereHEALPix<r_8> fits_sph2("myfile.fits");
    95                                            // load a delegated object
    96                                            // from FITS file
    97    SphereHEALPix<r_8> sphere2=(SphereHEALPix<r_8>)fits_sph2;
    98                                            // casting the delegated object
    99                                            // into a SOPHYA object
    100 \endverbatim
    101 
    102 */
    103 
    104 class FitsIOHandler {
    105 
    106 
    107   public:
    108 
    109 virtual ~FitsIOHandler() {}
    110 /*!
    111 this method is called from inherited objects :
    112 
    113 opens a file 'flnm'
    114 
    115 gets parameters in extension-header (hdunum)
    116 
    117 calls the method 'ReadFromFits' from the inherited  object
    118 
    119 
    120 */
    121   void   Read(char flnm[],int hdunum= 0);
    122 /*!
    123 this method is called from inherited objects :
    124 
    125 for writing a new object in a new fits-extension :
    126 
    127 ???
    128 
    129  at the end of
    130 
    131 the existing file (flnm), if OldFile=true.
    132 
    133 If OldFile=false, an exception occurs
     28 class FitsIOHandler {
     29
     30
     31 public:
     32
     33   virtual ~FitsIOHandler() {}
     34   void   Read(char flnm[],int hdunum= 0);
     35   void   Write(char flnm[]) ;
     36   void   Read(FitsInFile& ifts, int hdunum=0);
     37   void   Write(FitsOutFile& ofts) ;
     38
     39
     40 protected:
     41 
     42   virtual void    ReadFromFits(FitsInFile& is)=0;           
     43   virtual void    WriteToFits(FitsOutFile& os) =0;   
     44       
     45   friend class FitsInFile;
     46   friend class FitsOutFile;
     47  };
     48
     49
     50//! Class (virtual) for managing FITS format files
     51 class FitsFile {
     52
     53 public:
     54
     55   FitsFile() { InitNull(); };
     56   virtual ~FitsFile();
     57   static string GetErrStatus(int status);
     58   inline  int   statusF() const { return fits_status_;}
    13459 
    135 By convention, primary header does not contain fits-image data : i.e.
    136 all data are fits-extensions. The first relevant header will have hdunum=2.
    137 For switching off this convention use the method :
    138 
    139 firstImageOnPrimaryHeader() (see below)
    140 
    141 In that case do not forget to precise hdunum=1 when reading data on primary header.
    142 
    143 calls the method 'WriteToFits' from the inherited  object
    144 
    145 */
    146   void   Write(char flnm[]) ;
    147 
    148   /*!
    149 Read the data on extension hdunum (or primary header, if hdunum=1) from FitsInFIle. With default value for hdunum, one reads the next extension, with respect to the current position.
    150    */
    151   void   Read(FitsInFile& ifts, int hdunum=0);
    152   void   Write(FitsOutFile& ofts) ;
    153 
    154 
    155   protected: 
    156   virtual void    ReadFromFits(FitsInFile& is)=0;           
    157   virtual void    WriteToFits(FitsOutFile& os) =0;           
    158   friend class FitsInFile;
    159   friend class FitsOutFile;
    160   private :
    161   };
    162 
    163 
    164 
    165 class FitsFile
    166 {
    167 
    168 public:
    169 
    170 FitsFile()
    171 {
    172  InitNull();
    173 };
    174   virtual ~FitsFile();
    175 
    176   static string GetErrStatus(int status);
    177 
    178 
    179 
    180 
    181 inline  int statusF() const { return fits_status_;}
    182  
    183 
    184 protected:
    185 
    186    void ResetStatus(int& status) ;
     60
     61 protected:
     62
     63   void         ResetStatus(int& status) ;
    18764   static  void printerror(int&) ;
    18865   static  void printerror(int&,char* texte) ;
    189    inline void InitNull() {fptr_ = NULL; hdutype_= 0; hdunum_ = 1;
     66   inline void  InitNull() {fptr_ = NULL; hdutype_= 0; hdunum_ = 1;
    19067   fits_status_ = 0;}
    19168
    192   //! pointer to the FITS file, defined in fitsio.h
    193   fitsfile *fptr_;
    194  
    195   //!  image or bintable ?
    196   int hdutype_;
    197 
    198 //! index of header to be read/written
    199   int hdunum_; 
    200 
    201   //! last status returned by fitsio library. updated only by several methods
    202   int fits_status_;
    203 
    204    
    205 };
    206 
     69   fitsfile *fptr_;     /**<  pointer to the FITS file, defined in fitsio.h */
     70   int hdutype_;        /**<  image or bintable ? */
     71   int hdunum_;         /**<   index of header to be read/written */
     72   int fits_status_;    /**< last status returned by fitsio library. updated only by several methods */
     73
     74 };
     75
     76//! Class for saving  SOPHYA objects on FITS Format Files (uses cfitsio lib)
    20777
    20878 class FitsInFile : public  FitsFile {
     
    21080 public:
    21181   FitsInFile();
    212    //   FitsInFile(char flnm[], int hdunum=0);
    21382   FitsInFile(char flnm[]);
    21483   ~FitsInFile() { ; };
    21584
    216 
    217 //////////////////////////////////////////////////////////
    218 //     methods with general purpose
    219 ///////////////////////////////////////
    220 
    221 
    222 
    223   static int NbBlocks(char flnm[]);
    224   static void GetBlockType(char flnm[], int hdunum, string& typeOfExtension, int& naxis, vector<int>& naxisn, string& dataType, DVList& dvl  );
    225 
    226 
    227 
    228   //  void ReadFInit(char flnm[],int hdunum=0);
    229   void ReadFInit(int hdunum);
     85   static int  NbBlocks(char flnm[]);
     86   static void GetBlockType(char flnm[], int hdunum, string& typeOfExtension, int& naxis, vector<int>& naxisn, string& dataType, DVList& dvl  );
     87   void        ReadFInit(int hdunum);
    23088 
    231   /*! return a reference on a DVList containing the keywords from FITS file
    232    */
    233 inline const DVList& DVListFromFits() const { return dvl_;}
    234 
    235 /* get the keywords of primary header in a DVList */
    236 DVList  DVListFromPrimaryHeader() const;
    237 
    238 void moveToFollowingHeader();
     89  /*! \return a reference on a DVList containing the keywords from FITS file */
     90  inline const DVList& DVListFromFits() const { return dvl_;}
     91
     92  DVList  DVListFromPrimaryHeader() const;
     93  void    moveToFollowingHeader();
     94
     95
    23996
    24097
     
    24299       ///////   methods for managing extensions ////////////////
    243100       //////////////////////////////////////////////////////////
    244 
    245101
    246102
     
    250106
    251107
    252 /*! return true if the current header  corresponds to a FITS image extension */
     108
     109/*! \return true if the current header  corresponds to a FITS image extension */
    253110inline bool IsFitsImage() const { return (hdutype_ == IMAGE_HDU);}
    254111
    255112
    256113
    257   /*! number of dimensions of an image extension : NAXIS parameter (in FITS notations)
    258    */
     114  /*! \return number of dimensions of an image extension : NAXIS parameter (in FITS notations)   */
    259115inline int nbDimOfImage() const {return naxis_;}
    260 /*! a reference on a vector containing sizes of the NAXIS dimensions : NAXIS1, NAXIS2, NAXIS3 wtc.
    261  */
     116
     117/*! \return a reference on a vector containing sizes of the NAXIS dimensions : NAXIS1, NAXIS2, NAXIS3 etc. */
    262118 inline const vector<int>& dimOfImageAxes() const { return naxisn_;}
    263 /*!
    264  total number of data in the current IMAGE extension
    265  */
     119
     120
     121/*! \return total number of data in the current IMAGE extension */
    266122inline int nbOfImageData() const { return nbData_; }
    267123
     
    275131
    276132
    277 /*! return true if the current header  corresponds to a FITS ASCII or BINTABLE extension */
     133/*! \return true if the current header  corresponds to a FITS ASCII or BINTABLE extension */
    278134inline bool IsFitsTable() const {return (hdutype_ == ASCII_TBL || hdutype_ == BINARY_TBL);}
    279135
    280136
    281 static  void GetBinTabParameters(fitsfile* fileptr, int& nbcols, int& nrows,
     137 static  void GetBinTabParameters(fitsfile* fileptr, int& nbcols, int& nrows,
    282138                                  vector<int>& repeat,
    283139                                  vector<string>& noms,
    284140                                  vector<char>& types,   
    285141                                  vector<int>&  taille_des_chaines);
    286 
    287 
    288 
    289   /*! return a character denoting data type of column number 'nocol' in a BINTABLE :
    290 
    291 D : double
    292 
    293 E : float
    294 
    295 I : integer
    296    
    297 S : character string
    298 
    299   */
    300   char    ColTypeFromFits(int nocol) const;
    301   /*! name of the column number 'nocol' of the current BINTABLE extension
    302    */
    303   string  ColNameFromFits(int nocol) const;
    304 
    305   /*! number of characters of each data  for the column number 'nocol' (if char* typed) of the current BINTABLE extension
    306    */
    307   int     ColStringLengthFromFits(int nocol) const;
    308 
    309 
    310 
    311   /*!
    312 Get the NoLine-th 'line'  from the current BINTABLE extension on FITS file,
    313   */
    314   void GetBinTabLine(int NoLine, double* ddata, float* fdata, int* idata, char
     142 char   ColTypeFromFits(int nocol) const;
     143 string ColNameFromFits(int nocol) const;
     144 int    ColStringLengthFromFits(int nocol) const;
     145 void   GetBinTabLine(int NoLine, double* ddata, float* fdata, int* idata, char
    315146** cdata) ;
    316   /*!
    317 Get the NoLine-th 'line'  from the current BINTABLE extension on FITS file,
    318   */
    319   void GetBinTabLine(long NoLine, BnTblLine& ligne) ;
    320 
    321   /*!
    322 Get the NoLine-th 'line'  from the current BINTABLE extension on FITS file,
    323   */
    324   void GetBinTabLine(int NoLine, float* fdata) ;
    325 
    326 /*!
    327 fill the array 'valeurs' with double data from the current BINTABLE extension on FITS file, from column number 'NoCol'
    328 
    329 \param <nentries>  number of data to be read
    330   */
    331   void GetBinTabFCol(double* valeurs, int nentries, int NoCol) const;
    332 
    333   /*! same as previous method with float data */
    334   void GetBinTabFCol(float* valeurs, int nentries, int NoCol) const;
    335   /*! same as previous method with int data */
    336   void GetBinTabFCol(int* valeurs, int nentries,  int NoCol) const;
    337   /*! same as previous method with char* data */
    338   void GetBinTabFCol(char** valeurs,int nentries, int NoCol) const;
    339   // Write elements into the FITS data array
     147 void   GetBinTabLine(long NoLine, BnTblLine& ligne) ;
     148 void   GetBinTabLine(int NoLine, float* fdata) ;
     149 void   GetBinTabFCol(double* valeurs, int nentries, int NoCol) const;
     150 void   GetBinTabFCol(float* valeurs, int nentries, int NoCol) const;
     151 void   GetBinTabFCol(int* valeurs, int nentries,  int NoCol) const;
     152 void   GetBinTabFCol(char** valeurs,int nentries, int NoCol) const;
    340153
    341154/////////////////////////////////////////////////////////////
     
    343156////////////////////////////////////////////////////////
    344157
    345  /*! return number of columns (return 1 if IMAGE) */
    346158  int     NbColsFromFits() const;
    347   /*! number of data in the current IMAGE extension on FITS file, or number
    348  of data of column number 'nocol' of the current BINTABLE extension
    349   */
    350159  int     NentriesFromFits(int nocol) const;
    351160
    352161
    353 /*!
    354 fill the array 'map' with double data from the current extension on FITS file.
    355 If the extension is BINTABLE, the first column is provided.
    356 
    357 \param <nentries>  number of data to be read
    358   */
    359162  void    GetSingleColumn(double* map, int nentries) const;
    360163
    361   /*! same as above with float data */
    362164  void    GetSingleColumn(float*  map, int nentries) const;
    363165
    364   /*! same as above with int data */
    365166  void    GetSingleColumn(int* map, int nentries) const;
    366167
     
    371172static  void GetImageParameters (fitsfile* fileptr,int& bitpix,int& naxis,vector<int>& naxisn);
    372173
    373 
    374  
    375 
    376   //! fits-Image parameter
    377   int bitpix_;
    378 
    379   //! fits-Image parameter
    380   int naxis_;
    381 
    382   //! fits-Image parameters : sizes of dimensions
    383   vector<int> naxisn_;
    384 
    385   //! fits-Image parameter: number of data
    386   int nbData_;
    387 
    388   //! Bintable parameter
    389   int nrows_;
    390 
    391   //! Bintable parameter
    392   vector<int> repeat_;
    393 
    394   //! Bintable parameter
    395   int nbcols_;
    396 
    397   //! Bintable parameter: column names
    398   vector<string> noms_;
    399 
    400   //! Bintable parameters: types of columns (D: double, E: float, I: integers,  A: char*)
    401   vector<char> types_;   
    402 
    403   //! Bintable parameters:   length of the char* variables                 
    404   vector<int>  taille_des_chaines_;
    405 
    406   //! DVList for transferring keywords
    407   DVList dvl_;
    408 
    409 
    410  };
    411 
     174 int bitpix_;          /**< fits-Image parameter */
     175 int naxis_;           /**< fits-Image parameter */
     176 vector<int> naxisn_;  /**< fits-Image parameters : sizes of dimensions */
     177 int nbData_;          /*< fits-Image parameter: number of data */
     178 int nrows_;           /**< Bintable parameter */
     179 vector<int> repeat_;  /**< Bintable parameter */
     180 int nbcols_;          /**< Bintable parameter */
     181 vector<string> noms_; /**< Bintable parameter: column names */
     182 vector<char> types_;  /**< Bintable parameters: types of columns (D: double, E: float, I: integers,  A: char*) */
     183 DVList dvl_;          /**< DVList for transferring keywords */
     184 vector<int>  taille_des_chaines_; /**< Bintable parameters:   length of the char* variables */
     185
     186 };
     187
     188//! Class for loading  SOPHYA objects from FITS Format Files (uses cfitsio lib)
    412189
    413190 class FitsOutFile : public  FitsFile {
     
    417194
    418195   FitsOutFile();
    419    /*!
    420 \param <WriteMode>  enum , WriteMode = clear -> if alreadyy exists, the file will be overwritten (else created) ; WriteMode = append -> further objects will be appended to the file if it exists (else : file created). WriteMode = unknown -> file created if does not exist, else : exception. (the last situation is the default)
    421 
    422    */
    423196   FitsOutFile(char flnm[], WriteMode wrm = unknown );
    424197   ~FitsOutFile() { ;};
     
    437210
    438211   inline void firstImageOnPrimaryHeader() {imageOnPrimary_=true;}
    439 
    440   /*! create an IMAGE header on FITS file.
    441 \param <type> type of data (see method ColTypeFromFits)
    442 \param <nbdim>  number of dimensions : 1D, 2D, 3D etc. = NAXIS
    443 \param <naxisn>  array containind sizes of the different dimensions
    444   */
    445   void makeHeaderImageOnFits(char type, int nbdim, int* naxisn,  DVList &dvl) ;
    446 
    447 
    448   /*! write double data from array 'map'on an IMAGE extension
    449 \param <nbData>  number of data to be written
    450 
    451    */
    452   void PutImageToFits( int nbData, double* map) const;
    453 
    454   /*! same as previous method with float data */
    455   void PutImageToFits(int nbData, float* map ) const;
    456 
    457   /*! same as previous method with int data */
    458   void PutImageToFits(int nbData, int* map) const;
     212   void makeHeaderImageOnFits(char type, int nbdim, int* naxisn, DVList &dvl) ;
     213   void PutImageToFits( int nbData, double* map) const;
     214   void PutImageToFits(int nbData, float* map ) const;
     215   void PutImageToFits(int nbData, int* map) const;
    459216
    460217
     
    466223
    467224
    468   /*! create an BINTABLE header on FITS file.
    469 \param <fieldType> array conta
    470 ining characters denoting types of the different column (see method ColTypeFromFits)
    471 \param <Noms>  array of the names of columns
    472 \param <nentries>  number of data of each column
    473 \param <tfields> number of columns
    474 \param <dvl> a SOPHYA DVList containing keywords to be appended
    475 \param <extname> keyword EXTNAME for FITS file
    476 \param <taille_des_chaines> vector containing the number of characters of  data  for each char* typed column, with order of appearance in 'fieldType'
    477    */
    478225   void makeHeaderBntblOnFits ( string fieldType, vector<string> Noms, int nentries, int tfields, DVList &dvl, string extname,  vector<int> taille_des_chaines) ;
    479 
    480   /*! write double data from array 'donnees ' on column number 'nocol' of a BINTABLE  extension.
    481 \param <nentries>  number of data to be written
    482 
    483    */
    484   void PutColToFits(int nocol, int nentries, double* donnees) const;
    485 
    486   /*! same as previous method with float data */
    487   void PutColToFits(int nocol, int nentries, float* donnees) const;
    488 
    489   /*! same as previous method with int data */
    490   void PutColToFits(int nocol, int nentries, int* donnees) const;
    491 
    492   /*! same as previous method with char* data */
    493   void PutColToFits(int nocol, int nentries, char** donnees) const;
    494 
    495   void PutBinTabLine(long NoLine,  BnTblLine& ligne) const;
     226   void PutColToFits(int nocol, int nentries, double* donnees) const;
     227   void PutColToFits(int nocol, int nentries, float* donnees) const;
     228   void PutColToFits(int nocol, int nentries, int* donnees) const;
     229   void PutColToFits(int nocol, int nentries, char** donnees) const;
     230   void PutBinTabLine(long NoLine,  BnTblLine& ligne) const;
    496231
    497232
     
    501236
    502237
    503 /* Put keywords from a DVList into the primary header of the fits-file */
    504238void  DVListIntoPrimaryHeader(DVList& dvl) const;
    505239
     
    515249 };
    516250
     251 struct BnTblLine
     252 {
     253   BnTblLine() {}
     254   void setFormat(int dc, int fc, int ic, int cc, vector<string> names);
     255   bool sameFormat(const BnTblLine& btl) const;
     256
     257   void Print();
     258
     259   vector<double> ddata_;
     260   vector<float>  fdata_;
     261   vector<int>    idata_;
     262   vector<string>  cdata_;
     263   vector<string> ColName_;
     264
     265 };
    517266
    518267
  • trunk/SophyaLib/Samba/lambdaBuilder.cc

    r729 r1218  
    11#include "lambdaBuilder.h"
    22#include "nbconst.h"
     3
     4
     5/*!
     6  \class SOPHYA::Legendre
     7generate Legendre polynomials : use in two steps :
     8
     9a) instanciate Legendre(\f$x\f$, \f$lmax\f$) ; \f$x\f$ is the value for wich Legendre polynomials will be required (usually equal to \f$\cos \theta\f$) and \f$lmax\f$ is the MAXIMUM value of the order of polynomials wich will be required in the following code (all polynomials, from \f$l=0 to lmax\f$, are computed once for all by an iterative formula).
     10
     11b) get the value of Legendre polynomial for a particular value of \f$l\f$ by calling the method getPl.
     12
     13*/
    314
    415
     
    1223  array_init(lmax);
    1324}
     25   
     26/*! \fn void SOPHYA::Legendre::array_init(int_4 lmax)
     27
     28compute all \f$P_l(x,l_{max})\f$ for \f$l=1,l_{max}\f$
     29*/
    1430void Legendre::array_init(int_4 lmax)
    1531{
     
    2844TVector<r_8>* LambdaLMBuilder::normal_l_     = NULL;
    2945
     46
     47
     48/*! \class SOPHYA::LambdaLMBuilder
     49
     50
     51This class generate the coefficients :
     52\f[
     53            \lambda_l^m=\sqrt{\frac{2l+1}{4\pi}\frac{(l-m)!}{(l+m)!}}
     54            P_l^m(\cos{\theta})
     55\f]
     56where \f$P_l^m\f$ are the associated Legendre polynomials.  The above coefficients contain the theta-dependance of spheric harmonics :
     57\f[
     58            Y_{lm}(\cos{\theta})=\lambda_l^m(\cos{\theta}) e^{im\phi}.
     59\f]
     60
     61Each object has a fixed theta (radians), and maximum l and m to be calculated
     62(lmax and mmax).
     63 use the class in two steps :
     64a) instanciate  LambdaLMBuilder(\f$\theta\f$, \f$lmax\f$, \f$mmax\f$) ;  \f$lmax\f$ and \f$mmax\f$ are  MAXIMUM values for which \f$\lambda_l^m\f$ will be required in the following code (all coefficients, from \f$l=0 to lmax\f$, are computed once for all by an iterative formula).
     65b) get the values of coefficients for  particular values of \f$l\f$ and \f$m\f$ by calling the method lamlm.
     66*/
     67
     68
    3069LambdaLMBuilder::LambdaLMBuilder(r_8 theta,int_4 lmax, int_4 mmax)
    3170    {
     
    84123
    85124
     125/*! \fn void  SOPHYA::LambdaLMBuilder::updateArrayRecurrence(int_4 lmax)
     126
     127 compute a static array of coefficients independant from theta (common to all instances of the LambdaBuilder Class
     128*/
    86129void  LambdaLMBuilder::updateArrayRecurrence(int_4 lmax)
    87130   {
     
    98141   }
    99142
    100 
     143/*! \fn void  SOPHYA::LambdaLMBuilder::updateArrayLamNorm()
     144
     145 compute  static arrays of coefficients independant from theta (common to all instances of the derived  classes
     146*/
    101147void  LambdaLMBuilder::updateArrayLamNorm()
    102148     {
     
    119165
    120166
     167/*! \class SOPHYA::LambdaWXBuilder
     168
     169This class generates the coefficients :
     170\f[
     171            _{w}\lambda_l^m=-2\sqrt{\frac{2(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}} G^+_{lm}
     172\f]
     173\f[
     174            _{x}\lambda_l^m=-2\sqrt{\frac{2(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}}G^-_{lm}
     175\f]
     176where
     177\f[G^+_{lm}(\cos{\theta})=-\left( \frac{l-m^2}{\sin^2{\theta}}+\frac{1}{2}l\left(l-1\right)\right)P_l^m(\cos{\theta})+\left(l+m\right)\frac{\cos{\theta}}{\sin^2{\theta}}P^m_{l-1}(\cos{\theta})
     178\f]
     179and
     180\f[G^-_{lm}(\cos{\theta})=\frac{m}{\sin^2{\theta}}\left(\left(l-1\right)\cos{\theta}P^m_l(\cos{\theta})-\left(l+m\right)P^m_{l-1}(\cos{\theta})\right)
     181\f]
     182 \f$P_l^m\f$ are the associated Legendre polynomials.
     183
     184The coefficients express the theta-dependance of the \f$W_{lm}(\cos{\theta})\f$ and \f$X_{lm}(\cos{\theta})\f$ functions :
     185\f[W_{lm}(\cos{\theta}) = \sqrt{\frac{(l+2)!}{2(l-2)!}}_w\lambda_l^m(\cos{\theta})e^{im\phi}
     186\f]
     187\f[X_{lm}(\cos{\theta}) = -i\sqrt{\frac{(l+2)!}{2(l-2)!}}_x\lambda_l^m(\cos{\theta})e^{im\phi}
     188\f]
     189 where \f$W_{lm}(\cos{\theta})\f$ and \f$X_{lm}(\cos{\theta})\f$ are defined as :
     190
     191\f[
     192W_{lm}(\cos{\theta})=-\frac{1}{2}\sqrt{\frac{(l+2)!}{(l-2)!}}\left(
     193_{+2}Y_l^m(\cos{\theta})+_{-2}Y_l^m(\cos{\theta})\right)
     194\f]
     195\f[X_{lm}(\cos{\theta})=-\frac{i}{2}\sqrt{\frac{(l+2)!}{(l-2)!}}\left(
     196_{+2}Y_l^m(\cos{\theta})-_{-2}Y_l^m(\cos{\theta})\right)
     197\f]
     198
     199*/
    121200
    122201
     
    173252   }
    174253
     254/*!   \class SOPHYA::LambdaPMBuilder
     255
     256This class generates the coefficients
     257\f[
     258            _{\pm}\lambda_l^m=2\sqrt{\frac{(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}}\left( G^+_{lm} \mp G^-_{lm}\right)
     259\f]
     260where
     261\f[G^+_{lm}(\cos{\theta})=-\left( \frac{l-m^2}{\sin^2{\theta}}+\frac{1}{2}l\left(l-1\right)\right)P_l^m(\cos{\theta})+\left(l+m\right)\frac{\cos{\theta}}{\sin^2{\theta}}P^m_{l-1}(\cos{\theta})
     262\f]
     263and
     264\f[G^-_{lm}(\cos{\theta})=\frac{m}{\sin^2{\theta}}\left(\left(l-1\right)\cos{\theta}P^m_l(\cos{\theta})-\left(l+m\right)P^m_{l-1}(\cos{\theta})\right)
     265\f]
     266and \f$P_l^m\f$ are the associated Legendre polynomials.
     267The coefficients express the theta-dependance of the  spin-2 spherical harmonics :
     268\f[_{\pm2}Y_l^m(\cos{\theta})=_\pm\lambda_l^m(\cos{\theta})e^{im\phi}
     269\f]
     270*/
    175271
    176272LambdaPMBuilder::LambdaPMBuilder(r_8 theta, int_4 lmax, int_4 mmax) : LambdaLMBuilder(theta, lmax, mmax)
  • trunk/SophyaLib/Samba/lambdaBuilder.h

    r864 r1218  
    1010namespace SOPHYA {
    1111
    12 /*!
    13 generate Legendre polynomials : use in two steps :
    14 
    15 a) instanciate Legendre(\f$x\f$, \f$lmax\f$) ; \f$x\f$ is the value for wich Legendre polynomials will be required (usually equal to \f$\cos \theta\f$) and \f$lmax\f$ is the MAXIMUM value of the order of polynomials wich will be required in the following code (all polynomials, from \f$l=0 to lmax\f$, are computed once for all by an iterative formula).
    16 
    17 b) get the value of Legendre polynomial for a particular value of \f$l\f$ by calling the method getPl.
    18 */
     12/*!  classe pour les polynomes de legendre*/
    1913class Legendre {
    2014
     
    3226
    3327 private :
    34    /*! compute all \f$P_l(x,l_{max})\f$ for \f$l=1,l_{max}\f$ */
    3528  void array_init(int_4 lmax);
    3629
     
    4235
    4336
    44 /*!
    45 This class generate the coefficients :
    46 \f[
    47             \lambda_l^m=\sqrt{\frac{2l+1}{4\pi}\frac{(l-m)!}{(l+m)!}}
    48             P_l^m(\cos{\theta})
    49 \f]
    50 where \f$P_l^m\f$ are the associated Legendre polynomials.  The above coefficients contain the theta-dependance of spheric harmonics :
    51 \f[
    52             Y_{lm}(\cos{\theta})=\lambda_l^m(\cos{\theta}) e^{im\phi}.
    53 \f]
    54 
    55 Each object has a fixed theta (radians), and maximum l and m to be calculated
    56 (lmax and mmax).
    57  use the class in two steps :
    58 a) instanciate  LambdaLMBuilder(\f$\theta\f$, \f$lmax\f$, \f$mmax\f$) ;  \f$lmax\f$ and \f$mmax\f$ are  MAXIMUM values for which \f$\lambda_l^m\f$ will be required in the following code (all coefficients, from \f$l=0 to lmax\f$, are computed once for all by an iterative formula).
    59 b) get the values of coefficients for  particular values of \f$l\f$ and \f$m\f$ by calling the method lamlm.
    60 */
    6137 class LambdaLMBuilder {
    6238
     
    7147
    7248 private:
    73 /*! compute a static array of coefficients independant from theta (common to all instances of the LambdaBuilder Class */
    7449 void updateArrayRecurrence(int_4 lmax);
    7550
     
    8257 protected :
    8358
    84 /*! compute  static arrays of coefficients independant from theta (common to all instances of the derived  classes */
    8559 void updateArrayLamNorm();
    8660
     
    9670
    9771
    98 /*!
    99 
    100 This class generates the coefficients :
    101 \f[
    102             _{w}\lambda_l^m=-2\sqrt{\frac{2(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}} G^+_{lm}
    103 \f]
    104 \f[
    105             _{x}\lambda_l^m=-2\sqrt{\frac{2(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}}G^-_{lm}
    106 \f]
    107 where
    108 \f[G^+_{lm}(\cos{\theta})=-\left( \frac{l-m^2}{\sin^2{\theta}}+\frac{1}{2}l\left(l-1\right)\right)P_l^m(\cos{\theta})+\left(l+m\right)\frac{\cos{\theta}}{\sin^2{\theta}}P^m_{l-1}(\cos{\theta})
    109 \f]
    110 and
    111 \f[G^-_{lm}(\cos{\theta})=\frac{m}{\sin^2{\theta}}\left(\left(l-1\right)\cos{\theta}P^m_l(\cos{\theta})-\left(l+m\right)P^m_{l-1}(\cos{\theta})\right)
    112 \f]
    113  \f$P_l^m\f$ are the associated Legendre polynomials.
    114 
    115 The coefficients express the theta-dependance of the \f$W_{lm}(\cos{\theta})\f$ and \f$X_{lm}(\cos{\theta})\f$ functions :
    116 \f[W_{lm}(\cos{\theta}) = \sqrt{\frac{(l+2)!}{2(l-2)!}}_w\lambda_l^m(\cos{\theta})e^{im\phi}
    117 \f]
    118 \f[X_{lm}(\cos{\theta}) = -i\sqrt{\frac{(l+2)!}{2(l-2)!}}_x\lambda_l^m(\cos{\theta})e^{im\phi}
    119 \f]
    120  where \f$W_{lm}(\cos{\theta})\f$ and \f$X_{lm}(\cos{\theta})\f$ are defined as :
    121 
    122 \f[
    123 W_{lm}(\cos{\theta})=-\frac{1}{2}\sqrt{\frac{(l+2)!}{(l-2)!}}\left(
    124 _{+2}Y_l^m(\cos{\theta})+_{-2}Y_l^m(\cos{\theta})\right)
    125 \f]
    126 \f[X_{lm}(\cos{\theta})=-\frac{i}{2}\sqrt{\frac{(l+2)!}{(l-2)!}}\left(
    127 _{+2}Y_l^m(\cos{\theta})-_{-2}Y_l^m(\cos{\theta})\right)
    128 \f]
    129 
    130 */
    13172class LambdaWXBuilder : public LambdaLMBuilder
    13273{
     
    15697};
    15798
    158 /*!
    159 
    160 This class generates the coefficients
    161 \f[
    162             _{\pm}\lambda_l^m=2\sqrt{\frac{(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}}\left( G^+_{lm} \mp G^-_{lm}\right)
    163 \f]
    164 where
    165 \f[G^+_{lm}(\cos{\theta})=-\left( \frac{l-m^2}{\sin^2{\theta}}+\frac{1}{2}l\left(l-1\right)\right)P_l^m(\cos{\theta})+\left(l+m\right)\frac{\cos{\theta}}{\sin^2{\theta}}P^m_{l-1}(\cos{\theta})
    166 \f]
    167 and
    168 \f[G^-_{lm}(\cos{\theta})=\frac{m}{\sin^2{\theta}}\left(\left(l-1\right)\cos{\theta}P^m_l(\cos{\theta})-\left(l+m\right)P^m_{l-1}(\cos{\theta})\right)
    169 \f]
    170 and \f$P_l^m\f$ are the associated Legendre polynomials.
    171 The coefficients express the theta-dependance of the  spin-2 spherical harmonics :
    172 \f[_{\pm2}Y_l^m(\cos{\theta})=_\pm\lambda_l^m(\cos{\theta})e^{im\phi}
    173 \f]
    174 */
    17599class LambdaPMBuilder : public LambdaLMBuilder
    176100{
  • trunk/SophyaLib/Samba/sphericaltransformserver.cc

    r833 r1218  
    88#include "nbmath.h"
    99
    10 
     10/*! \class SOPHYA::SphericalTransformServer
     11
     12 Class for performing analysis and synthesis of sky maps using spin-0 or spin-2 spherical harmonics.
     13
     14Maps must be SOPHYA SphericalMaps (SphereGorski or SphereThetaPhi).
     15
     16Temperature and polarization (Stokes parameters) can be developped on spherical harmonics :
     17\f[
     18\frac{\Delta T}{T}(\hat{n})=\sum_{lm}a_{lm}^TY_l^m(\hat{n})
     19\f]
     20\f[
     21Q(\hat{n})=\frac{1}{\sqrt{2}}\sum_{lm}N_l\left(a_{lm}^EW_{lm}(\hat{n})+a_{lm}^BX_{lm}(\hat{n})\right)
     22\f]
     23\f[
     24U(\hat{n})=-\frac{1}{\sqrt{2}}\sum_{lm}N_l\left(a_{lm}^EX_{lm}(\hat{n})-a_{lm}^BW_{lm}(\hat{n})\right)
     25\f]
     26\f[
     27\left(Q \pm iU\right)(\hat{n})=\sum_{lm}a_{\pm 2lm}\, _{\pm 2}Y_l^m(\hat{n})
     28\f]
     29
     30\f[
     31Y_l^m(\hat{n})=\lambda_l^m(\theta)e^{im\phi}
     32\f]
     33\f[
     34_{\pm}Y_l^m(\hat{n})=_{\pm}\lambda_l^m(\theta)e^{im\phi}
     35\f]
     36\f[
     37W_{lm}(\hat{n})=\frac{1}{N_l}\,_{w}\lambda_l^m(\theta)e^{im\phi}
     38\f]
     39\f[
     40X_{lm}(\hat{n})=\frac{-i}{N_l}\,_{x}\lambda_l^m(\theta)e^{im\phi}
     41\f]
     42
     43(see LambdaLMBuilder, LambdaPMBuilder, LambdaWXBuilder classes)
     44
     45power spectra :
     46
     47\f[
     48C_l^T=\frac{1}{2l+1}\sum_{m=0}^{+ \infty }\left|a_{lm}^T\right|^2=\langle\left|a_{lm}^T\right|^2\rangle
     49\f]
     50\f[
     51C_l^E=\frac{1}{2l+1}\sum_{m=0}^{+\infty}\left|a_{lm}^E\right|^2=\langle\left|a_{lm}^E\right|^2\rangle
     52\f]
     53\f[
     54C_l^B=\frac{1}{2l+1}\sum_{m=0}^{+\infty}\left|a_{lm}^B\right|^2=\langle\left|a_{lm}^B\right|^2\rangle
     55\f]
     56
     57\arg
     58\b Synthesis : Get temperature and polarization maps  from \f$a_{lm}\f$ coefficients or from power spectra, (methods GenerateFrom...).
     59
     60\b Temperature:
     61\f[
     62\frac{\Delta T}{T}(\hat{n})=\sum_{lm}a_{lm}^TY_l^m(\hat{n}) = \sum_{-\infty}^{+\infty}b_m(\theta)e^{im\phi}
     63\f]
     64
     65with
     66\f[
     67b_m(\theta)=\sum_{l=\left|m\right|}^{+\infty}a_{lm}^T\lambda_l^m(\theta)
     68\f]
     69
     70\b Polarisation
     71\f[
     72Q \pm iU = \sum_{-\infty}^{+\infty}b_m^{\pm}(\theta)e^{im\phi}
     73\f]
     74
     75where :
     76\f[
     77b_m^{\pm}(\theta) = \sum_{l=\left|m\right|}^{+\infty}a_{\pm 2lm}\,_{\pm}\lambda_l^m(\theta)
     78\f]
     79
     80or :
     81\f[
     82Q  = \sum_{-\infty}^{+\infty}b_m^{Q}(\theta)e^{im\phi}
     83\f]
     84\f[
     85U  = \sum_{-\infty}^{+\infty}b_m^{U}(\theta)e^{im\phi}
     86\f]
     87
     88where:
     89\f[
     90b_m^{Q}(\theta) = \frac{1}{\sqrt{2}}\sum_{l=\left|m\right|}^{+\infty}\left(a_{lm}^E\,_{w}\lambda_l^m(\theta)-ia_{lm}^B\,_{x}\lambda_l^m(\theta)\right)
     91\f]
     92\f[
     93b_m^{U}(\theta) = \frac{1}{\sqrt{2}}\sum_{l=\left|m\right|}^{+\infty}\left(ia_{lm}^E\,_{x}\lambda_l^m(\theta)+a_{lm}^B\,_{w}\lambda_l^m(\theta)\right)
     94\f]
     95
     96Since the pixelization provides "slices" with constant \f$\theta\f$ and \f$\phi\f$ equally distributed  on \f$2\pi\f$  \f$\frac{\Delta T}{T}\f$, \f$Q\f$,\f$U\f$  can be computed by FFT.
     97
     98
     99\arg
     100\b Analysis :  Get \f$a_{lm}\f$ coefficients or  power spectra from temperature and polarization maps   (methods DecomposeTo...).
     101
     102\b Temperature:
     103\f[
     104a_{lm}^T=\int\frac{\Delta T}{T}(\hat{n})Y_l^{m*}(\hat{n})d\hat{n}
     105\f]
     106
     107approximated as :
     108\f[
     109a_{lm}^T=\sum_{\theta_k}\omega_kC_m(\theta_k)\lambda_l^m(\theta_k)
     110\f]
     111where :
     112\f[
     113C_m (\theta _k)=\sum_{\phi _{k\prime}}\frac{\Delta T}{T}(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
     114\f]
     115Since the pixelization provides "slices" with constant \f$\theta\f$ and \f$\phi\f$ equally distributed  on \f$2\pi\f$ (\f$\omega_k\f$ is the solid angle of each pixel of the slice \f$\theta_k\f$) \f$C_m\f$ can be computed by FFT.
     116
     117\b polarisation:
     118
     119\f[
     120a_{\pm 2lm}=\sum_{\theta_k}\omega_kC_m^{\pm}(\theta_k)\,_{\pm}\lambda_l^m(\theta_k)
     121\f]
     122where :
     123\f[
     124C_m^{\pm} (\theta _k)=\sum_{\phi _{k\prime}}\left(Q \pm iU\right)(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
     125\f]
     126or :
     127
     128\f[
     129a_{lm}^E=\frac{1}{\sqrt{2}}\sum_{\theta_k}\omega_k\left(C_m^{Q}(\theta_k)\,_{w}\lambda_l^m(\theta_k)-iC_m^{U}(\theta_k)\,_{x}\lambda_l^m(\theta_k)\right)
     130\f]
     131\f[
     132a_{lm}^B=\frac{1}{\sqrt{2}}\sum_{\theta_k}\omega_k\left(iC_m^{Q}(\theta_k)\,_{x}\lambda_l^m(\theta_k)+C_m^{U}(\theta_k)\,_{w}\lambda_l^m(\theta_k)\right)
     133\f]
     134
     135where :
     136\f[
     137C_m^{Q} (\theta _k)=\sum_{\phi _{k\prime}}Q(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
     138\f]
     139\f[
     140C_m^{U} (\theta _k)=\sum_{\phi _{k\prime}}U(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
     141\f]
     142
     143 */
     144
     145 /*! \fn void SOPHYA::SphericalTransformServer::GenerateFromAlm( SphericalMap<T>& map, int_4 pixelSizeIndex, const Alm<T>& alm) const
     146
     147 synthesis of a temperature  map from  Alm coefficients
     148*/
    11149template<class T>
    12150void SphericalTransformServer<T>::GenerateFromAlm( SphericalMap<T>& map, int_4 pixelSizeIndex, const Alm<T>& alm) const
     
    125263
    126264
     265  /*! \fn TVector< complex<T> >  SOPHYA::SphericalTransformServer::fourierSynthesisFromB(const Bm<complex<T> >& b_m,  int_4 nph, r_8 phi0) const
     266
     267\return a vector with nph elements  which are  sums :\f$\sum_{m=-mmax}^{mmax}b_m(\theta)e^{im\varphi}\f$ for nph values of \f$\varphi\f$ regularly distributed in \f$[0,\pi]\f$ ( calculated by FFT)
     268
     269  The object b_m (\f$b_m\f$) of the class Bm is a special vector which index goes from -mmax to mmax.
     270  */
    127271template<class T>
    128272TVector< complex<T> >  SphericalTransformServer<T>::fourierSynthesisFromB(const Bm<complex<T> >& b_m,  int_4 nph, r_8 phi0) const
     
    205349
    206350//********************************************
     351/*! \fn TVector<T>  SOPHYA::SphericalTransformServer::RfourierSynthesisFromB(const Bm<complex<T> >& b_m,  int_4 nph, r_8 phi0) const
     352
     353same as fourierSynthesisFromB, but return a real vector, taking into account the fact that b(-m) is conjugate of b(m) */
    207354template<class T>
    208355TVector<T>  SphericalTransformServer<T>::RfourierSynthesisFromB(const Bm<complex<T> >& b_m,  int_4 nph, r_8 phi0) const
     
    285432//*******************************************
    286433
     434 /*! \fn  Alm<T> SOPHYA::SphericalTransformServer::DecomposeToAlm(const SphericalMap<T>& map, int_4 nlmax, r_8 cos_theta_cut) const
     435
     436\return the Alm coefficients from analysis of a temperature map.
     437
     438    \param<nlmax> : maximum value of the l index
     439
     440     \param<cos_theta_cut> : cosinus of the symmetric cut EULER angle theta : cos_theta_cut=0 means no cut ; cos_theta_cut=1 all the sphere is cut.
     441  */
    287442template<class T>
    288443 Alm<T> SphericalTransformServer<T>::DecomposeToAlm(const SphericalMap<T>& map, int_4 nlmax, r_8 cos_theta_cut) const
     
    346501  return alm;
    347502}
     503  /*! \fn TVector< complex<T> > SOPHYA::SphericalTransformServer::CFromFourierAnalysis(int_4 nmmax, const TVector<complex<T> >datain, r_8 phi0) const
     504
     505\return a vector with mmax elements  which are  sums :
     506\f$\sum_{k=0}^{nphi}datain(\theta,\varphi_k)e^{im\varphi_k}\f$ for (mmax+1) values of \f$m\f$ from 0 to mmax.
     507   */
    348508template<class T>
    349509TVector< complex<T> > SphericalTransformServer<T>::CFromFourierAnalysis(int_4 nmmax, const TVector<complex<T> >datain, r_8 phi0) const
     
    384544
    385545//&&&&&&&&& nouvelle version
     546/* \fn TVector< complex<T> > SOPHYA::SphericalTransformServer::CFromFourierAnalysis(int_4 nmmax, const TVector<T> datain, r_8 phi0) const
     547
     548same as previous one, but with a "datain" which is real (not complex) */
    386549template<class T>
    387550TVector< complex<T> > SphericalTransformServer<T>::CFromFourierAnalysis(int_4 nmmax, const TVector<T> datain, r_8 phi0) const
     
    445608}
    446609
     610 /*! \fn void SOPHYA::SphericalTransformServer::GenerateFromAlm(SphericalMap<T>& mapq,
     611                                               SphericalMap<T>& mapu,
     612                                               int_4 pixelSizeIndex,
     613                                               const Alm<T>& alme,
     614                                               const Alm<T>& almb) const
     615
     616synthesis of a polarization map from  Alm coefficients. The spheres mapq and mapu contain respectively the Stokes parameters. */
    447617template<class T>
    448618void SphericalTransformServer<T>::GenerateFromAlm(SphericalMap<T>& mapq,
     
    521691
    522692
     693 /*! \fn void SOPHYA::SphericalTransformServer::DecomposeToAlm(const SphericalMap<T>& mapq,
     694                                              const SphericalMap<T>& mapu,
     695                                              Alm<T>& alme,
     696                                              Alm<T>& almb,
     697                                              int_4 nlmax,
     698                                              r_8 cos_theta_cut) const
     699
     700analysis of a polarization map into Alm coefficients.
     701
     702 The spheres \c mapq and \c mapu contain respectively the Stokes parameters.
     703
     704 \c a2lme and \c a2lmb will receive respectively electric and magnetic Alm's
     705    nlmax : maximum value of the l index
     706
     707 \c cos_theta_cut : cosinus of the symmetric cut EULER angle theta : cos_theta_cut=0 means no cut ; cos_theta_cut=1 all the sphere is cut.
     708 */
    523709template<class T>
    524710void SphericalTransformServer<T>::DecomposeToAlm(const SphericalMap<T>& mapq,
     
    578764
    579765
     766 /*! \fn void SOPHYA::SphericalTransformServer::almFromWX(int_4 nlmax, int_4 nmmax,
     767                                         r_8 phi0, r_8 domega,
     768                                         r_8 theta,
     769                                         const TVector<T>& dataq,
     770                                         const TVector<T>& datau,
     771                                         Alm<T>& alme,
     772                                         Alm<T>& almb) const
     773
     774Compute polarized Alm's as :
     775\f[
     776a_{lm}^E=\frac{1}{\sqrt{2}}\sum_{slices}{\omega_{pix}\left(\,_{w}\lambda_l^m\tilde{Q}-i\,_{x}\lambda_l^m\tilde{U}\right)}
     777\f]
     778\f[
     779a_{lm}^B=\frac{1}{\sqrt{2}}\sum_{slices}{\omega_{pix}\left(i\,_{x}\lambda_l^m\tilde{Q}+\,_{w}\lambda_l^m\tilde{U}\right)}
     780\f]
     781
     782where \f$\tilde{Q}\f$ and \f$\tilde{U}\f$ are C-coefficients computed by FFT (method CFromFourierAnalysis, called by present method) from the Stokes parameters.
     783
     784\f$\omega_{pix}\f$ are solid angle of each pixel.
     785
     786dataq, datau : Stokes parameters.
     787
     788  */
    580789template<class T>
    581790void SphericalTransformServer<T>::almFromWX(int_4 nlmax, int_4 nmmax,
     
    628837
    629838
    630 template<class T>
    631 void SphericalTransformServer<T>::almFromPM(int_4 nph, int_4 nlmax, int_4 nmmax,
     839 /*! \fn void SOPHYA::SphericalTransformServer::almFromPM(int_4 nph, int_4 nlmax,
     840                                         int_4 nmmax,
     841                                         r_8 phi0, r_8 domega, 
     842                                         r_8 theta,
     843                                         const TVector<T>& dataq,
     844                                         const TVector<T>& datau,
     845                                         Alm<T>& alme,
     846                                         Alm<T>& almb) const
     847
     848Compute polarized Alm's as :
     849\f[
     850a_{lm}^E=-\frac{1}{2}\sum_{slices}{\omega_{pix}\left(\,_{+}\lambda_l^m\tilde{P^+}+\,_{-}\lambda_l^m\tilde{P^-}\right)}
     851\f]
     852\f[
     853a_{lm}^B=\frac{i}{2}\sum_{slices}{\omega_{pix}\left(\,_{+}\lambda_l^m\tilde{P^+}-\,_{-}\lambda_l^m\tilde{P^-}\right)}
     854\f]
     855
     856where \f$\tilde{P^{\pm}}=\tilde{Q}\pm\tilde{U}\f$  computed by FFT (method CFromFourierAnalysis, called by present method) from the Stokes parameters,\f$Q\f$ and \f$U\f$ .
     857
     858\f$\omega_{pix}\f$ are solid angle of each pixel.
     859
     860dataq, datau : Stokes parameters.
     861
     862  */
     863template<class T>
     864void SphericalTransformServer<T>::almFromPM(int_4 nph, int_4 nlmax,
     865                                         int_4 nmmax,
    632866                                         r_8 phi0, r_8 domega, 
    633867                                         r_8 theta,
     
    673907
    674908
     909/*! \fn void SOPHYA::SphericalTransformServer::mapFromWX(int_4 nlmax, int_4 nmmax,
     910                                         SphericalMap<T>& mapq,
     911                                         SphericalMap<T>& mapu,
     912                                         const Alm<T>& alme,
     913                                         const Alm<T>& almb) const
     914
     915synthesis of Stokes parameters following formulae :
     916
     917\f[
     918Q=\sum_{m=-mmax}^{mmax}b_m^qe^{im\varphi}
     919\f]
     920\f[
     921U=\sum_{m=-mmax}^{mmax}b_m^ue^{im\varphi}
     922\f]
     923
     924computed by FFT (method fourierSynthesisFromB called by the present one)
     925
     926with :
     927
     928\f[
     929b_m^q=-\frac{1}{\sqrt{2}}\sum_{l=|m|}^{lmax}{\left(\,_{w}\lambda_l^ma_{lm}^E-i\,_{x}\lambda_l^ma_{lm}^B\right) }
     930\f]
     931\f[
     932b_m^u=\frac{1}{\sqrt{2}}\sum_{l=|m|}^{lmax}{\left(i\,_{x}\lambda_l^ma_{lm}^E+\,_{w}\lambda_l^ma_{lm}^B\right) }
     933\f]
     934 */
    675935template<class T>
    676936void SphericalTransformServer<T>::mapFromWX(int_4 nlmax, int_4 nmmax,
     
    7441004    }
    7451005}
     1006/*! \fn void SOPHYA::SphericalTransformServer::mapFromPM(int_4 nlmax, int_4 nmmax,
     1007                                         SphericalMap<T>& mapq,
     1008                                         SphericalMap<T>& mapu,
     1009                                         const Alm<T>& alme,
     1010                                         const Alm<T>& almb) const
     1011
     1012synthesis of polarizations following formulae :
     1013
     1014\f[
     1015P^+ = \sum_{m=-mmax}^{mmax} {b_m^+e^{im\varphi} }
     1016\f]
     1017\f[
     1018P^- = \sum_{m=-mmax}^{mmax} {b_m^-e^{im\varphi} }
     1019\f]
     1020
     1021computed by FFT (method fourierSynthesisFromB called by the present one)
     1022
     1023with :
     1024
     1025\f[
     1026b_m^+=-\sum_{l=|m|}^{lmax}{\,_{+}\lambda_l^m \left( a_{lm}^E+ia_{lm}^B \right) }
     1027\f]
     1028\f[
     1029b_m^-=-\sum_{l=|m|}^{lmax}{\,_{+}\lambda_l^m \left( a_{lm}^E-ia_{lm}^B \right) }
     1030\f]
     1031 */
    7461032template<class T>
    7471033void SphericalTransformServer<T>::mapFromPM(int_4 nlmax, int_4 nmmax,
     
    8091095
    8101096
     1097  /*! \fn void SOPHYA::SphericalTransformServer::GenerateFromCl(SphericalMap<T>& sphq,
     1098                                              SphericalMap<T>& sphu,
     1099                                              int_4 pixelSizeIndex,
     1100                                              const TVector<T>& Cle,
     1101                                              const TVector<T>& Clb,
     1102                                              const r_8 fwhm) const
     1103
     1104synthesis of a polarization  map from  power spectra electric-Cl and magnetic-Cl (Alm's are generated randomly, following a gaussian distribution).
     1105  \param fwhm FWHM in arcmin for random generation of Alm's (eg. 5)
     1106*/
    8111107template<class T>
    8121108void SphericalTransformServer<T>::GenerateFromCl(SphericalMap<T>& sphq,
     
    8331129  GenerateFromAlm(sphq,sphu,pixelSizeIndex,a2lme,a2lmb);
    8341130}
     1131 /*! \fn void SOPHYA::SphericalTransformServer::GenerateFromCl(SphericalMap<T>& sph,
     1132                                                 int_4 pixelSizeIndex,
     1133                                              const TVector<T>& Cl,
     1134                                                 const r_8 fwhm)  const
     1135
     1136synthesis of a temperature  map from  power spectrum Cl (Alm's are generated randomly, following a gaussian distribution). */
    8351137template<class T>
    8361138void SphericalTransformServer<T>::GenerateFromCl(SphericalMap<T>& sph,
     
    8461148
    8471149
     1150/*! \fn TVector<T>  SOPHYA::SphericalTransformServer::DecomposeToCl(const SphericalMap<T>& sph, int_4 nlmax, r_8 cos_theta_cut) const
     1151
     1152\return power spectrum from analysis of a temperature map.
     1153
     1154     \param<nlmax> : maximum value of the l index
     1155
     1156     \param<cos_theta_cut> : cosinus of the symmetric cut EULER angle theta : cos_theta_cut=0 means no cut ; cos_theta_cut=1 all the sphere is cut.
     1157  */
    8481158template <class T>
    8491159TVector<T>  SphericalTransformServer<T>::DecomposeToCl(const SphericalMap<T>& sph, int_4 nlmax, r_8 cos_theta_cut) const
  • trunk/SophyaLib/Samba/sphericaltransformserver.h

    r866 r1218  
    1111namespace SOPHYA {
    1212
    13 //
    14 /*! Class for performing analysis and synthesis of sky maps using spin-0 or spin-2 spherical harmonics.
    1513
    16 Maps must be SOPHYA SphericalMaps (SphereGorski or SphereThetaPhi).
    17 
    18 Temperature and polarization (Stokes parameters) can be developped on spherical harmonics :
    19 \f[
    20 \frac{\Delta T}{T}(\hat{n})=\sum_{lm}a_{lm}^TY_l^m(\hat{n})
    21 \f]
    22 \f[
    23 Q(\hat{n})=\frac{1}{\sqrt{2}}\sum_{lm}N_l\left(a_{lm}^EW_{lm}(\hat{n})+a_{lm}^BX_{lm}(\hat{n})\right)
    24 \f]
    25 \f[
    26 U(\hat{n})=-\frac{1}{\sqrt{2}}\sum_{lm}N_l\left(a_{lm}^EX_{lm}(\hat{n})-a_{lm}^BW_{lm}(\hat{n})\right)
    27 \f]
    28 \f[
    29 \left(Q \pm iU\right)(\hat{n})=\sum_{lm}a_{\pm 2lm}\, _{\pm 2}Y_l^m(\hat{n})
    30 \f]
    31 
    32 \f[
    33 Y_l^m(\hat{n})=\lambda_l^m(\theta)e^{im\phi}
    34 \f]
    35 \f[
    36 _{\pm}Y_l^m(\hat{n})=_{\pm}\lambda_l^m(\theta)e^{im\phi}
    37 \f]
    38 \f[
    39 W_{lm}(\hat{n})=\frac{1}{N_l}\,_{w}\lambda_l^m(\theta)e^{im\phi}
    40 \f]
    41 \f[
    42 X_{lm}(\hat{n})=\frac{-i}{N_l}\,_{x}\lambda_l^m(\theta)e^{im\phi}
    43 \f]
    44 
    45 (see LambdaLMBuilder, LambdaPMBuilder, LambdaWXBuilder classes)
    46 
    47 power spectra :
    48 
    49 \f[
    50 C_l^T=\frac{1}{2l+1}\sum_{m=0}^{+ \infty }\left|a_{lm}^T\right|^2=\langle\left|a_{lm}^T\right|^2\rangle
    51 \f]
    52 \f[
    53 C_l^E=\frac{1}{2l+1}\sum_{m=0}^{+\infty}\left|a_{lm}^E\right|^2=\langle\left|a_{lm}^E\right|^2\rangle
    54 \f]
    55 \f[
    56 C_l^B=\frac{1}{2l+1}\sum_{m=0}^{+\infty}\left|a_{lm}^B\right|^2=\langle\left|a_{lm}^B\right|^2\rangle
    57 \f]
    58 
    59 \arg
    60 \b Synthesis : Get temperature and polarization maps  from \f$a_{lm}\f$ coefficients or from power spectra, (methods GenerateFrom...).
    61 
    62 \b Temperature:
    63 \f[
    64 \frac{\Delta T}{T}(\hat{n})=\sum_{lm}a_{lm}^TY_l^m(\hat{n}) = \sum_{-\infty}^{+\infty}b_m(\theta)e^{im\phi}
    65 \f]
    66 
    67 with
    68 \f[
    69 b_m(\theta)=\sum_{l=\left|m\right|}^{+\infty}a_{lm}^T\lambda_l^m(\theta)
    70 \f]
    71 
    72 \b Polarisation
    73 \f[
    74 Q \pm iU = \sum_{-\infty}^{+\infty}b_m^{\pm}(\theta)e^{im\phi}
    75 \f]
    76 
    77 where :
    78 \f[
    79 b_m^{\pm}(\theta) = \sum_{l=\left|m\right|}^{+\infty}a_{\pm 2lm}\,_{\pm}\lambda_l^m(\theta)
    80 \f]
    81 
    82 or :
    83 \f[
    84 Q  = \sum_{-\infty}^{+\infty}b_m^{Q}(\theta)e^{im\phi}
    85 \f]
    86 \f[
    87 U  = \sum_{-\infty}^{+\infty}b_m^{U}(\theta)e^{im\phi}
    88 \f]
    89 
    90 where:
    91 \f[
    92 b_m^{Q}(\theta) = \frac{1}{\sqrt{2}}\sum_{l=\left|m\right|}^{+\infty}\left(a_{lm}^E\,_{w}\lambda_l^m(\theta)-ia_{lm}^B\,_{x}\lambda_l^m(\theta)\right)
    93 \f]
    94 \f[
    95 b_m^{U}(\theta) = \frac{1}{\sqrt{2}}\sum_{l=\left|m\right|}^{+\infty}\left(ia_{lm}^E\,_{x}\lambda_l^m(\theta)+a_{lm}^B\,_{w}\lambda_l^m(\theta)\right)
    96 \f]
    97 
    98 Since the pixelization provides "slices" with constant \f$\theta\f$ and \f$\phi\f$ equally distributed  on \f$2\pi\f$  \f$\frac{\Delta T}{T}\f$, \f$Q\f$,\f$U\f$  can be computed by FFT.
    99 
    100 
    101 \arg
    102 \b Analysis :  Get \f$a_{lm}\f$ coefficients or  power spectra from temperature and polarization maps   (methods DecomposeTo...).
    103 
    104 \b Temperature:
    105 \f[
    106 a_{lm}^T=\int\frac{\Delta T}{T}(\hat{n})Y_l^{m*}(\hat{n})d\hat{n}
    107 \f]
    108 
    109 approximated as :
    110 \f[
    111 a_{lm}^T=\sum_{\theta_k}\omega_kC_m(\theta_k)\lambda_l^m(\theta_k)
    112 \f]
    113 where :
    114 \f[
    115 C_m (\theta _k)=\sum_{\phi _{k\prime}}\frac{\Delta T}{T}(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
    116 \f]
    117 Since the pixelization provides "slices" with constant \f$\theta\f$ and \f$\phi\f$ equally distributed  on \f$2\pi\f$ (\f$\omega_k\f$ is the solid angle of each pixel of the slice \f$\theta_k\f$) \f$C_m\f$ can be computed by FFT.
    118 
    119 \b polarisation:
    120 
    121 \f[
    122 a_{\pm 2lm}=\sum_{\theta_k}\omega_kC_m^{\pm}(\theta_k)\,_{\pm}\lambda_l^m(\theta_k)
    123 \f]
    124 where :
    125 \f[
    126 C_m^{\pm} (\theta _k)=\sum_{\phi _{k\prime}}\left(Q \pm iU\right)(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
    127 \f]
    128 or :
    129 
    130 \f[
    131 a_{lm}^E=\frac{1}{\sqrt{2}}\sum_{\theta_k}\omega_k\left(C_m^{Q}(\theta_k)\,_{w}\lambda_l^m(\theta_k)-iC_m^{U}(\theta_k)\,_{x}\lambda_l^m(\theta_k)\right)
    132 \f]
    133 \f[
    134 a_{lm}^B=\frac{1}{\sqrt{2}}\sum_{\theta_k}\omega_k\left(iC_m^{Q}(\theta_k)\,_{x}\lambda_l^m(\theta_k)+C_m^{U}(\theta_k)\,_{w}\lambda_l^m(\theta_k)\right)
    135 \f]
    136 
    137 where :
    138 \f[
    139 C_m^{Q} (\theta _k)=\sum_{\phi _{k\prime}}Q(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
    140 \f]
    141 \f[
    142 C_m^{U} (\theta _k)=\sum_{\phi _{k\prime}}U(\theta _k,\phi_{k\prime})e^{-im\phi _{k\prime}}
    143 \f]
    144 
    145  */
    14614template <class T>
    14715class SphericalTransformServer
     
    15018 public:
    15119
    152  SphericalTransformServer() 
    153 {
    154     fftIntfPtr_=new FFTPackServer;
    155     fftIntfPtr_->setNormalize(false);
    156 };
    157  ~SphericalTransformServer(){ if (fftIntfPtr_!=NULL) delete fftIntfPtr_;};
    158 /*!
     20  SphericalTransformServer() 
     21    {
     22      fftIntfPtr_=new FFTPackServer;
     23      fftIntfPtr_->setNormalize(false);
     24    };
     25  ~SphericalTransformServer(){ if (fftIntfPtr_!=NULL) delete fftIntfPtr_;};
     26 
     27 /*!
    15928 Set a fft server. The constructor sets a default fft server (fft-pack). So it is not necessary to call this method for a standard use.
    16029*/
    161  void SetFFTServer(FFTServerInterface* srv=NULL)
    162 {
    163   if (fftIntfPtr_!=NULL) delete fftIntfPtr_;
    164   fftIntfPtr_=srv;
    165   fftIntfPtr_->setNormalize(false);
    166 }
    167  /*! synthesis of a temperature  map from  Alm coefficients */
    168  void GenerateFromAlm( SphericalMap<T>& map, int_4 pixelSizeIndex, const Alm<T>& alm) const;
    169  /*! synthesis of a polarization map from  Alm coefficients. The spheres mapq and mapu contain respectively the Stokes parameters. */
    170  void GenerateFromAlm(SphericalMap<T>& mapq, SphericalMap<T>& mapu, int_4 pixelSizeIndex, const Alm<T>& alme,    const Alm<T>& almb) const;
     30  void SetFFTServer(FFTServerInterface* srv=NULL)
     31    {
     32      if (fftIntfPtr_!=NULL) delete fftIntfPtr_;
     33      fftIntfPtr_=srv;
     34      fftIntfPtr_->setNormalize(false);
     35    }
    17136
    172  /*! synthesis of a temperature  map from  power spectrum Cl (Alm's are generated randomly, following a gaussian distribution). */
    173  void GenerateFromCl(SphericalMap<T>& sph, int_4 pixelSizeIndex,
     37  void GenerateFromAlm( SphericalMap<T>& map, int_4 pixelSizeIndex, const Alm<T>& alm) const;
     38  void GenerateFromAlm(SphericalMap<T>& mapq, SphericalMap<T>& mapu, int_4 pixelSizeIndex, const Alm<T>& alme,    const Alm<T>& almb) const;
     39
     40  void GenerateFromCl(SphericalMap<T>& sph, int_4 pixelSizeIndex,
    17441                     const TVector<T>& Cl, const r_8 fwhm) const;
    175   /*! synthesis of a polarization  map from  power spectra electric-Cl and magnetic-Cl (Alm's are generated randomly, following a gaussian distribution).
    176   \param fwhm FWHM in arcmin for random generation of Alm's (eg. 5)
    177 
    178 */
    179  void GenerateFromCl(SphericalMap<T>& sphq, SphericalMap<T>& sphu,
     42  void GenerateFromCl(SphericalMap<T>& sphq, SphericalMap<T>& sphu,
    18043                     int_4 pixelSizeIndex,
    18144                     const TVector<T>& Cle, const TVector<T>& Clb,
    18245                    const r_8 fwhm) const;
    183  /*!return the Alm coefficients from analysis of a temperature map.
    184 
    185     \param<nlmax> : maximum value of the l index
    186 
    187      \param<cos_theta_cut> : cosinus of the symmetric cut EULER angle theta : cos_theta_cut=0 means no cut ; cos_theta_cut=1 all the sphere is cut.
    188   */
    18946
    19047
    191 Alm<T> DecomposeToAlm(const SphericalMap<T>& map, int_4 nlmax, r_8 cos_theta_cut) const;
    192  /*!analysis of a polarization map into Alm coefficients.
     48  Alm<T> DecomposeToAlm(const SphericalMap<T>& map, int_4 nlmax, r_8 cos_theta_cut) const;
    19349
    194  The spheres \c mapq and \c mapu contain respectively the Stokes parameters.
    195 
    196  \c a2lme and \c a2lmb will receive respectively electric and magnetic Alm's
    197     nlmax : maximum value of the l index
    198 
    199  \c cos_theta_cut : cosinus of the symmetric cut EULER angle theta : cos_theta_cut=0 means no cut ; cos_theta_cut=1 all the sphere is cut.
    200  */
    201 
    202  void DecomposeToAlm(const SphericalMap<T>& mapq, const SphericalMap<T>& mapu,
     50  void DecomposeToAlm(const SphericalMap<T>& mapq, const SphericalMap<T>& mapu,
    20351                     Alm<T>& a2lme, Alm<T>& a2lmb,
    20452                     int_4 nlmax, r_8 cos_theta_cut) const;
    20553
    206 /*!return power spectrum from analysis of a temperature map.
    207 
    208      \param<nlmax> : maximum value of the l index
    209 
    210      \param<cos_theta_cut> : cosinus of the symmetric cut EULER angle theta : cos_theta_cut=0 means no cut ; cos_theta_cut=1 all the sphere is cut.
    211   */
    212  TVector<T>  DecomposeToCl(const SphericalMap<T>& sph, 
     54  TVector<T>  DecomposeToCl(const SphericalMap<T>& sph, 
    21355                           int_4 nlmax, r_8 cos_theta_cut) const;
    21456
    21557
    216   private:
    217   /*! return a vector with nph elements  which are  sums :\f$\sum_{m=-mmax}^{mmax}b_m(\theta)e^{im\varphi}\f$ for nph values of \f$\varphi\f$ regularly distributed in \f$[0,\pi]\f$ ( calculated by FFT)
    218 
    219   The object b_m (\f$b_m\f$) of the class Bm is a special vector which index goes from -mmax to mmax.
    220   */
    221  TVector< complex<T> >  fourierSynthesisFromB(const Bm<complex<T> >& b_m, 
     58 private:
     59  TVector< complex<T> >  fourierSynthesisFromB(const Bm<complex<T> >& b_m, 
    22260                                  int_4 nph, r_8 phi0) const;
    223 /*! same as fourierSynthesisFromB, but return a real vector, taking into account the fact that b(-m) is conjugate of b(m) */
    224  TVector<T>  RfourierSynthesisFromB(const Bm<complex<T> >& b_m, 
     61  TVector<T>  RfourierSynthesisFromB(const Bm<complex<T> >& b_m, 
    22562                                  int_4 nph, r_8 phi0) const;
    22663
    227   /*! return a vector with mmax elements  which are  sums :
    228 \f$\sum_{k=0}^{nphi}datain(\theta,\varphi_k)e^{im\varphi_k}\f$ for (mmax+1) values of \f$m\f$ from 0 to mmax.
    229    */
    230  TVector< complex<T> > CFromFourierAnalysis(int_4 mmax,
     64  TVector< complex<T> > CFromFourierAnalysis(int_4 mmax,
    23165                                  const TVector<complex<T> > datain,
    23266                                  r_8 phi0) const;
    233 /* same as previous one, but with a "datain" which is real (not complex) */
    234  TVector< complex<T> > CFromFourierAnalysis(int_4 mmax,
     67  TVector< complex<T> > CFromFourierAnalysis(int_4 mmax,
    23568                         const TVector<T> datain, 
    23669                                  r_8 phi0) const;
    23770
    238  /*!
    239 Compute polarized Alm's as :
    240 \f[
    241 a_{lm}^E=\frac{1}{\sqrt{2}}\sum_{slices}{\omega_{pix}\left(\,_{w}\lambda_l^m\tilde{Q}-i\,_{x}\lambda_l^m\tilde{U}\right)}
    242 \f]
    243 \f[
    244 a_{lm}^B=\frac{1}{\sqrt{2}}\sum_{slices}{\omega_{pix}\left(i\,_{x}\lambda_l^m\tilde{Q}+\,_{w}\lambda_l^m\tilde{U}\right)}
    245 \f]
    246 
    247 where \f$\tilde{Q}\f$ and \f$\tilde{U}\f$ are C-coefficients computed by FFT (method CFromFourierAnalysis, called by present method) from the Stokes parameters.
    248 
    249 \f$\omega_{pix}\f$ are solid angle of each pixel.
    250 
    251 dataq, datau : Stokes parameters.
    252 
    253   */
    254 void almFromWX(int_4 nlmax, int_4 nmmax, r_8 phi0,
     71  void almFromWX(int_4 nlmax, int_4 nmmax, r_8 phi0,
    25572               r_8 domega, r_8 theta,
    25673               const TVector<T>& dataq, const TVector<T>& datau,
    25774               Alm<T>& alme, Alm<T>& almb) const;
    258  /*!
    259 Compute polarized Alm's as :
    260 \f[
    261 a_{lm}^E=-\frac{1}{2}\sum_{slices}{\omega_{pix}\left(\,_{+}\lambda_l^m\tilde{P^+}+\,_{-}\lambda_l^m\tilde{P^-}\right)}
    262 \f]
    263 \f[
    264 a_{lm}^B=\frac{i}{2}\sum_{slices}{\omega_{pix}\left(\,_{+}\lambda_l^m\tilde{P^+}-\,_{-}\lambda_l^m\tilde{P^-}\right)}
    265 \f]
    266 
    267 where \f$\tilde{P^{\pm}}=\tilde{Q}\pm\tilde{U}\f$  computed by FFT (method CFromFourierAnalysis, called by present method) from the Stokes parameters,\f$Q\f$ and \f$U\f$ .
    268 
    269 \f$\omega_{pix}\f$ are solid angle of each pixel.
    270 
    271 dataq, datau : Stokes parameters.
    272 
    273   */
    274 void almFromPM(int_4 nph, int_4 nlmax, int_4 nmmax,
     75  void almFromPM(int_4 nph, int_4 nlmax, int_4 nmmax,
    27576               r_8 phi0, r_8 domega, r_8 theta,
    27677               const TVector<T>& dataq, const TVector<T>& datau,
    27778               Alm<T>& alme, Alm<T>& almb) const;
    27879
    279 /*! synthesis of Stokes parameters following formulae :
    280 
    281 \f[
    282 Q=\sum_{m=-mmax}^{mmax}b_m^qe^{im\varphi}
    283 \f]
    284 \f[
    285 U=\sum_{m=-mmax}^{mmax}b_m^ue^{im\varphi}
    286 \f]
    287 
    288 computed by FFT (method fourierSynthesisFromB called by the present one)
    289 
    290 with :
    291 
    292 \f[
    293 b_m^q=-\frac{1}{\sqrt{2}}\sum_{l=|m|}^{lmax}{\left(\,_{w}\lambda_l^ma_{lm}^E-i\,_{x}\lambda_l^ma_{lm}^B\right) }
    294 \f]
    295 \f[
    296 b_m^u=\frac{1}{\sqrt{2}}\sum_{l=|m|}^{lmax}{\left(i\,_{x}\lambda_l^ma_{lm}^E+\,_{w}\lambda_l^ma_{lm}^B\right) }
    297 \f]
    298  */
    299 void mapFromWX(int_4 nlmax, int_4 nmmax,
     80  void mapFromWX(int_4 nlmax, int_4 nmmax,
    30081               SphericalMap<T>& mapq, SphericalMap<T>& mapu,
    30182               const Alm<T>& alme, const Alm<T>& almb) const;
    30283
    303 /*! synthesis of polarizations following formulae :
    30484
    305 \f[
    306 P^+ = \sum_{m=-mmax}^{mmax} {b_m^+e^{im\varphi} }
    307 \f]
    308 \f[
    309 P^- = \sum_{m=-mmax}^{mmax} {b_m^-e^{im\varphi} }
    310 \f]
    311 
    312 computed by FFT (method fourierSynthesisFromB called by the present one)
    313 
    314 with :
    315 
    316 \f[
    317 b_m^+=-\sum_{l=|m|}^{lmax}{\,_{+}\lambda_l^m \left( a_{lm}^E+ia_{lm}^B \right) }
    318 \f]
    319 \f[
    320 b_m^-=-\sum_{l=|m|}^{lmax}{\,_{+}\lambda_l^m \left( a_{lm}^E-ia_{lm}^B \right) }
    321 \f]
    322  */
    323 
    324 void mapFromPM(int_4 nlmax, int_4 nmmax,
     85  void mapFromPM(int_4 nlmax, int_4 nmmax,
    32586               SphericalMap<T>& mapq, SphericalMap<T>& mapu,
    32687               const Alm<T>& alme, const Alm<T>& almb) const;
Note: See TracChangeset for help on using the changeset viewer.