Changeset 3943 in Sophya for trunk

Timestamp:

Feb 1, 2011, 9:30:20 AM (15 years ago)

Author:

campagne

Message:

release.txt

Location:

trunk/AddOn/TAcq

Files:

• branap.cc (modified) (7 diffs)
• branap.h (modified) (2 diffs)
• brproc.cc (modified) (16 diffs)
• brproc.h (modified) (7 diffs)
• brprocGain.cc (added)
• brprocGain.h (added)
• makefile (modified) (3 diffs)
• specmfib.cc (modified) (5 diffs)

trunk/AddOn/TAcq/branap.cc

@@ -26 +26 @@
   nzones_=nzon;
   npaqinzone_=npaqz;
-  fgdatafft_=false;  fgsinglechannel_=false;
+  fgdatafft_=false;  fgsinglechannel_=false; fgforcesinglechan_=false;
   prtlevel_=0;
   prtmodulo_=50000;
@@ -37 +37 @@
 
   vmin_=0.;  vmax_=9e99;
+  nbands_=0; bandfirst_ = bandlast_ = 0;
+  fgdt_ = false;
+
+  fgfreqfilter_=false; //JEC 1/2/11
+
   gainfile_="";
 }
@@ -94 +99 @@
       ka++;
     }
+    else if (strcmp(arg[ka],"-singlechan")==0) {
+      fgsinglechannel_=true; fgforcesinglechan_ = true; 
+      ka++;
+    }
     else if (strcmp(arg[ka],"-varcut")==0) {
       sscanf(arg[ka+1],"%lg,%lg",&vmin_,&vmax_);
       ka+=2;
     }
+    else if (strcmp(arg[ka],"-nband")==0) {
+      sscanf(arg[ka+1],"%d,%d,%d",&nbands_,&bandfirst_,&bandlast_);
+      ka+=2;
+    }
+    else if (strcmp(arg[ka],"-filldt")==0) {
+      fgdt_=true;
+      ka++;
+    }
+    else if (strcmp(arg[ka],"-freqfilter")==0) {
+      fgfreqfilter_=true;
+      ka++;
+    }    
     else if (strcmp(arg[ka],"-gain")==0) {
       gainfile_=arg[ka+1];
@@ -147 +168 @@
   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" 
-       << "                 [-gain filename] [-varcut min,max] [-tspwin wsz,extsz,nfiles] \n" 
+       << "                 [-prt lev,modulo] [-nvcal n] [-nthr n] [-nosfc] [-singlechan]\n"
+       << "                 [-freqfilter]\n"
+       << "                 [-gain filename] [-varcut min,max] [-nband nband,first,last] \n"
+       << "                 [-filldt] [-tspwin wsz,extsz,nfiles] \n" 
        << "         -in Imin,Imax,Istep InPath FiberList [InPath2 FiberList2 InPath3 FiberList3 ...] \n" << endl; 
   if (fgshort) {
@@ -168 +191 @@
        << " -nthr n : number of threads for parallel execution in BRVisibilityCalculator (default=1) \n" 
        << " -nosfc : Don't force reading with SAME FrameCounter \n" 
+       << " -singlechan : Force one channel per fiber \n"
        << " -varcut min,max : min-max cut on variance \n" 
+       << " -nband nband,first,last: numbers of freq. bands and first and last bands used for cuts \n"
+       << " -filldt : force data table filling \n"
+       << " -freqfilter : force median filtering on the frequencies \n"
        << " -gain filename : spectral response fits file name \n" 
        << " -tspwin wsz,extsz,nfiles : spectra time (paquet no) window (ex: -tspwin 120,4,5) \n"
@@ -185 +212 @@
   char flnm[1024];
   sprintf(flnm,"%s/signal%d.fits",dirlist_[0].c_str(),imin_);
-  return DecodeMiniFitsHeader(flnm,paqsize_, npaqinfile_,fgdatafft_, fgsinglechannel_); 
+  int rc = DecodeMiniFitsHeader(flnm,paqsize_, npaqinfile_,fgdatafft_, fgsinglechannel_); 
+  if(fgforcesinglechan_) fgsinglechannel_= true;
+  return rc;
 }
 
@@ -199 +228 @@
   cout << " Action=" << action_ << "  NMean=" << nmean_ << " NBloc=" << nbloc_ << endl;
   cout << " FreqMin= " << freqmin_ << " FreqMax= " << freqmax_ << " NBinFreq= " << nbinfreq_ << endl;
-  cout << " GainFileName=" << gainfile_ << " VarianceMin= " << vmin_ << " VarianceMax= " << vmax_ << endl;
+  cout<<  ((fgdt_)?" Fill DadaTable ":" NO DataTable") << endl;
+  cout << " GainFileName=" << gainfile_ << " Variance: Min= " << vmin_ << " Max= " << vmax_ 
+       << " Bands: N=" <<  nbands_ << " First=" << bandfirst_ << " Last=" << bandlast_
+       << endl;
   cout << " Spectra TimeWindow (Nb.Paquets) Size=" << spec_win_sz_ << " ExtensionSize=" << spw_ext_sz_ 
        << " MaxNbFile=" << nbmax_specwfiles_ << endl; 

trunk/AddOn/TAcq/branap.h

@@ -41 +41 @@
   bool fgdatafft_;  // true -> donnee FFT (provenant du firmware FFT)
   bool fgsinglechannel_; // true -> un seul canal par fibre (par defaut=2 canaux/fibres)
+  bool fgforcesinglechan_; // true force par la ligne de commande
   int prtlevel_;
   long int prtmodulo_; 
@@ -51 +52 @@
 
   double vmin_,vmax_;   // coupure min-max sur la variance 
+  int nbands_,bandfirst_,bandlast_; //bandes en frequences pour mean, variance et cut
+  bool fgdt_; //if true fill datatable
+
+  bool fgfreqfilter_; // if true force median filtering on the frequencies
+
   string gainfile_;     // nom du fichier de gain (reponse spectrale)
 };

trunk/AddOn/TAcq/brproc.cc

@@ -43 +43 @@
     fgdatafft_(fgdatafft), fgsinglechannel_(fgsinglechan), 
     nbpaq4mean_(fgsinglechan?memgr_.NbFibres():2*memgr_.NbFibres()), 
-    nbadpaq_(fgsinglechan?memgr_.NbFibres():2*memgr_.NbFibres()) 
+    nbadpaq_(fgsinglechan?memgr_.NbFibres():2*memgr_.NbFibres())
 {
   setNameId("meanSpecCalc",1);
-  BRPaquet paq(memgr_.PaqSize());
+
+  uint_4 nb_octets_entrop = 0; //this value is valid for Dec. 2010 data at Nancay
+  const char* venvp = NULL; 
+  venvp=getenv("BRANA_NBYTECUT");
+  if (venvp!=NULL){
+    nb_octets_entrop = atoi(venvp);
+    cout << "BRMeanSpecCalculator : BRANA_NBYTECUT : " << nb_octets_entrop << endl;
+  }
+
+  BRPaquet paq(memgr_.PaqSize()-nb_octets_entrop);
+
   if (fgsinglechannel_) {
     mspecmtx_.SetSize(memgr_.NbFibres(), paq.DataSize()/2); 
@@ -70 +80 @@
   }
 
+  //
+  
+
   // Definition des tailles de fenetres de spectres, etc ...
   SetSpectraWindowSize();
-  SetMaxNbSepcWinFiles();
+  SetMaxNbSpecWinFiles();
   nbtot_specwin_=0;
   SetVarianceLimits();
+  SetNumberOfBands();
 
   ofsdtp_=NULL;
   dtp_=NULL;
-  xnt_=new double[nchan*2];  // CHECK : ATTENTION la taille depend de nombre de colonne du NTuple !
+  
+  //  cout << "(JEC) creation tuple de " << nchan*(1+numberOfBands_) << " doubles " << endl;
+  //  xnt_=new double[nchan*(1+numberOfBands_)];  // CHECK : ATTENTION la taille depend de nombre de colonne du NTuple !
+  xnt_=NULL;
 }
 
@@ -84 +101 @@
 BRMeanSpecCalculator::~BRMeanSpecCalculator()
 {
+  cout << " ---------------- BRMeanSpecCalculator()_Finalizing -------------------- " << endl;
   uint_8 npqm=0;
   for(size_t i=0; i<nbpaq4mean_.size(); i++)  npqm+=nbpaq4mean_[i];
   if (npqm>nmean_*nbpaq4mean_.size()/10)  SaveMeanSpectra();
-  cout << " ---------------- BRMeanSpecCalculator()_Finalizing -------------------- " << endl;
   for(size_t i=0; i<nbadpaq_.size(); i++) {
     cout << " Channel " << i << "  NBadPaq=" << nbadpaq_[i] << " / TotNbPaq=" << totnbpaq_ << endl;
@@ -127 +144 @@
 void BRMeanSpecCalculator::DefineDataTable()
 {
+  cout << "(JEC) BRMeanSpecCalculator::DefineDataTable START" << endl;
   string dtfile="!"+outpath_+"/dtspec.fits";
   ofsdtp_ = new FitsInOutFile(dtfile,FitsInOutFile::Fits_Create);
@@ -136 +154 @@
     dtp_->AddFloatColumn(cnom);    
   }
+  for(int i=0; i<nchan; i++) {
+    for(int j=0;j<numberOfBands_;j++){
+      sprintf(cnom,"varnormb%d%d",i,j);
+      dtp_->AddFloatColumn(cnom);    
+    }
+  }
   /*
   for(int i=0; i<nchan; i++) {
@@ -143 +167 @@
   */
   //  xnt_=new double[nchan*2];  CHECK : faut-il reallouer ?
+  cout << "(JEC) creation tuple de " << nchan*(1+numberOfBands_) << " doubles " << endl;
+  xnt_=new double[nchan*(1+numberOfBands_)];
+  cout << "(JEC) BRMeanSpecCalculator::DefineDataTable END" << endl;
+}
+  
+/* --Methode-- */
+void BRMeanSpecCalculator::SetNumberOfBands(int numberOfBands, int ibandfirst, int ibandlast)
+{
+  if (numberOfBands < 1) numberOfBands = 1;
+  numberOfBands_ = numberOfBands;
+  if (ibandfirst < 0 )ibandfirst = 0;
+  if (ibandlast >= numberOfBands_ ) ibandlast = numberOfBands_-1;
+  ibandfirst_=ibandfirst; ibandlast_=ibandlast;
+
+  cout << "(JEC): SetNumberOfBands (END) : "
+       << numberOfBands_ << " "
+       << ibandfirst_ << " "
+       << ibandlast_ << endl;
+
 }
 
@@ -154 +197 @@
   cout << " MeanGain=" << sgain_.Sum()/sgain_.Size() << " ApplyGains=" 
        << ((fg_apply_gains_)?"true":"false") << endl;
-}
-
-static inline r_4 Zmod2(complex<r_4> z) 
-{ return (z.real()*z.real()+z.imag()*z.imag()); }
+  if( (spec_window_.SizeX()!= sgain_.NCols()) || (spec_window_.SizeY()!= sgain_.NRows()) ){
+    cout << " ReadGainFitsFile: BAD Gain Matrix sizes " << endl;
+    sgain_.Show();
+    spec_window_.Show();
+    throw ParmError("ReadGainFitsFile: BAD Gain Matrix sizes");
+  }
+}
+
+//JEC
+//static inline r_4 Zmod2(complex<r_4> z) 
+//{ return (z.real()*z.real()+z.imag()*z.imag()); }
 
 
@@ -169 +219 @@
   // des spectres le cas echeant ProcSpecWin()
 
+ 
   int_8 nbpaqdec = (int_8)totnbpaq_-(int_8)GetSpecWinOverlapSize();
   if ((nbpaqdec>0)&&(nbpaqdec%GetSpectraWindowSize()==0)) {
@@ -208 +259 @@
   if (fg_apply_gains_) {   // Application des gains, si demande 
     sa_size_t kz=PaqNumToArrayIndex(totnbpaq_);
-    for(sa_size_t i=0; i<spec_window_.SizeY(); i++) 
-      (spec_window_(Range::all(), Range(i), Range(kz))).Div(sgain_.Row(i));
+    for(sa_size_t i=0; i<spec_window_.SizeY(); i++) {
+      (spec_window_(Range::all(), Range(i), Range(kz)).CompactAllDimensions()).Div(sgain_.Row(i).CompactAllDimensions());
+    }
   }
 
@@ -231 +283 @@
     }
   }
-  nbtot_specwin_++;
-  return;
   // On appelle la routine de nettoyage qui doit flagger les mauvais paquets 
   FlagBadPackets(numpaqstart, numpaqend);
@@ -249 +299 @@
       // Calcul de spectres moyennes et variance 
       for(sa_size_t f=1; f<spec.Size(); f++)  {   // boucle sur les frequences 
-        spec(f) += spec_window_(f,i,kz);
+        spec(f)  += spec_window_(f,i,kz);
         sspec(f) += spec_window_(f,i,kz)*spec_window_(f,i,kz);  
       }
@@ -270 +320 @@
     for(sa_size_t i=0; i<spec_window_.SizeY(); i++)  {
       double mean, sigma;
-      sa_size_t kz=PaqNumToArrayIndex(totnbpaq_); 
+      ////////BUG      sa_size_t kz=PaqNumToArrayIndex(totnbpaq_); 
       double variance=0.;
       variance=spec_window_(Range(1,Range::lastIndex()), Range(i), Range(kz)).Sum();
-      //      xnt_[i]=variance;
+      xnt_[i]=variance;
+      //Compute nomalized variance in bands freq. 
+      sa_size_t fMin;
+      sa_size_t fMax;
+      int bandW = spec_window_.SizeX()/numberOfBands_;
+      vector<double> varNomBinned(numberOfBands_);
+      for (sa_size_t j=0; j<numberOfBands_; j++){
+        fMin = j*bandW;
+        fMax =fMin+bandW-1;
+        varNomBinned[j]=spec_window_(Range(fMin,fMax), Range(i), Range(kz)).Sum();
+        varNomBinned[j]/=(r_4)bandW;
+        xnt_[spec_window_.SizeY()+i*numberOfBands_+j] =  varNomBinned[j];
+      }//eof
+
       clnflg_(i,kz)=0;
-      if (variance<varmin_) { clnflg_(i,kz)=1;  nbadpaq_[i]++; }
-      else if (variance>varmax_) { clnflg_(i,kz)=2;  nbadpaq_[i]++; }
+      for (sa_size_t j=ibandfirst_; j<=ibandlast_; j++){
+        //        cout << "(jec) var["<<j<<"] =" << varNomBinned[j]
+        //     << " min " << varmin_
+        //     << " max " << varmax_ << endl;
+        if(varNomBinned[j]<varmin_)      { clnflg_(i,kz)=10+j;  nbadpaq_[i]++; break;}
+        else if(varNomBinned[j]>varmax_) { clnflg_(i,kz)=100+j; nbadpaq_[i]++; break;}
+      }
+      //cout << "clnflg_("<<i<<","<<kz<<"): " << clnflg_(i,kz)  << endl;
     }
     if (dtp_)  dtp_->AddRow(xnt_);
@@ -343 +412 @@
   : BRBaseProcessor(memgr), fgsinglechannel_(fgsinglechannel), totnbfftpaq_(0)
 {
-  BRPaquet paq(memgr_.PaqSize());
+  uint_4 nb_octets_entrop = 0; //this value is valid for Dec. 2010 data at Nancay
+  const char* venvp = NULL; 
+  venvp=getenv("BRANA_NBYTECUT");
+  if (venvp!=NULL){
+    nb_octets_entrop = atoi(venvp);
+    cout << "BRFFTCalculator  : BRANA_NBYTECUT : " << nb_octets_entrop << endl;
+  }
+
+  BRPaquet paq(memgr_.PaqSize()-nb_octets_entrop);
+  //JEC END
+  //  BRPaquet paq(memgr_.PaqSize());
   setNameId("FFTCalc",2);
   ffts_.SetInDataSize((fgsinglechannel_)?paq.DataSize():paq.DataSize()/2);
@@ -633 +712 @@
         }
         for(sa_size_t j=0; j<spectreV1.Size(); j++) 
-          spectreV1(j) += Zmod2(cfour1(j));
+          spectreV1(j) += BRMeanSpecCalculator::Zmod2(cfour1(j));
         memcpy(procbuff+i*procpaqsz, cfour1.Data(), sizeof(complex<r_4>)*cfour1.Size());
         if (fgtimfreq) {   // Remplissage tableau temps-frequence
           for(sa_size_t c=1; c<timfreqV1.NCols(); c++) {
             for(sa_size_t j=c*nfsmap_; j<(c+1)*nfsmap_; j++) 
-              timfreqV1(nzm, c) += Zmod2(cfour1(j));
+              timfreqV1(nzm, c) += BRMeanSpecCalculator::Zmod2(cfour1(j));
           }
         }
         for(sa_size_t j=0; j<spectreV2.Size(); j++) 
-          spectreV2(j) += Zmod2(cfour2(j));  // Zmod2(zp2[j]); 
+          spectreV2(j) += BRMeanSpecCalculator::Zmod2(cfour2(j));  // BRMeanSpecCalculator::Zmod2(zp2[j]); 
         memcpy(procbuff+i*procpaqsz+procpaqsz/2, cfour2.Data(), sizeof(complex<r_4>)*cfour2.Size());
         if (fgtimfreq) {   // Remplissage tableau temps-frequence
           for(sa_size_t c=1; c<timfreqV2.NCols(); c++) {
             for(sa_size_t j=c*nfsmap_; j<(c+1)*nfsmap_; j++) 
-              timfreqV2(nzm,c) += Zmod2(cfour2(j));
+              timfreqV2(nzm,c) += BRMeanSpecCalculator::Zmod2(cfour2(j));
           }
         }
@@ -667 +746 @@
         if (nmnt==0)  { xnt[0]=paq.FrameCounter();  xnt[1]=paq.TimeTag(); }
         for(sa_size_t j=2700; j<2800; j++) { 
-          ms1 += Zmod2(cfour1(j)); ms2 += Zmod2(cfour2(j)); 
+          ms1 += BRMeanSpecCalculator::Zmod2(cfour1(j)); ms2 += BRMeanSpecCalculator::Zmod2(cfour2(j)); 
           complex<r_4> zvis =  cfour1(j)*conj(cfour2(j));
-          ms12 += Zmod2(zvis);   ms12re += zvis.real();  ms12im += zvis.imag();
+          ms12 += BRMeanSpecCalculator::Zmod2(zvis);   ms12re += zvis.real();  ms12im += zvis.imag();
           ms12phi+= atan2(zvis.imag(),zvis.real());
         } 
@@ -742 +821 @@
           int_4 sspvmaxidx[3] = {-1,-1,-1};
           for(int jji=1;jji<visiV12.Size()-1;jji++) {
-            r_4 zmv2 = Zmod2(visiV12(jji)); 
+            r_4 zmv2 = BRMeanSpecCalculator::Zmod2(visiV12(jji)); 
             if (zmv2>sspvmax[2]) { sspvmax[2]=zmv2; sspvmaxidx[2]=jji; }
           }

trunk/AddOn/TAcq/brproc.h

@@ -37 +37 @@
 //---------------------------------------------------------------------
 
-class BRMeanSpecCalculator : public BRBaseProcessor
+class BRMeanSpecCalculator : public BRBaseProcessor 
 {
 public:
@@ -43 +43 @@
                        bool fgdatafft=true, bool fgsinglechan=false);
   virtual ~BRMeanSpecCalculator(); // Fait le nettoyage final , sauvegarde des objets ...
+
+ //JEC 19/1/11 
+  static inline r_4 Zmod2(complex<r_4> z){ return (z.real()*z.real()+z.imag()*z.imag()); }
 
 // Definit la taille de la fenetre temporelle (winsz = nb de paquets) des spectres a traiter
@@ -55 +58 @@
 // Retourne l'index (selon Z) dans le tableau des spectres en fonction du numero de paquet
   inline sa_size_t PaqNumToArrayIndex(uint_8 numpaq)  { return numpaq%spec_window_.SizeZ() ; }
+  inline int GetNumberOfBands() { return numberOfBands_; }
 
 // Pour definir le nombre maximum de fichiers de fenetres de spectres a sauver sur fichiers
-  inline void SetMaxNbSepcWinFiles(uint_4 nmax=0)  {  nmaxfiles_specw_=nmax; }
+  inline void SetMaxNbSpecWinFiles(uint_4 nmax=0)  {  nmaxfiles_specw_=nmax; }
 
 // Pour definir le datatable a remplir - pas de DataTable rempli si pas appele 
@@ -64 +68 @@
   inline void SetVarianceLimits(double vmin=0., double vmax=9.e99)
     { varmin_=vmin;  varmax_=vmax; }
+  
+  //POur definir le nbre de bandes de freequences pour calcul mean, variance, cut
+  //ibandlast incluse
+  void SetNumberOfBands(int numberOfBands = 1, int ibandfirst = 0, int ibandlast = 99999);
+
 
 // lecture des gains en fonction de la frequence depuis un fichier FITS
@@ -77 +86 @@
   virtual void SaveSpectraWindow();  // Pour sauver les spectres de la fenetre temporel ds un fichier 
 
-
+   
   uint_4 nmean_;  // Nombre de spectres pour le calcul des moyennes 
   string outpath_;  // directory pour fichiers de sortie
@@ -97 +106 @@
   uint_4 nbtot_specwin_;    // nombre total de fenetre de spectres traites 
   uint_4 nmaxfiles_specw_;   // nombre maximum de fichier de fenetre_spectre 
+
   
   TMatrix< r_4 > sgain_;   // gain en fonction de la frequence  
@@ -102 +112 @@
 
   double varmin_, varmax_;  // limites pour la coupure sur la variance
+  int ibandfirst_, ibandlast_; //limites pour les numero de bandes a utiliser pour les cut en varaince  
+  int numberOfBands_; // nombre de bandes de frequences pour calcul mean, variance et cut
+
   vector<uint_8> nbadpaq_;
 

trunk/AddOn/TAcq/makefile

@@ -5 +5 @@
 
 PKGHLIST = bracqvers.h branap.h brparam.h racqumem.h brpaqu.h brtypes.h \
-           brdiskw.h brfitsrd.h minifits.h brviscalc.h brbaseproc.h brproc.h \
+           brdiskw.h brfitsrd.h minifits.h brviscalc.h brbaseproc.h brproc.h brprocGain.h\
            mfacq.h  tmtacq.h racqueth.h  racquproc.h racqurw.h pciewrap.h  \
            swrapsock.h burawstream.h sockrawstream.h tmtacq.h
@@ -11 +11 @@
 PKGOLIST = $(OBJ)branap.o  $(OBJ)brparam.o $(OBJ)racqumem.o $(OBJ)brpaqu.o \
            $(OBJ)brfitsrd.o $(OBJ)brdiskw.o  $(OBJ)minifits.o \
-           $(OBJ)brbaseproc.o $(OBJ)brproc.o $(OBJ)brviscalc.o \
+           $(OBJ)brbaseproc.o $(OBJ)brproc.o $(OBJ)brprocGain.o $(OBJ)brviscalc.o \
            $(OBJ)racqueth.o $(OBJ)racquproc.o $(OBJ)racqurw.o  $(OBJ)pciewrap.o \
            $(OBJ)swrapsock.o $(OBJ)sockrawstream.o $(OBJ)burawstream.o
@@ -115 +115 @@
         $(CXXCOMPILE) -c  -o $(OBJ)brviscalc.o brviscalc.cc 
 
+$(OBJ)brprocGain.o : brprocGain.cc $(PKGHLIST)
+        $(CXXCOMPILE) -c  -o $(OBJ)brprocGain.o brprocGain.cc 
+
 $(OBJ)brproc.o : brproc.cc $(PKGHLIST)
         $(CXXCOMPILE) -c  -o $(OBJ)brproc.o brproc.cc 

trunk/AddOn/TAcq/specmfib.cc

@@ -26 +26 @@
 #include "brfitsrd.h"
 #include "brproc.h"
+
+//JEC 19/1/11 START
+#include "brprocGain.h"
+//JEC END
+
 #include "brdiskw.h"
 
@@ -50 +55 @@
     if (rc) return rc;
     par.Print(cout);
-    if ((par.action_!="cube3d")&&(par.action_!="mspec")&&(par.action_!="bproc")) {
+    if ((par.action_!="cube3d")
+        &&(par.action_!="mspec")
+        &&(par.action_!="gain") //JEC 19/1/11 add gain action
+        &&(par.action_!="bproc") 
+        ) {
       cout << " !!! specmfib.cc  BAD action = " << par.action_ << " possible values: mspec,cube3d,bproc" << endl;
       return 5;
@@ -73 +82 @@
     BRMeanSpecCalculator procms(mmgr, par.outpath_, par.nmean_, par.fgdatafft_, par.fgsinglechannel_);
     procms.SetVarianceLimits(par.vmin_, par.vmax_);
-    if (par.gainfile_.length()>0) procms.ReadGainFitsFile(par.gainfile_);
+    if (par.nbands_>0) procms.SetNumberOfBands(par.nbands_,par.bandfirst_,par.bandlast_);
+    //JEC 27/1/11 see below    if (par.gainfile_.length()>0) procms.ReadGainFitsFile(par.gainfile_);
     procms.SetPrintLevel(par.prtlevel_,par.prtmodulo_);
-    // Il faut appeler procms.DefineDataTable() si on veut remplir le DataTable
+    if (par.fgdt_) procms.DefineDataTable();
+
+    BRGainCalculator* procgain_p=NULL;
+    if (par.action_=="gain") {
+      procgain_p = 
+        new BRGainCalculator(mmgr, par.outpath_, par.nmean_, par.fgdatafft_, par.fgsinglechannel_);
+      procgain_p->SetPrintLevel(par.prtlevel_,par.prtmodulo_);
+    }
+
 
     BRFFTCalculator procfft(mmgr, par.fgsinglechannel_);
-    if (!par.fgdatafft_) procms.SetMemZAction(MemZA_ProcB);
+    if (!par.fgdatafft_) { 
+      procms.SetMemZAction(MemZA_ProcB);
+      if(procgain_p) procgain_p->SetMemZAction(MemZA_ProcB);
+    }
+
     if (par.spec_win_sz_>0) { 
       procms.SetSpectraWindowSize(par.spec_win_sz_, par.spw_ext_sz_);
-      procms.SetMaxNbSepcWinFiles(par.nbmax_specwfiles_);
+      procms.SetMaxNbSpecWinFiles(par.nbmax_specwfiles_);
+      if(procgain_p){
+        procgain_p->SetSpectraWindowSize(par.spec_win_sz_, par.spw_ext_sz_);
+        procgain_p->SetMaxNbSpecWinFiles(par.nbmax_specwfiles_);
+      }
     }
+    //JEC 27/1/11 should be done after SetSpectraWindowSize to perform x-checking
+    if (par.gainfile_.length()>0) procms.ReadGainFitsFile(par.gainfile_);
+
     FitsCubeWriter wrt(mmgr, par.outpath_, par.nbloc_);
 
@@ -92 +121 @@
     }
 
-    cout << " specmfib: Starting threads (reader meanSpecCalculator ...) ... " << endl;
+    cout << " specmfib: Starting threads (reader + "<< par.action_ << ")" << endl;
     reader.start();
     if (par.action_ == "cube3d")   wrt.start();
     else {  // Calcul spectre moyenne 
       if (!par.fgdatafft_)  procfft.start(); 
+      //JEC 19/1/11 START add the gain process alternative
       if (par.action_=="bproc")  bproc_p->start();
+      else if (par.action_=="gain") procgain_p->start();
       else procms.start();
+      //JEC END
     }
     usleep(200000);
@@ -105 +137 @@
     else {
       if (!par.fgdatafft_)  procfft.join(); 
+
+      //JEC 19/1/11 START add the gain process alternative
       if (par.action_=="bproc")  bproc_p->join();
+      else if (par.action_=="gain") procgain_p->join();
       else procms.join();
+      //JEC END
     }
     mmgr.Print(cout);
     if (bproc_p) delete bproc_p;
+    if (procgain_p) delete procgain_p;
     cout << resu ; 
   }