Changeset 1516 in Sophya for trunk/ArchTOIPipe


Ignore:
Timestamp:
Jun 12, 2001, 6:00:35 PM (24 years ago)
Author:
cmv
Message:

Re-shape de toi2map cmv 12/6/01

Location:
trunk/ArchTOIPipe
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • trunk/ArchTOIPipe/Kernel/flagtoidef.h

    r1502 r1516  
     1#ifndef FLAGTOIDEF_H
     2#define FLAGTOIDEF_H
    13
    24enum FlagToiDef {
    3 // Out of Range Valuet of Range Value
     5// Out of Range Value of Range Value
    46  FlgToiOut     =  (unsigned long long)1 << 0,
    57// Spike-like sample (ex: glitch)
    6   FlgToiSpike   =  (unsigned long long)1 << 2,
     8  FlgToiSpike   =  (unsigned long long)1 << 1,
    79// Discontinuity-like sample (ex: step)
    8   FlgToiDisc    =  (unsigned long long)1 << 3,
     10  FlgToiDisc    =  (unsigned long long)1 << 2,
    911// Large-Sigma-like sample
    10   FlgToiSigma   =  (unsigned long long)1 << 4,
     12  FlgToiSigma   =  (unsigned long long)1 << 3,
    1113// Cleaned sample (ex: clean after spike,step...)
    12   FlgToiClean   =  (unsigned long long)1 << 5,
     14  FlgToiClean   =  (unsigned long long)1 << 4,
    1315// Sample Value has been interpolated/Replaced ...
    14   FlgToiInterp  =  (unsigned long long)1 << 6,
     16  FlgToiInterp  =  (unsigned long long)1 << 5,
     17};
     18
    1519// et on peut continuer ainsi jusqu'a 1<<63 cad le bit 64 ieme
    1620// Pour imprimer la valeur: cout<<(unsigned long long)FlgToi...<<endl;
    17 }
    1821
     22#endif
  • trunk/ArchTOIPipe/ProcWSophya/toi2map.cc

    r1498 r1516  
    44#include "ctimer.h"
    55#include "toi2map.h"
    6 #include "xastropack.h"
    76
    87////////////////////////////////////////////////////////////////////////
    98TOI2Map::TOI2Map(SphereHEALPix<r_8>* sph,SphereHEALPix<r_8>* wsph)
    10   : mSph(sph), mWSph(wsph), mWSphInternal(false), mTypCoor(false), fTypCoor(false), mActualYear(2001.)
     9  : mSph(sph), mWSph(wsph), mWSphInternal(false)
    1110{
     11 SetEquinox();
     12 SetCoorIn();
     13 SetCoorOut();
     14
    1215 if(mSph->NbPixels()<1) {
    1316  cout<<"TOI2Map::TOI2Map() Bad number of pixels in sphere mSph "
     
    4447  declareInput("Coord1In");     // input index 0
    4548  declareInput("Coord2In");     // input index 1
    46   declareInput("BoloIn");      // input index 2
     49  declareInput("BoloIn");       // input index 2
    4750}
    4851
     
    6164  throw ParmError("TOI2Map::run() Output TOI (Coord1In or Coord2In or BoloIn) not connected!");
    6265}
     66if( !(mTypCoorIn&TypCoordEq || mTypCoorIn&TypCoordGal) ) {
     67  cout<<"TOI2Map::run() - Input Coordinates not Eq or Gal! "<<endl;
     68  throw ParmError("TOI2Map::run() - Input Coordinates not Eq or Gal!");
     69}
     70if( !(mTypCoorOut&TypCoordEq || mTypCoorOut&TypCoordGal) ) {
     71  cout<<"TOI2Map::run() - Output Coordinates not Eq or Gal! "<<endl;
     72  throw ParmError("TOI2Map::run() - Output Coordinates not Eq or Gal!");
     73}
    6374
    6475//---------------------------------------------------------
     
    6980uint_4 mNSnFill=0, mNpixFill=0, NFill[NFILL];
    7081for(ii=0;ii<NFILL;ii++) NFill[ii]=0;
     82double mjd = MJDfrYear(mActualYear);
    7183
    7284// Remplissage des spheres
    7385for(int s=snb;s<=sne;s++) {
    74   int_8 fgbolo = 0; 
     86  int_8 fgbolo = 0;
    7587  double bolo;
     88  //              Equatoriales   /   Galactiques
     89  // coord1,2 =   alpha,delta    /   gLon,gLat
     90  double coord1 = getData(0,s);
     91  double coord2 = getData(1,s);
    7692
    77   double coord1 = getData(0,s); // gLat ou delta entre [-90,90] en degres
    78   double coord2 = getData(1,s); // gLon entre [0,360[ en degres ou alpha entre [0,24[ en heures
    7993  getData(2,s,bolo,fgbolo);
     94  if(bolo<-32767.) continue; // Bidouille (A CORRIGER QUAND FGBOLO DISPO)
    8095
    81   if(coord2<-90. || coord2>90.) fgbolo=1;
    82   if((coord1<0.) || (!mTypCoor && coord1>=24.) || (mTypCoor && coord1>=360.) ) {fgbolo=1;
    83   cout << "!!!!!!!!" <<coord1 << endl;
     96  // sphere phi   entre [0,2*Pi] en radians
     97  // sphere theta entre [0,Pi]   en radians
     98  double phi,theta;
     99  CoordConvertToStd(mTypCoorIn,coord1,coord2);
     100  if(mTypCoorIn&TypCoordEq && mTypCoorOut&TypCoordGal) { // Eq -> Gal
     101    EqtoGal(mjd,coord1,coord2,&coord1,&coord2);
     102    phi   = coord1 * M_PI/180.;
     103  } else if(mTypCoorIn&TypCoordGal && mTypCoorOut&TypCoordEq) { // Gal -> Eq
     104    GaltoEq(mjd,coord1,coord2,&coord1,&coord2);
     105    phi   = coord1 * M_PI/12.;
     106  } else if(mTypCoorOut&TypCoordGal) { // Gal -> Gal
     107    phi   = coord1 * M_PI/180.;
     108  } else if(mTypCoorOut&TypCoordEq) { // Eq -> Eq
     109    phi   = coord1 * M_PI/12.;
    84110  }
    85  
    86  
    87  
    88   if(bolo<-32767.) fgbolo=1;   // Bidouille Archeops
    89  
    90   if(!fgbolo) {
    91     // sphere phi entre [0,2*Pi] en radians
    92     // sphere phi entre [0,2*Pi] en radians
    93     // sphere theta entre [0,Pi] en radians
    94     double phi,theta;
    95     if(fTypCoor && !mTypCoor) { //on a alpha,delta et on veut l,b
    96       double mjd = MJDfrYear(mActualYear);
    97       EqtoGal(mjd,coord1,coord2,&coord1,&coord2);
    98     }
     111  theta = (90.-coord2) * M_PI/180.;
     112  if(phi<0.   || phi>=2*M_PI) continue;
     113  if(theta<0. || theta>=M_PI) continue;
    99114
    100     if(!fTypCoor && mTypCoor) { //on a l,b et on veut alpha,delta
    101       double mjd = MJDfrYear(mActualYear);
    102       GaltoEq(mjd,coord1,coord2,&coord1,&coord2);
    103     }
    104    
    105     if(fTypCoor) phi = coord1 * M_PI/180.;
    106     else  phi = coord1 * M_PI/12.;
    107 
    108     theta = (90.-coord2)*M_PI/180.;
    109 
    110     int_4 ipix = mSph->PixIndexSph(theta,phi);
    111     (*mSph)(ipix) += bolo;   
    112     ((*mWSph)(ipix))++;
    113     mNSnFill++;
    114   }
     115  int_4 ipix = mSph->PixIndexSph(theta,phi);
     116  (*mSph)(ipix) += bolo;   
     117  ((*mWSph)(ipix))++;
     118  mNSnFill++;
    115119}
    116120
     
    129133     <<"  mNpixFill="<<mNpixFill
    130134     <<"  mNSnFill="<<mNSnFill<<endl
    131      <<" --> FracSky="<<mNpixFill*100./(double)mSph->NbPixels()<<"%"<<endl
    132      <<"NFill["<<NFILL<<"] = "<<endl;
    133  for(ii=0;ii<NFILL;ii++) cout<<NFill[ii]<<" ";
     135     <<"  --> FracSky="<<mNpixFill*100./(double)mSph->NbPixels()<<"%"
     136     <<"  NFill["<<NFILL<<"] ="<<endl;
     137 for(ii=0;ii<NFILL;ii++) {cout<<NFill[ii]<<" "; if(ii%10==9) cout<<endl;}
    134138 cout<<endl;
    135139
  • trunk/ArchTOIPipe/ProcWSophya/toi2map.h

    r1498 r1516  
    55#include "toiprocessor.h"
    66#include "spherehealpix.h"
     7#include "flagtoidef.h"
     8#include "xastropack.h"
    79
    810//-- Un projecteur de TOI sur une sphere Healpix
    9 // Lecture de 3 TOI alpha,delta,boloMuV
    10 // Sortie rien
     11// Lecture de 3 TOI coord1,coord2,boloMuV
     12// Sortie pas de TOI, une (2) sphere(s) Healpix
    1113//
    1214// Structure generale :
    13 // (les alpha[0,360[,delta[-90,90] sont en degres decimaux)
    14 //                         |---->Sphere
    15 //                         |
    16 //                    -----------
     15//                          |---->Sphere
     16//                          |
     17//                     -----------
    1718//   toi Coord1In ---> |         |
    1819//   toi Coord2In ---> | TOI2Map |
    19 //   toi BoloIn  ---> |         |
    20 //                    -----------
     20//   toi BoloIn   ---> |         |
     21//                     -----------
     22// Gestion du type de coordonnees :
     23// Coord1In,Coord2In : soit Equatoriales (Alpha,Delta)
     24//                     soit Galactiques  (GLong,GLat)
     25// Sortie sur une sphere en coordonnees Equatoriales ou Galactiques
    2126
    2227class TOI2Map : public TOIProcessor {
     
    2530  virtual       ~TOI2Map();
    2631
    27   virtual void  init(void); 
     32  virtual void  init(void);
    2833  virtual void  run(void);
    29   inline  void  SetCoorGal(bool mfg=false,bool ffg=false,double actualyear=2000.)
    30                           {mTypCoor = mfg; fTypCoor = ffg; mActualYear = actualyear;}
    31  
     34
     35  // Coordonnees donnees en entree et en sortie
     36  inline void SetEquinox(double actualyear=2000.)
     37              {mActualYear = actualyear;}
     38  inline void SetCoorIn(TypAstroCoord mfg=TypCoordGalStd)
     39              {mTypCoorIn = mfg;}
     40  inline void SetCoorOut(TypAstroCoord mfg=TypCoordGalStd)
     41              {mTypCoorOut = mfg;}
     42
    3243protected:
    3344  SphereHEALPix<r_8>* mSph;
    3445  SphereHEALPix<r_8>* mWSph;
    3546  bool mWSphInternal;
    36   bool mTypCoor;
    37   bool fTypCoor;
     47  TypAstroCoord mTypCoorIn, mTypCoorOut;
    3848  double mActualYear;
    3949};
  • trunk/ArchTOIPipe/SMakefile

    r1496 r1516  
    1515            $(AOBJ)simtoipr.o $(AOBJ)map2toi.o $(AOBJ)toi2map.o $(AOBJ)nooppr.o
    1616EXEOLIST := $(AOBJ)fits2asc.o $(AOBJ)tsttoi.o $(AOBJ)tsttoi2.o $(AOBJ)simtst.o \
    17             $(AOBJ)tstrztoi.o $(AOBJ)rztoi.o \
     17            $(AOBJ)tstrztoi.o $(AOBJ)rztoi.o $(AOBJ)sfftc.o $(AOBJ)mesovh.o \
    1818            $(AOBJ)fits2ascii.o $(AOBJ)tgenw.o $(AOBJ)tstmap2toi.o $(AOBJ)tsttoi2map.o
    1919EXELIST := $(AOBJ)sfftc $(AOBJ)simtst $(AOBJ)mesovh $(AOBJ)tstrztoi $(AOBJ)fits2asc $(AOBJ)tsttoi \
  • trunk/ArchTOIPipe/TestPipes/tsttoi2map.cc

    r1498 r1516  
    1919 cout<<"tsttoi2map [-h] [-p lp] [-s samplemin,samplemax] [-w data_window_size]"<<endl
    2020     <<"           [-a label_coord1] [-d label_coord2] [-b label_bolomuv]"<<endl
    21      <<"           [-n nlat] coord_ini coord_fin fitsin_bolo fitsin_point fitsphout fitsphwout"<<endl
    22      <<" coord_ini[_fin] = G for Galactic coordinates"<<endl
    23      <<"                 = Q for equatorial coordinates"<<endl
    24      <<" coord1 = alpha or longitude ; coord2 = delta or latitude"<<endl;
    25  return;
     21     <<"           [-n nlat] [-i c,h] [-o c,h]"<<endl
     22     <<"           fitsin_point fitsin_bolo fitsphout [fitsphwout]"<<endl
     23     <<" -p lp : print level (def=0)"<<endl
     24     <<" -s samplemin,samplemax : sample range to be treated (def=all)"<<endl
     25     <<" -w data_window_size : window size for pipe (def=8192)"<<endl
     26     <<" -a label_coord1 : label fits for alpha/gLong (def=coord1)"<<endl
     27     <<" -d label_coord2 : label fits for delta/gLat (def=coord2)"<<endl
     28     <<"          coord1 = alpha or gLong ; coord2 = delta or gLat"<<endl
     29     <<" -b label_bolomuv : label fits for bolo value (def=boloMuV)"<<endl
     30     <<" -n nlat : nlat for Healpix sphere (def=128)"<<endl
     31     <<" -i c,h : coord1 caracteristics (c=G/E h=H/D) (def=G,D)"<<endl
     32     <<" -o c,h : idem -i for coord2"<<endl
     33     <<" fitsin_point : fits file for pointing"<<endl
     34     <<" fitsin_bolo : fits file for bolo values"<<endl
     35     <<" fitsphout : fits file for output Healpix sphere"<<endl
     36     <<" fitsphwout : fits file for output Healpix nFilled sphere (def=no)"<<endl;
     37}
     38
     39unsigned long typecoord(char typc=' ',char hd=' ');
     40unsigned long typecoord(char typc,char hd)
     41// typc : G=galactiques, E=equatoriales, autres=galactiques
     42// hd : H=heure, D=degre, autres=(heure si typc==E, degre si typc==G)
     43{
     44  if(typc!='G' && typc!='E') typc='G';
     45  if(hd!='H' && hd!='D') {if(typc=='E') hd='H'; else hd='D';}
     46 unsigned long rc=TypCoordUndef;
     47  if(typc=='G') rc |= TypCoordGal;
     48    else        rc |= TypCoordEq;
     49  if(hd=='D')   rc |= TypCoordDD;
     50    else        rc |= TypCoordHD;
     51  return rc;
    2652}
    2753
     
    3662char *label_coord1 = "coord1", *label_coord2 = "coord2", *label_bolomuv = "boloMuV";
    3763long sdeb,sfin;
    38 int c;
    39 while((c = getopt(narg,arg,"hGp:s:w:a:d:b:n:")) != -1) {
     64string fitsphwout = "";
     65unsigned long tcoorin=typecoord(), tcoorout=typecoord();
     66int c; char t=' ',h=' ';
     67while((c = getopt(narg,arg,"hp:s:w:a:d:b:n:i:o:")) != -1) {
    4068  switch (c) {
    4169  case 's' :
     
    4371    cout<<"Requested Samples from "<<sdeb<<" , "<<sfin<<endl;
    4472    if(sfin>=sdeb) mgr->setRequestedSample(sdeb,sfin);
    45     else {cout<<"Bad sample interval "<<endl; exit(-2);}
     73    else {cout<<"Bad sample interval "<<endl; exit(2);}
    4674    break;
    4775  case 'w' :
     
    6795    if(nlat<0) nlat=128;
    6896    break;
     97  case 'i' :
     98    sscanf(optarg,"%c,%c",&t,&h);
     99    tcoorin=typecoord(t,h);
     100    break;
     101  case 'o' :
     102    sscanf(optarg,"%c,%c",&t,&h);
     103    tcoorout=typecoord(t,h);
     104    break;
    69105  case 'h' :
    70     usage(); exit(-1);
    71     break;
    72106  default:
    73     usage(); exit(-1);
     107    usage(); exit(1);
     108    break;
    74109  }
    75110}
    76 if(optind+3>=narg) {usage(); exit(-2);}
    77 bool mcoorgal=true;
    78 bool fcoorgal=true;
     111if(optind+2>=narg) {usage(); exit(3);}
    79112 
    80 if ( strcmp(arg[optind],"Q") == 0) mcoorgal=false;
    81 if ( strcmp(arg[optind+1],"Q") == 0) fcoorgal=false;
    82  
    83 char * fitsin_bolo = arg[optind+2];
    84 char * fitsin_point = arg[optind+3];
    85 string const fitsphout = arg[optind+4];
    86 string fitsphwout = "";
    87 if(optind+5<narg) fitsphwout = arg[optind+5];
    88 
     113char * fitsin_point = arg[optind];
     114char * fitsin_bolo = arg[optind+1];
     115string const fitsphout = arg[optind+2];
     116if(optind+3<narg) fitsphwout = arg[optind+3];
     117
     118{
     119unsigned long tg,te,hd,dd;
    89120cout<<">>>> tsttoi2map:"<<endl
    90     <<"Fits Infile(snum,bolomuv)= "<<fitsin_bolo<<endl
    91     <<"Fits Infile(snum,coord1,coord2)= "<<fitsin_point<<endl
    92     <<"  ...label_coord1 "<<label_coord1<<"  ,  label_coord2 "<<label_coord2<<endl
    93     <<"  with coordinates Gal "<<mcoorgal<<endl
    94     <<"  ...label_bolomuv "<<label_bolomuv<<endl
    95     <<"Fits Sphere Healpix"<<fitsphout<<endl
     121    <<"Pipe Window Size "<<width<<endl
     122    <<"Fits Infile Bolo "<<fitsin_bolo<<endl
     123    <<"  ...label_bolomuv "<<label_bolomuv<<endl;
     124tg = tcoorin&TypCoordGal; te = tcoorin&TypCoordEq;
     125hd = tcoorin&TypCoordHD;  dd = tcoorin&TypCoordDD;
     126cout<<"  ...label_coord1 "<<label_coord1<<endl
     127    <<"  ...label_coord2 "<<label_coord2<<endl
     128    <<"  ...... Gal="<<tg<<" Eq="<<te<<"   hour="<<hd<<" deg="<<dd<<endl;
     129tg = tcoorout&TypCoordGal; te = tcoorout&TypCoordEq;
     130hd = tcoorout&TypCoordHD;  dd = tcoorout&TypCoordDD;
     131cout<<"Fits Healpix Sphere "<<fitsphout<<endl
    96132    <<"  ...nlat "<<nlat<<endl
    97     <<"  with coordinates Gal "<<fcoorgal<<endl
    98     <<"Fits Sphere Healpix Error "<<fitsphwout<<endl;
     133    <<"  ...... Gal="<<tg<<" Eq="<<te<<"   hour="<<hd<<" deg="<<dd<<endl;
     134cout<<"Fits Healpix Weight Sphere "<<fitsphwout<<endl;
     135}
    99136
    100137SophyaInit();
     
    109146 int ncolb = rfitsb.getNOut();
    110147 cout<<"Number of columns in fits Infile_bolo : "<<ncolb<<endl;
    111  if(ncolb<1) exit(-3);
     148 if(ncolb<1) exit(-4);
    112149
    113150 FITSTOIReader rfitsp(fitsin_point);
    114151 int ncolp = rfitsp.getNOut();
    115152 cout<<"Number of columns in fits Infile_point : "<<ncolp<<endl;
    116  if(ncolp<2) exit(-3);
     153 if(ncolp<2) exit(-5);
    117154
    118155 // Creation de la sphere Healpix
     
    126163   wsph = new SphereHEALPix<r_8>;
    127164   cout<<"SphereHEALPix Weight: Type de map : "<<wsph->TypeOfMap()<<endl
    128        <<"               Nombre de pixels : "<<wsph->NbPixels()<<endl;
     165       <<"              Nombre de pixels : "<<wsph->NbPixels()<<endl;
    129166 }
    130167
     
    132169 TOI2Map toi2m(sph,wsph);
    133170 cout<<"TOI2Map created"<<endl;
    134  toi2m.SetCoorGal(mcoorgal,fcoorgal,2000.);  // equinoxe de ref. 2000. (pour archtoi)
     171 toi2m.SetEquinox(2000.);
     172 toi2m.SetCoorIn((TypAstroCoord) tcoorin);
     173 toi2m.SetCoorOut((TypAstroCoord) tcoorout);
    135174
    136175 // Definition des tuyaux
    137176 TOISeqBuffered * toicoord1in = new TOISeqBuffered("toi_coord1_in",width);
    138  if(lp) toicoord1in->setDebugLevel(1);
     177 // toicoord1in->setDebugLevel(1);
    139178 rfitsp.addOutput(label_coord1,toicoord1in);
    140179 toi2m.addInput("Coord1In",toicoord1in);
    141180
    142181 TOISeqBuffered * toicoord2in = new TOISeqBuffered("toi_coord2_in",width);
    143  if(lp) toicoord2in->setDebugLevel(1);
     182 // toicoord2in->setDebugLevel(1);
    144183 rfitsp.addOutput(label_coord2,toicoord2in);
    145184 toi2m.addInput("Coord2In",toicoord2in);
    146185 
    147186 TOISeqBuffered * toibolin = new TOISeqBuffered("toi_bolo_in",width);
    148  if(lp) toibolin->setDebugLevel(1);
     187 // toibolin->setDebugLevel(1);
    149188 rfitsb.addOutput(label_bolomuv,toibolin);
    150189 toi2m.addInput("BoloIn",toibolin);
Note: See TracChangeset for help on using the changeset viewer.