Changeset 552 for BAORadio/AmasNancay/trunk/mergeRawOnOff.cc

Timestamp:

Oct 5, 2011, 3:13:36 PM (13 years ago)

Author:

campagne

Message:

improvment (jec)

File:

• BAORadio/AmasNancay/trunk/mergeRawOnOff.cc (modified) (5 diffs)

BAORadio/AmasNancay/trunk/mergeRawOnOff.cc

@@ -55 +55 @@
 //--------------------------------------------------------------
 //Utility functions
+
+  class median_of_empty_list_exception:public std::exception{
+    virtual const char* what() const throw() {
+    return "Attempt to take the median of an empty list of numbers.  "
+      "The median of an empty list is undefined.";
+  }
+  };
+  template<class RandAccessIter>
+    double median(RandAccessIter begin, RandAccessIter end) 
+    throw(median_of_empty_list_exception){
+    if(begin == end){ throw median_of_empty_list_exception(); }
+    std::size_t size = end - begin;
+    std::size_t middleIdx = size/2;
+    RandAccessIter target = begin + middleIdx;
+    std::nth_element(begin, target, end);
+    
+    if(size % 2 != 0){ //Odd number of elements
+      return *target;
+    }else{            //Even number of elements
+      double a = *target;
+      RandAccessIter targetNeighbor= target-1;
+      std::nth_element(begin, targetNeighbor, end);
+      return (a+*targetNeighbor)/2.0;
+    }
+  }
+
 //--------------------------------------------------------------
 char *regexp (const char *string, const char *patrn, int *begin, int *end) {    
@@ -291 +317 @@
 
       if ((sa_size_t)nSpectra < TOTAL_NUM_CYCLES ) {
-        //STORE THE spectra
-        //      tableOfSpectra(Range::all(),Range::all(),Range(nSpectra)).Show();
-        //      aSpec.Show();
         tableOfSpectra(Range::all(),Range::all(),Range(nSpectra)) = aSpec;
-//      aSpec.Show();
-//      cout << "deb 1" << endl;
-//      tmp = aSpec;
-//      cout << "deb 2" << endl;
       } else {
         stringstream tmp;
@@ -310 +329 @@
   //Resize
 
-  nSpectra--;//For further selection in Range
-  if (para.debuglev_>90){
-    cout << "Cycles < 0," <<  nSpectra << ">" << endl;
-  }
  
-
   TMatrix<r_4> medianOfSpectra(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ);
   //Compute the median for each freq. and Channel
   for (sa_size_t iCh=0; iCh<NUMBER_OF_CHANNELS; iCh++){
     for (sa_size_t freq =0; freq<NUMBER_OF_FREQ; freq++){
-      TVector<r_4> tmp0(tableOfSpectra(Range(freq),Range(iCh),Range(0,nSpectra)).CompactAllDimensions());
+      TVector<r_4> tmp0(tableOfSpectra(Range(freq),Range(iCh),Range(0,nSpectra-1)).CompactAllDimensions());
       vector<r_4> tmp;
       tmp0.FillTo(tmp);
-      nth_element(tmp.begin(),tmp.begin()+tmp.size()/2,tmp.end());
-      medianOfSpectra(iCh,freq) = *(tmp.begin()+tmp.size()/2);
+      medianOfSpectra(iCh,freq) = median(tmp.begin(),tmp.end());
     }
   }
@@ -338 +351 @@
     tmp << nSpectra;
     string fileName = para.outPath_+"/medianDiffOnOffRaw_"+StringToLower(para.sourceName_)+"-"+tmp.str()+"Cycles.ppf";
-    cout << "Save median based on " <<  (nSpectra+1) << " cycles " << endl;
+    cout << "Save median based on " <<  nSpectra << " cycles " << endl;
     POutPersist fos(fileName);
 
@@ -347 +360 @@
     tag = "sigmamedred";
     fos << PPFNameTag(tag) << sigmaRedMtx;
+    tag = "cycleVsfreq";
+    
+    TArray<r_4> tarr(tableOfSpectra(Range::all(),Range::all(),Range(0,nSpectra-1)));
+    fos << PPFNameTag(tag) << tarr;
   }
 }