Changeset 560

Timestamp:

Oct 7, 2011, 3:25:15 PM (13 years ago)

Author:

campagne

Message:

start merging with calibration options (jec)

File:

• BAORadio/AmasNancay/trunk/mergeAnaFiles.cc (modified) (12 diffs)

BAORadio/AmasNancay/trunk/mergeAnaFiles.cc

@@ -41 +41 @@
 const sa_size_t NUMBER_OF_FREQ     = 8192;
 const sa_size_t TOTAL_NUM_CYCLES   = 500;
-const r_4    LOWER_FREQUENCY = 1250.0; //MHz
-const r_4    TOTAL_BANDWIDTH = 250.0; //MHz
+const r_4    LOWER_FREQUENCY       = 1250.0; //MHz
+const r_4    TOTAL_BANDWIDTH       = 250.0;  //MHz
 //-----------------------------------------------
 //Input parameters
@@ -52 +52 @@
   string ppfFile_;    //generic name of the input files
   int nSliceInFreq_;  //used by reduceSpectra() fnc
+  string typeOfCalib_;//type of calibration to be done
 } para;
 //--------------------------------------------------------------
@@ -207 +208 @@
 }
 //-------------------------------------------------------
-//Compute the mean of Raw spectra and also the mean/sigma of rebinned spectra 
+//Compute the mean of Diff ON-OFF BAO-calibrated spectra and also the mean/sigma of rebinned spectra 
+//Used like: 
+//
+void meanCalibBAODiffOnOffCycles() throw(string) {
+
+  list<string> listOfFiles;
+  string directoryName;
+  directoryName = para.inPath_ + "/" + para.sourceName_;
+
+  //Make the listing of the directory
+  listOfFiles = ListOfFileInDir(directoryName,para.ppfFile_);
+  
+  list<string>::const_iterator iFile, iFileEnd, iSpecOff, iSpecOffEnd, iSpecOn, iSpecOnEnd;
+  iFileEnd = listOfFiles.end();
+  
+  //Loop on files
+  for (iFile = listOfFiles.begin(); iFile != iFileEnd; ++iFile) {
+    if (para.debuglev_>90){
+      cout << "load file <" << *iFile << ">" << endl;
+    }
+    PInPersist fin(*iFile);
+    vector<string> vec = fin.GetNameTags();
+
+    if (para.typeOfCalib_ == "perRun") {
+      ///////////////////
+      //make the calibration of the mean of all Off and On of the run and perform the difference
+      ///////////////////
+      //
+      //Compute the mean of the Off
+      //
+      list<string> listOfSpectraOff;
+      //Keep only required PPF objects
+      std::remove_copy_if(
+                          vec.begin(), vec.end(), back_inserter(listOfSpectra),
+                          not1(match)
+                          );
+      
+      listOfSpectraOff.sort(stringCompare);
+      iSpecOffEnd = listOfSpectraOff.end();
+      TMatrix<r_4> meanOfSpectraOff(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ);
+      uint_4 nSpectraOff=0;
+      //Loop of spectra matrix
+      for (iSpecOff = listOfSpectraOff.begin(); iSpecOff !=iSpecOffEnd;  ++iSpecOff){
+        if (para.debuglev_>90){
+          cout << " spactra <" << *iSpecOff << ">" << endl;
+        }
+        TMatrix<r_4> aSpec(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ);
+        fin.GetObject(aSpec,*iSpecOff);
+        //How to see if the GetObject is ok?? Ask Reza
+        nSpectraOff++;
+        meanOfSpectraOff+=aSpec;
+      }//end loop Off
+      //Normalisation
+      if(nSpectraOff>0)meanOfSpectraOff/=(r_4)(nSpectraOff);
+      //
+      //Compute the mean of the On
+      //
+      list<string> listOfSpectraOn;
+      //Keep only required PPF objects
+      std::remove_copy_if(
+                          vec.begin(), vec.end(), back_inserter(listOfSpectra),
+                          not1(match)
+                          );
+      
+      listOfSpectraOn.sort(stringCompare);
+      iSpecOnEnd = listOfSpectraOn.end();
+      TMatrix<r_4> meanOfSpectraOn(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ);
+      uint_4 nSpectraOn=0;
+      //Loop of spectra matrix
+      for (iSpecOn = listOfSpectraOn.begin(); iSpecOn !=iSpecOnEnd;  ++iSpecOn){
+        if (para.debuglev_>90){
+          cout << " spactra <" << *iSpecOn << ">" << endl;
+        }
+        TMatrix<r_4> aSpec(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ);
+        fin.GetObject(aSpec,*iSpecOn);
+        //How to see if the GetObject is ok?? Ask Reza
+        nSpectraOn++;
+        meanOfSpectraOn+=aSpec;
+      }//end loop On
+      //Normalisation
+      if(nSpectraOn>0)meanOfSpectraOn/=(r_4)(nSpectraOn);
+      //
+      //Get the calibration file
+      //
+      
+
+
+    } else if {para.typeOfCalib_ == "perCycle") {
+      //perform the calibration of the OFF and ON per cycle, then make the mean and take the diff
+    } else {
+      string msg="FATAL (meanCalibBAODiffOnOffCycles); unknown calibration mode " 
+        + para.typeOfCalib_ ;
+      throw(msg);
+    }
+
+    
+    }//eo loop on spectra in a file
+  }//eo loop on files
+  
+}
+//-------------------------------------------------------
+//Compute the mean of Diff ON-OFF Raw spectra and also the mean/sigma of rebinned spectra 
 //Used like:
 //
-void meanOnCycles() throw(string) {
+void meanRawDiffOnOffCycles() throw(string) {
   list<string> listOfFiles;
   string directoryName;
@@ -238 +340 @@
                         );
    
+    listOfSpectra.sort(stringCompare);
     iSpecEnd = listOfSpectra.end();
-    listOfSpectra.sort(stringCompare);
     //Loop of spectra matrix
     for (iSpec = listOfSpectra.begin(); iSpec !=iSpecEnd;  ++iSpec){
@@ -277 +379 @@
 }
 //-------------------------------------------------------
-//Compute the median of Raw spectra and also the mean/sigma of rebinned spectra 
+//Compute the median of Diff ON-OFF Raw spectra and also the mean/sigma of rebinned spectra 
 //Used like:
 //
-void medianOnCycles() throw(string) {
+void medianRawDiffOnOffCycles() throw(string) {
   list<string> listOfFiles;
   string directoryName;
@@ -310 +412 @@
                         );
    
+    listOfSpectra.sort(stringCompare);
     iSpecEnd = listOfSpectra.end();
-    listOfSpectra.sort(stringCompare);
     //Loop of spectra matrix
     for (iSpec = listOfSpectra.begin(); iSpec !=iSpecEnd && (sa_size_t)nSpectra < TOTAL_NUM_CYCLES ;  ++iSpec){
@@ -375 +477 @@
   }
   TVector<r_4>sigmaMed(NUMBER_OF_CHANNELS);
-  sa_size_t redf1320=(1320.0-LOWER_FREQUENCY)/TOTAL_BANDWIDTH*para.nSliceInFreq_;
-  sa_size_t redf1380=(1380.0-LOWER_FREQUENCY)/TOTAL_BANDWIDTH*para.nSliceInFreq_;
+  sa_size_t redf1320=(sa_size_t)((1320.0-LOWER_FREQUENCY)/TOTAL_BANDWIDTH*para.nSliceInFreq_);
+  sa_size_t redf1380=(sa_size_t)((1380.0-LOWER_FREQUENCY)/TOTAL_BANDWIDTH*para.nSliceInFreq_);
 
   for (sa_size_t iCh=0; iCh<NUMBER_OF_CHANNELS; ++iCh){
@@ -442 +544 @@
   string ppfFile;
   string nSliceInFreq = "32";
- 
+  string typeOfCalib="perRun";
   
   //  bool okarg=false;
@@ -448 +550 @@
   while (ka<narg) {
     if (strcmp(arg[ka],"-h")==0) {
-      cout << "Usage:  -act [meanDiffOnOff]|medianDiffOnOff\n"
+      cout << "Usage:  -act [meanRawDiffOnOff]|medianRawDiffOnOff|meanCalibBAODiffOnOff\n"
+           << " -calibopt [perRun]|perCcyle (cf. in case of calibration to be done)"
            << " -src [Abell85]\n -inPath [.]|<top_root_dir of the ppf file>\n" 
            << " (ex. /sps/baoradio/AmasNancay/JEC/\n " 
@@ -466 +569 @@
       ka+=2;
     }
+    else if (strcmp(arg[ka],"-calibopt")==0) {
+      typeOfCalib=arg[ka+1];
+      ka+=2;
+    }    
     else if (strcmp(arg[ka],"-inPath")==0) {
       inputPath=arg[ka+1];
@@ -507 +614 @@
 
   if ( "" == ppfFile ) {
-    cerr << "mergeRawOnOff.cc: you have forgotten ppfFile option"
+    cerr << "mergeAnaFiles.cc: you have forgotten ppfFile option"
          << endl;
     return 999;
@@ -519 +626 @@
 //     printf("->%s<-\n(b=%d e=%d)\n",match,b,e);
 
-    if ( action == "meanDiffOnOff" ) {
-      meanOnCycles();
-    } else if (action == "medianDiffOnOff") {
-      medianOnCycles();
+    if ( action == "meanRawDiffOnOff" ) {
+      meanRawDiffOnOffCycles();
+    } else if (action == "medianRawDiffOnOff") {
+      medianRawDiffOnOffCycles();
+    } else if (action == "meanCalibBAODiffOnOff") {
+      meanCalibBAODiffOnOffCycles();
     } else {
       msg = "Unknown action " + action;