Changeset 3886 in Sophya

Timestamp:

Sep 24, 2010, 12:31:21 AM (15 years ago)

Author:

ansari

Message:

Codage remplissage n_tuple avec les variances, Reza 23/09/2010

Location:

trunk/AddOn/TAcq

Files:

• branap.cc (modified) (5 diffs)
• branap.h (modified) (1 diff)
• brproc.cc (modified) (9 diffs)
• brproc.h (modified) (3 diffs)
• specmfib.cc (modified) (1 diff)

trunk/AddOn/TAcq/branap.cc

@@ -31 +31 @@
   nbcalgrp_=1;
   nthreads_=1;
+
+  vmin_=0.;  vmax_=9e99;
 }
 
@@ -86 +88 @@
       rdsamefc_=false;
       ka++;
+    }
+    else if (strcmp(arg[ka],"-varcut")==0) {
+      sscanf(arg[ka+1],"%lg,%lg",&vmin_,&vmax_);
+      ka+=2;
     }
     else if (strcmp(arg[ka],"-in")==0) {
@@ -126 +132 @@
   cout << " Usage:  prgname [-act ACT] [-out OutPath] [-nmean NMean] [-zones NZones,nPaqinZone] \n" 
        << "                 [-nbloc NBloc] [-freq NumFreqMin,NumFreqMax,NBinFreq] \n" 
-       << "                 [-prt lev,modulo] [-nvcal n] [-nthr n] [-nosfc] \n" 
+       << "                 [-prt lev,modulo] [-nvcal n] [-nthr n] [-nosfc] [-varcut min,max] \n" 
        << "         -in Imin,Imax,Istep InPath FiberList [InPath2 FiberList2 InPath3 FiberList3 ...] \n" << endl; 
   if (fgshort) {
@@ -145 +151 @@
        << " -nthr n : number of threads for parallel execution in BRVisibilityCalculator (default=1) \n" 
        << " -nosfc : Don't force reading with SAME FrameCounter \n" 
+       << " -varcut min,max : min-max cut on variance \n" 
        << " -in : input files/directory definition : \n"
        << "       Imin,Imax,Istep: fits files signalII.fits Imin<=II<=Imax Istep=increment \n" 
@@ -172 +179 @@
   cout << " Action=" << action_ << "  NMean=" << nmean_ << " NBloc=" << nbloc_ << endl;
   cout << " FreqMin= " << freqmin_ << " FreqMax= " << freqmax_ << " NBinFreq= " << nbinfreq_ << endl;
+  cout << " VarianceMin= " << vmin_ << " VarianceMax= " << vmax_ << endl;
   cout << " PaqSize=" << paqsize_ << "  -  NZones=" << nzones_ << " NPaqZone=" << npaqinzone_ 
        << " PrtLevel=" << prtlevel_ << " PrtCntModulo=" <<  prtmodulo_ << endl;

trunk/AddOn/TAcq/branap.h

@@ -45 +45 @@
   uint_4 nbcalgrp_;    // Nb d'objets/threads dans BRVisCalcGroup 
   uint_4 nthreads_;  // Nb de threads pour l'execution parallele ds BRVisibilityCalculator
+
+  double vmin_,vmax_;   // coupure min-max sur la variance 
 };
 

trunk/AddOn/TAcq/brproc.cc

@@ -20 +20 @@
 #include "ctimer.h"
 #include "fftpserver.h"
-#include "fftwserver.h"
+#include "fitsarrhand.h"
 
 #include "FFTW/fftw3.h"
@@ -38 +38 @@
                                            bool fgdatafft, bool fgsinglechan)
   : BRBaseProcessor(memgr), outpath_(outpath), nmean_(nmean),
-    fgdatafft_(fgdatafft), fgsinglechannel_(fgsinglechan)
+    fgdatafft_(fgdatafft), fgsinglechannel_(fgsinglechan), 
+    clnflg_(fgsinglechan?memgr_.NbFibres():2*memgr_.NbFibres()),
+    nbadpaq_(fgsinglechan?memgr_.NbFibres():2*memgr_.NbFibres()) 
 {
   setNameId("meanSpecCalc",1);
@@ -55 +57 @@
   nbpaq4mean_=0;
   totnbpaq_=0;
+
+  SetVarianceLimits();
+
+  ofsdtp_=NULL;
+  dtp_=NULL;
+  string dtfile="!"+outpath_+"/dtspec.fits";
+  ofsdtp_ = new FitsInOutFile(dtfile,FitsInOutFile::Fits_Create);
+  dtp_ = new SwFitsDataTable(*ofsdtp_,1024,true);
+  int nchan=(fgsinglechannel_?memgr_.NbFibres():2*memgr_.NbFibres());
+  char cnom[32];
+  for(int i=0; i<nchan; i++) {
+    sprintf(cnom,"var%d",i);
+    dtp_->AddFloatColumn(cnom);    
+  }
+  xnt_=new double[nchan];
 }
 
@@ -61 +78 @@
 {
   if (nbpaq4mean_>1)  SaveSpectra();
+  cout << " ---------------- BRMeanSpecCalculator()_Finalizing -------------------- " << endl;
+  for(size_t i=0; i<nbadpaq_.size(); i++) {
+    cout << " Channel " << i << "  NBadPaq=" << nbadpaq_[i] << " / TotNbPaq=" << totnbpaq_ << endl;
+  }
+  if (dtp_) {
+    cout << *dtp_; 
+    delete dtp_;
+    delete ofsdtp_;
+    delete xnt_;
+  }
+  cout << " ------------------------------------------------------------------------ " << endl;
 }
 
@@ -70 +98 @@
 int BRMeanSpecCalculator::Process()
 {
-      
   if (nbpaq4mean_==nmean_)  SaveSpectra();
+  FlagBadPackets();
     
   if (fgdatafft_) {  // Donnees firmware FFT 
     for(sa_size_t i=0; i<(size_t)mspecmtx_.NRows(); i++) {
+      if (clnflg_[i])  { nbadpaq_[i]++;  continue; }   // si le paquet a ete flagge mauvais ( clnflg_[i] <> 0 ) 
       TwoByteComplex* zp=NULL;
       if (fgsinglechannel_) {
@@ -94 +123 @@
   else {  // Donnees RAW qui ont du etre processe par BRFFTCalculator
     for(sa_size_t i=0; i<(size_t)mspecmtx_.NRows(); i++) {
+      if (clnflg_[i])  { nbadpaq_[i]++;  continue; }   // si le paquet a ete flagge mauvais ( clnflg_[i] <> 0 ) 
       complex<ODT>* zp=NULL;
       if (fgsinglechannel_) {
@@ -116 +146 @@
 
 /* --Methode-- */
+void BRMeanSpecCalculator::FlagBadPackets()
+{
+  if (fgdatafft_) {  // Donnees firmware FFT 
+    for(sa_size_t i=0; i<(size_t)mspecmtx_.NRows(); i++) {
+      TwoByteComplex* zp=NULL;
+      if (fgsinglechannel_) {
+        zp=vpaq_[i].Data1C();
+      }
+      else {
+        zp=vpaq_[i/2].Data1C();
+        if (i%2==1)  zp=vpaq_[i/2].Data2C();
+      }
+      double variance=0.;
+      for(sa_size_t f=1; f<mspecmtx_.NCols(); f++)   variance+=(double)zp[f].module2F();
+      xnt_[i]=variance;
+      clnflg_[i]=0;
+      if (variance<varmin_) clnflg_[i]=1;
+      else if (variance>varmax_) clnflg_[i]=2;
+    }
+  }
+  else {  // Donnees RAW qui ont du etre processe par BRFFTCalculator
+    for(sa_size_t i=0; i<(size_t)mspecmtx_.NRows(); i++) {
+      complex<ODT>* zp=NULL;
+      if (fgsinglechannel_) {
+        zp=reinterpret_cast< complex<ODT>* > (vprocpaq_[i]);
+      }
+      else {
+        zp=reinterpret_cast< complex<ODT>* > (vprocpaq_[i/2]);
+        if (i%2==1)  zp= reinterpret_cast< complex<ODT>* >(vprocpaq_[i/2]+memgr_.ProcPaqSize()/2) ;
+      }
+      double variance=0.;
+      for(sa_size_t f=1; f<mspecmtx_.NCols(); f++)   variance+=(double)Zmod2(zp[f]);
+      xnt_[i]=variance;
+      clnflg_[i]=0;
+      if (variance<varmin_) clnflg_[i]=1;
+      else if (variance>varmax_) clnflg_[i]=2;     
+    }
+  }
+  dtp_->AddRow(xnt_);
+  return;
+}
+
+/* --Methode-- */
 void BRMeanSpecCalculator::SaveSpectra()
 {
@@ -125 +198 @@
   char nfile[64];
   string flnm;
+  /*
   {
   sprintf(nfile,"mspecmtx%d.ppf",numfile_);
@@ -137 +211 @@
   po << sigspecmtx_;
   }
+  */
+  {
+  sprintf(nfile,"mspecmtx%d.fits",numfile_);
+  flnm="!"+outpath_+nfile;
+  FitsInOutFile fos(flnm,FitsInOutFile::Fits_Create);
+  fos << mspecmtx_;
+  }
+  {
+  sprintf(nfile,"sigspecmtx%d.fits",numfile_);
+  flnm="!"+outpath_+nfile;
+  FitsInOutFile fos(flnm,FitsInOutFile::Fits_Create);
+  fos << sigspecmtx_;
+  }
+
   cout << numfile_ << "-BRMeanSpecCalculator::SaveSpectra() NPaqProc="
        << totnbpaq_ << "  -> Mean/Sig spectra Matrix in " << flnm << " /sigspec...ppf" << endl;

trunk/AddOn/TAcq/brproc.h

@@ -22 +22 @@
 #include "ntuple.h"
 #include "datatable.h"
+#include "swfitsdtable.h"
 #include "FFTW/fftw3.h"
 
@@ -39 +40 @@
   virtual ~BRMeanSpecCalculator(); // Fait le nettoyage final , sauvegarde des objets ...
 
+  inline void SetVarianceLimits(double vmin=0., double vmax=9.e99)
+     { varmin_=vmin;  varmax_=vmax; }
+
 protected:
   virtual int Process(); 
+  virtual void FlagBadPackets();
   virtual void SaveSpectra();
+
 
   uint_4 nmean_;  // Nombre de spectres pour le calcul des moyennes 
@@ -54 +60 @@
   TMatrix< r_4 > sigspecmtx_;    // Matrice des sigmas des spectres  
 
+  vector<int> clnflg_;      // resultat de la routine "cleaner" 
+  double varmin_, varmax_;  // limites pour la coupure sur la variance
+  vector<uint_8> nbadpaq_;
+
+  FitsInOutFile* ofsdtp_;   // Output fits stream for datatable 
+  SwFitsDataTable* dtp_;    // DataTable 
+  double* xnt_; 
 }; 
 

trunk/AddOn/TAcq/specmfib.cc

@@ -72 +72 @@
 
     BRMeanSpecCalculator procms(mmgr, par.outpath_, par.nmean_, par.fgdatafft_, par.fgsinglechannel_);
+    procms.SetVarianceLimits(par.vmin_, par.vmax_);
     BRFFTCalculator procfft(mmgr, par.fgsinglechannel_);
     if (!par.fgdatafft_) procms.SetMemZAction(MemZA_ProcB);