source: Sophya/trunk/ArchTOIPipe/TestPipes/tsttoi2map.cc@ 2054

//#define TOISEQBUFFERED

#include <unistd.h>
#include <stdlib.h>
#include <stdio.h>
#include "toi.h"
#include "toiprocessor.h"
#include "fitstoirdr.h"
#include "fitstoiwtr.h"
#include "toimanager.h"
#ifdef TOISEQBUFFERED
#include "toiseqbuff.h"
#else
#include "toisegment.h"
#endif

#include "sambainit.h"
#include "toi2map.h"
#include "fitsspherehealpix.h"
#include "timing.h"

#include <stdexcept>

void usage(void);
void usage(void) {
 cout<<"tsttoi2map [-h] [-p lp] [-s samplemin,samplemax] [-w data_window_size]"<<endl
     <<"           [-a label_coord1] [-d label_coord2] [-b label_bolomuv]"<<endl
     <<"           [-n nlat] [-i c,h] [-o c,h]"<<endl
     <<"           [-m vmin] [-M vmax] [-f flag]"<<endl
     <<"           fitsin_point fitsin_bolo fitsphout [fitsphwout]"<<endl
     <<" -p lp : print level (def=0)"<<endl
     <<" -s samplemin,samplemax : sample range to be treated (def=all)"<<endl
     <<" -w data_window_size : window size for pipe (def=8192)"<<endl
     <<" -a label_coord1 : label fits for alpha/gLong (def=coord1)"<<endl
     <<" -d label_coord2 : label fits for delta/gLat (def=coord2)"<<endl
     <<"          coord1 = alpha or gLong ; coord2 = delta or gLat"<<endl
     <<" -b label_bolomuv : label fits for bolo value (def=boloMuV)"<<endl
     <<" -n nlat : nlat for Healpix sphere (def=128)"<<endl
     <<" -i cin : coordIn caracteristics (def=\"gdcdl\")"<<endl
     <<" -o cmap : idem -i for coordMap (def=\"g\")"<<endl
     <<" -e equi : equinoxe en annee (def=2000.)"<<endl
     <<" -m vmin : samples are good if sample value >= vmin"<<endl
     <<" -M vmax : samples are good if sample value <= vmax"<<endl
     <<" -f flag : samples are bad if match flag"<<endl
     <<" -N nbfiles : Number of fitsin_point files fitsin_point%d.fits 0..nb-1"<<endl
     <<" -I : sampleNum are implicit in fits files (def=no)"<<endl
     <<" fitsin_point : fits file for pointing"<<endl
     <<" fitsin_bolo : fits file for bolo values"<<endl
     <<" fitsphout : fits file for output Healpix sphere"<<endl
     <<" fitsphwout : fits file for output Healpix nFilled sphere (def=no)"<<endl;
}

////////////////////////////////////////////////////////////////
int main(int narg, char** arg) {

TOIManager* mgr = TOIManager::getManager();

//-- Decodage arguments
int lp = 0, width = 8192;
int nlat = 128;   // npixel = 12 * nlat^2
bool tflg=false, tmin=false, tmax=false;
r_8 vmin=-1.e30, vmax=1.e30; uint_8 badflg=0;
char *label_coord1 = "coord1", *label_coord2 = "coord2"
   , *label_bolomuv = "boloMuV";
double equi=2000.;
char *tcoorin="gdcdl", *tcoormap="g";
string fitsphwout = "";
long sdeb,sfin;
bool snimplicit = false;
bool fgprstat = true;
int nbpointfiles = 0;
int c;
while((c = getopt(narg,arg,"hIp:s:w:a:d:b:n:i:o:m:M:f:e:N:")) != -1) {
  switch (c) {
  case 's' :
    sscanf(optarg,"%ld,%ld",&sdeb,&sfin);
    cout<<"Requested Samples from "<<sdeb<<" , "<<sfin<<endl;
    if(sfin>=sdeb) mgr->setRequestedSample(sdeb,sfin);
    //else {cout<<"Bad sample interval "<<endl; exit(2);}
    break;
  case 'w' :
    sscanf(optarg,"%d",&width);
    if(width<=0) width=8192;
    cout<<"Data window size "<<width<<endl;
    break;
  case 'p' :
    sscanf(optarg,"%d",&lp);
    if(lp<0) lp=0;
    break;
  case 'a' :
    label_coord1 = optarg;
    break;
  case 'd' :
    label_coord2 = optarg;
    break;
  case 'b' :
    label_bolomuv = optarg;
    break;
  case 'n' :
    sscanf(optarg,"%d",&nlat);
    if(nlat<0) nlat=128;
    break;
  case 'i' :
    tcoorin=optarg;
    break;
  case 'o' :
    tcoormap=optarg;
    break;
  case 'e' :
    sscanf(optarg,"%lf",&equi);
    break;
  case 'm' :
    sscanf(optarg,"%lf",&vmin);
    tmin = true;
    break;
  case 'M' :
    sscanf(optarg,"%lf",&vmax);
    tmax = true;
    break;
  case 'f' :
    sscanf(optarg,"%ul",&badflg);
    tflg = true;
    break;
  case 'N' :
    nbpointfiles = atoi(optarg);
    break;
  case 'I' :
    snimplicit = true;
    break;
  case 'h' :
  default:
    usage(); exit(1);
    break;
  }
}
if(optind+2>=narg) {usage(); exit(3);}
 
char * fitsin_point = arg[optind];
char * fitsin_bolo = arg[optind+1];
string const fitsphout = arg[optind+2];
if(optind+3<narg) fitsphwout = arg[optind+3];

cout<<">>>> tsttoi2map:"<<endl
    <<"Pipe Window Size "<<width<<endl
    <<"Fits Infile Bolo "<<fitsin_bolo<<endl
    <<"  ...label_bolomuv "<<label_bolomuv<<endl;
cout<<"Fits Infile Pointing "<<fitsin_point<<" NbFiles=" <<nbpointfiles<<endl
    <<"  ...label_coord1 "<<label_coord1<<endl
    <<"  ...label_coord2 "<<label_coord2<<endl
    <<"  ...... ctype="<<tcoorin<<endl;
cout<<"Fits Healpix Sphere "<<fitsphout<<endl
    <<"  ...nlat "<<nlat<<endl;
cout<<"Fits Healpix Weight Sphere "<<fitsphwout<<endl
    <<"  ...... ctype="<<tcoormap<<endl;
cout<<"Equinoxe "<<equi<<" years"<<endl;

SophyaInit();
InitTim();

//--------------------------------------------------------------------
try {
//--------------------------------------------------------------------


 // FITS reader
 FITSTOIReader rfitsb(fitsin_bolo);
 if(snimplicit) rfitsb.setImplicitSN();
 int ncolb = rfitsb.getNOut();
 cout<<"Number of columns in fits Infile_bolo : "<<ncolb<<endl;
 if(ncolb<1) exit(-4);

 string pointfileI = fitsin_point;
 MuTyV numf=0;
 if(nbpointfiles > 0) 
   pointfileI += (string)numf + ".fits";
 cout << " Creating FITSTOIReader for pointing file " << pointfileI << endl;
 FITSTOIReader rfitsp(pointfileI);
 if(nbpointfiles > 0) 
   for(int kff=1; kff<nbpointfiles; kff++) {
     numf=kff;
     pointfileI = fitsin_point;
     pointfileI += (string)numf + ".fits";
     cout << " pointing FITSTOIReader.addFile( " << pointfileI << " )" << endl;
     rfitsp.addFile( pointfileI );
   }
 
 if(snimplicit) rfitsp.setImplicitSN();
 int ncolp = rfitsp.getNOut();
 cout<<"Number of columns in fits Infile_point : "<<ncolp<<endl;
 if(ncolp<2) exit(-5);

 // Creation de la sphere Healpix
 SphereHEALPix<r_8>* sph = new SphereHEALPix<r_8>(nlat);
 cout<<"SphereHEALPix: Type de map : "<<sph->TypeOfMap()<<endl
     <<"               Nombre de pixels : "<<sph->NbPixels()<<endl
     <<"               Nlat : "<<sph->SizeIndex()<<endl;

 // Creation de la sphere de poids Healpix
 SphereHEALPix<r_8>* wsph = NULL;
 if(fitsphwout.size()>0) {
   wsph = new SphereHEALPix<r_8>;
   cout<<"SphereHEALPix Weight Created"<<endl;
 }

 // TOI Processor
 TOI2Map toi2m(sph,wsph);
 cout<<"TOI2Map created"<<endl;
 toi2m.SetEquinox(equi);
 toi2m.SetCoorIn(tcoorin);
 toi2m.SetCoorMap(tcoormap);
 toi2m.SetTestFlag(tflg,badflg);
 toi2m.SetTestMin(tmin,vmin);
 toi2m.SetTestMax(tmax,vmax);
 toi2m.Print(cout);

 // Definition des tuyaux
#ifdef TOISEQBUFFERED
 cout<<">>>> Using TOISeqBuffered"<<endl;
 TOISeqBuffered * toicoord1in = new TOISeqBuffered("toi_coord1_in",width);
#else
 cout<<">>>> Using TOISegmented"<<endl;
 TOISegmented * toicoord1in = new TOISegmented("toi_coord1_in",width);
#endif
 // toicoord1in->setDebugLevel(1);
 rfitsp.addOutput(label_coord1,toicoord1in);
 toi2m.addInput("Coord1In",toicoord1in);

#ifdef TOISEQBUFFERED
 TOISeqBuffered * toicoord2in = new TOISeqBuffered("toi_coord2_in",width);
#else
 TOISegmented * toicoord2in = new TOISegmented("toi_coord2_in",width);
#endif
 // toicoord2in->setDebugLevel(1);
 rfitsp.addOutput(label_coord2,toicoord2in);
 toi2m.addInput("Coord2In",toicoord2in);
 
#ifdef TOISEQBUFFERED
 TOISeqBuffered * toibolin = new TOISeqBuffered("toi_bolo_in",width);
#else
 TOISegmented * toibolin = new TOISegmented("toi_bolo_in",width);
#endif
 // toibolin->setDebugLevel(1);
 rfitsb.addOutput(label_bolomuv,toibolin);
 toi2m.addInput("BoloIn",toibolin);

 // Run
 cout<<"----- FITSReaderTOI::PrintStatus() : -----"<<endl;
 rfitsp.PrintStatus(cout);
 rfitsb.PrintStatus(cout);

 PrtTim("starting threads");
 rfitsp.start();
 rfitsb.start();
 toi2m.start();

 // if(lp>1)
 //   for(int jjjj=0;jjjj<5;jjjj++) {
 //     cout<<*toicoord1in; 
 //     cout<<*toibolin; 
 //     sleep(2);
 //   }

 // Affichage de l'avancement des TOIProcessors 
 ProcSampleCounter<TOI2Map>  stats(toi2m);
 stats.InfoMessage() = "tsttoi2map/Info"; 
 stats.PrintStats(); 

 // Fin des traitements des TOIProcessors  
 mgr->joinAll();
 PrtTim("End threads");


 // Ecriture de la sphere Healpix sur fits
 {
 FitsOutFile sfits(fitsphout,FitsFile::clear);
 cout<<"tsttoi2map: Creating sphere fits file "<<fitsphout<<endl;
 sfits << *sph;
 }

 // Ecriture de la sphere Healpix sur fits
 if(wsph) {
   FitsOutFile swfits(fitsphwout,FitsFile::clear);
   cout<<"tsttoi2map: Creating sphere weight fits file "<<fitsphwout<<endl;
   swfits << *wsph;
 }

 // Nettoyage
 cout << "tsttoi2map: cleanup " << endl;
 delete sph;
 if(wsph) delete wsph;
 cout << "tsttoi2map: ----------- End of job -------------- " << endl;

//--------------------------------------------------------------------
} catch (PThrowable & exc) {
 cout<<"\ntsttoi2map: Catched Exception \n"<<(string)typeid(exc).name() 
     <<" - Msg= "<<exc.Msg()<<endl;
} catch (const std::exception & sex) {
 cout<<"\ntsttoi2map: Catched std::exception \n" 
     <<(string)typeid(sex).name()<<endl; 
} catch (...) {
 cout<<"\ntsttoi2map: some other exception was caught ! "<<endl;
}
//--------------------------------------------------------------------

exit(0);
}