source: Sophya/trunk/AddOn/TAcq/brviscalc.cc@ 4016

//----------------------------------------------------------------
// Projet BAORadio - (C) LAL/IRFU  2008-2010
// Classes de threads de calcul de visibilites pour BAORadio 
//----------------------------------------------------------------

#include <stdlib.h>
#include <string.h>
#include <unistd.h>
#include <fstream>
#include <signal.h>

#include "pexceptions.h"
#include "fioarr.h"
#include "fitsarrhand.h"
#include "matharr.h"
#include "arrctcast.h"
#include "timestamp.h"
#include "ctimer.h"

#include "brviscalc.h"


//---------------------------------------------------------------------
// Classe de traitement - calcul de visibilite pour n fibres 
//---------------------------------------------------------------------
/*!
   \class BRVisibilityCalculator
   \ingroup TAcq

   \brief Visibility computation class.
*/

/* --Methode-- */
BRVisibilityCalculator::BRVisibilityCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean, size_t nthr)
  : BRBaseProcessor(memgr), paralex_(*this, nthr), nparthr_(nthr), 
    outpath_(outpath), nbcalc_(1), calcid_(0), 
    vpdata_(2*memgr.NbFibres()), vpdatar_(2*memgr.NbFibres())
    // , dtfos_(outpath+"visdt.fits", Fits_Create), visdt_(dtfos_, 1024, true);
{  
  SetFFTData();
  SetNPaqIntervalMode(nmean);
  DefineRank(1,0);
  SetPPFOutput();

  npaqcumul_=0;
  threshold_tt_=0;
  moyfc_=moytt_=0.;
  first_fc_=first_tt_=0;

  uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
  chanids_.SetSize(2*memgr_.NbFibres());
  chanpairnumall_.SetSize(maxnpairs);
  chanpairsall_.SetSize(maxnpairs,2);
  for(size_t i=0; i<2*memgr_.NbFibres(); i++)   { 
    vpdata_[i]=NULL;   vpdatar_[i]=NULL;  
  }
  // Pas necessaire - appele par DefineRank-  SelectPairs();

  //  visdt_.AddFloatColumn("mfc");
  visdt_.AddDoubleColumn("mfc");
  visdt_.AddDoubleColumn("mtt");
  visdt_.AddIntegerColumn("jfreq");
  visdt_.AddIntegerColumn("numch");
  visdt_.AddFloatColumn("vre");
  visdt_.AddFloatColumn("vim");
  ActivateVisDTable(); 

  if (nmean_ < 1)  nmean_=memgr_.NbPaquets();
  if (nmean_ < 1)  nmean_=1;

  totnbpaq_=0;
  numfile_=0;
  nb_flop_=0.;

  fgallfibok=NULL;
  fgcktt_=false;
  setNameId("viscalc", 0); 
}

/* --Methode-- */
BRVisibilityCalculator::~BRVisibilityCalculator()
{
  if (totnbpaq_<1)   return; 
  if (fgvisdt_) {
    cout << " ~BRVisibilityCalculator - Visibility Datatable : " << endl;
    cout << visdt_;
  
    string filename;
    filename = outpath_+"visdt.ppf";
    if (nbcalc_>1) {
      char sbuff[32];
      sprintf(sbuff,"visdt_%d.ppf",(int)calcid_);
      filename = outpath_+sbuff;
    }
    POutPersist po(filename);
    po << visdt_;
  }
  if (calcid_ == 0) {    
    if (fgcktt_) {
      cout << " BRVisibilityCalculator -  Check TimeTag done: TotNPaqProc= " << totnbpaq_ << endl;
      for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  {
        cout << " BRTTCheck-Fiber[" << fib << "] NBadTT=" << vbadtt_[fib] << " NDiffTT>5="
             << vndiff5tt_[fib] << "  NotSameTT=" << vnsamett_[fib] << endl; 
      }
      POutPersist pott(outpath_+"ttfcmtx.ppf");
      pott << PPFNameTag("FC") << fcmtx_; 
      pott << PPFNameTag("TT") << ttmtx_;     
    }
  }
}

/* --Methode-- */
void BRVisibilityCalculator::DefineRank(uint_4 nbc, uint_4 cid, uint_4 pair1, uint_4 nbpairs, bool fgpimp)
{
  if ((nbc>6)||(cid>=nbc)) 
    throw ParmError("BRVisibilityCalculator::DefineRank() NbCalc > 6 !");
  uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
  if ((pair1>=maxnpairs)||(nbpairs<1))
    throw ParmError("BRVisibilityCalculator::DefineRank() (pair1>=maxnpairs)||(nbpairs<1) !");
  maxnpairs-=pair1;
  if (maxnpairs>nbpairs) maxnpairs=nbpairs;

  nbcalc_=nbc;
  calcid_=cid;
  if (nbcalc_>1) {
    uint_4 npairspth=maxnpairs/nbcalc_;  // nb de paires par thread
    if (calcid_==(nbcalc_-1)) 
      SelectPairs(calcid_*npairspth+pair1, maxnpairs-calcid_*npairspth, fgpimp);
    else 
      SelectPairs(calcid_*npairspth+pair1, npairspth, fgpimp);
    MemZaction mmzas[10]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF,
                          MemZA_ProcG,MemZA_ProcH,MemZA_ProcI,MemZA_ProcJ};
    SetMemZAction(mmzas[calcid_]);
    setNameId("viscalc_grp", calcid_); 
  }
  else SelectPairs(pair1,maxnpairs, fgpimp);
  return ;
}

/* --Methode-- */
uint_4 BRVisibilityCalculator::SelectPairs(uint_4 pair1, uint_4 nbpairs, bool fgpimp)
{
  BRPaquet paq(memgr_.PaqSize());
  uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();

  if (pair1 >= maxnpairs)  pair1=maxnpairs-1;
  if (nbpairs > maxnpairs-pair1)  nbpairs=maxnpairs-pair1;

  pairst_=pair1;
  nbpairs_=nbpairs;
  fgpimp_=fgpimp;

  uint_4 nbvisicomp=0;
  sa_size_t kpair=0;  
  for(size_t i=0; i<vpdata_.size(); i++) {
    for(size_t j=i; j<vpdata_.size(); j++) {
      kpair++;
      if (kpair<(pairst_+1))  continue;
      if (kpair>=(pairst_+nbpairs_+1))  break;
      if (fgpimp_&&(i!=j)&&((i+j)%2==0))  continue;  // calcul des visib avec numero pair-impair + autocorrel
      nbvisicomp++;
    }
  }

  vismtx_.SetSize(nbvisicomp, paq.DataSize()/4); 

  chanpairnum_.SetSize(nbvisicomp);
  chanpairs_.SetSize(nbvisicomp,2);

  return nbvisicomp;
}

/* --Methode-- */
int BRVisibilityCalculator::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
{
  if (freq2<freq1) {
    freq1=0; freq2=vismtx_.NCols()-1; nbfreq=1; 
  }
  jf1_=freq1; jf2_=freq2;
  if ((jf1_<1)||(jf1_>=vismtx_.NCols()))  jf1_=1;
  if ((jf2_<1)||(jf2_>=vismtx_.NCols())||(jf2_<jf1_))  jf2_=vismtx_.NCols()-1;
  if (nbfreq<1) nbfreq=(jf2_-jf1_);
  djf_=(jf2_-jf1_)/nbfreq;
  if (djf_<1) djf_=1;
  cout << " BRVisibilityCalculator::SelectFreqBinning/Info  JF1=" << jf1_ 
       << " JF2=" << jf2_ << " DJF=" << djf_ << endl;

}

/* --Methode-- */
void BRVisibilityCalculator::UpdateChanIds()
{

  for(size_t i=0; i<memgr_.NbFibres(); i++) { 
    chanids_(2*i)=memgr_.FiberId(i)*2-1; 
    chanids_(2*i+1)=memgr_.FiberId(i)*2; 
  } 
  sa_size_t k=0;     // numero de ligne dans la matrice des visibilites 
  for(size_t i=0; i<vpdata_.size(); i++) {
    for(size_t j=i; j<vpdata_.size(); j++) {
      chanpairnumall_(k)=chanids_(i)*CHANPAIRCONVFAC+chanids_(j);
      chanpairsall_(k,0)=chanids_(i);  chanpairsall_(k,1)=chanids_(j);     k++;
    }
  }
  sa_size_t kpair=0;  
  k=0;  // numero de ligne dans la matrice des visibilites 
  for(size_t i=0; i<vpdata_.size(); i++) {
    for(size_t j=i; j<vpdata_.size(); j++) {
      kpair++;
      if (kpair<(pairst_+1))  continue;
      if (kpair>=(pairst_+nbpairs_+1))  break;
      if (fgpimp_&&(i!=j)&&((i+j)%2==0))  continue;  // calcul des visib avec numero pair-impair + autocorrel
      chanpairnum_(k)=chanids_(i)*CHANPAIRCONVFAC+chanids_(j);
      chanpairs_(k,0)=chanids_(i);  chanpairs_(k,1)=chanids_(j);     k++;
    }
  }

  string filename;
  filename = (outpath_+"chanum.")+OutFileExtension();
  if (nbcalc_>1) {
    char sbuff[32];
    sprintf(sbuff,"chanum_%d.%s",(int)calcid_,OutFileExtension());
    filename = outpath_+sbuff;
  }
  if (fgfitsout_) {  // Ecriture au format FITS
    FitsInOutFile foc(filename, FitsInOutFile::Fits_Create);
    foc.SetNextExtensionName("chanids"); 
    foc << chanids_; 
    foc.SetNextExtensionName("chanpairs");
    foc << chanpairs_; 
    foc.SetNextExtensionName("chanpairnum");
    foc << chanpairnum_; 
    foc.SetNextExtensionName("chanpairsall");
    foc << chanpairsall_; 
    foc.SetNextExtensionName("chanpairnumall");
    foc << chanpairnumall_; 
  }
  else {  // Format PPF 
    POutPersist poc(filename);
    poc << PPFNameTag("chanids") << chanids_; 
    poc << PPFNameTag("chanpairs") << chanpairs_; 
    poc << PPFNameTag("chanpairnum") << chanpairnum_; 
    poc << PPFNameTag("chanpairsall") << chanpairsall_; 
    poc << PPFNameTag("chanpairnumall") << chanpairnumall_; 
  }
  cout << "BRVisibilityCalculator[" << calcid_ << "]::UpdateChanIds() Channel Ids/Pairs saved to file " 
       << filename << endl; 
  cout << " ... NbVisib=NbChanPairs=" << chanpairs_.NRows() << " ChannelPairs= " ;
  for(sa_size_t ir=0; ir<chanpairs_.NRows(); ir++) {
    if (ir%10==0)  cout << endl;
    cout << "(" << chanpairs_(ir,0) << "," << chanpairs_(ir,1) << ") ";
  }
  cout << endl;
  return;
}

/* --Methode-- */
int BRVisibilityCalculator::ActivateTimeTagCheck(uint_8 maxnpaq)
{
  mindeltatt_=memgr_.PaqSize()/2;
  if (mindeltatt_<1) mindeltatt_=1;
  fcmtx_.SetSize(memgr_.NbFibres(), maxnpaq);
  ttmtx_.SetSize(memgr_.NbFibres(), maxnpaq);
  vlasttt_.resize(memgr_.NbFibres(), 0);
  vbadtt_.resize(memgr_.NbFibres(), 0);
  vnsamett_.resize(memgr_.NbFibres(), 0);
  vndiff5tt_.resize(memgr_.NbFibres(), 0);

  fgcktt_=true;
  cout << " BRVisibilityCalculator::ActivateTimeTagCheck() - TT/Fc matrix NCols=" << maxnpaq 
       << " MinDeltaTT=" << mindeltatt_ << endl;

  return 0;
}


/* --Methode-- */
void BRVisibilityCalculator::run()
{
  cout << " BRVisibilityCalculator[" << calcid_ << "/" << nbcalc_ 
       << "]::run() - Starting " << " NFibers=" << memgr_.NbFibres() 
       << " NChan=" << 2*memgr_.NbFibres() << " NPairs=" << nbpairs_ << " First:" << pairst_ << endl;   

  if (fgnpaq_time) 
    cout << " BRVisibilityCalculator::run()/Info  NPaq=" << nmean_ << " Interval mode for visibility mtrix save" << endl;
  else 
    cout << " BRVisibilityCalculator::run()/Info  Time Tnterval (dt=" << dtimemean_ << ") for visibility mtrix save" << endl;

  if (nparthr_ < 2)  return BRBaseProcessor::run();
  // Execution multithread parallele
  setRC(1);     
  int rc=0;
  try {
    if ((nmean_%memgr_.NbPaquets())!=0) {
      uint_4 mnmean = (nmean_/memgr_.NbPaquets()+1)*memgr_.NbPaquets();
      cout << " BRVisibilityCalculator::run()/Info changing nmean=" << nmean_  << " to multiple of"
           << " memgr_.NbPaquets() -> " << mnmean << endl;
      nmean_=mnmean;
    }
    paralex_.SetParallelTask(*this);
    cout << " BRVisibilityCalculator::run()/Info : starting ParallelExecutor with nThreads=" 
         << paralex_.nThreads() << " ... " << endl;
    paralex_.start();

    fgallfibok = new bool[memgr_.NbPaquets()];

    size_t paqsz=memgr_.PaqSize();
    size_t procpaqsz=memgr_.ProcPaqSize();
    bool fgrun=true;
    while (fgrun) {
      if (stop_) break;
      if (memgr_.GetRunState() == MemZR_Stopped) break;
      int mid = memgr_.FindMemZoneId(mmact_);  // (MemZA_ProcA);
      //      Byte* buffg = memgr_.GetMemZone(mid);
      //      if (buffg == NULL) {
      if (mid < 0) {
        cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR FindMemZoneId(" 
             << (int)mmact_ << ") ->" << mid << ") -> NULL" << endl;
        setRC(7);      fgrun=false;             
        break;  
      }
      cts_=memgr_.GetAuxData(mid)->FillTime();   // get associated date/time (DATEOBS)
      if (totnbpaq_==0) UpdateChanIds();  // Appele ici pour etre sur que le thread de remplissage a mis l'info a jour.

      for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
        fbuff_[fib] = memgr_.GetMemZone(mid,fib);
        if (fbuff_[fib] == NULL) { // cela ne devrait pas arriver 
          cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetMemZone(" << mid << "," << fib << ") -> NULL" << endl;
          setRC(9);       fgrun=false;          
          break;
        }
        if ((procpaqsz>0)&&((fprocbuff_[fib]=memgr_.GetProcMemZone(mid,fib))==NULL)) {   // cela ne devrait pas arriver non plus
          cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetProcMemZone(" 
               << mid << "," << fib << ") -> NULL" << endl;
          setRC(9);       fgrun=false;          
          break;
        }
      }

      if ( CheckInterval4VisMtxSave() ) {
        if ((prtlev_>0)&&(numfile_%prtmodulo_==0)) {
          cout << numfile_ << "-BRVisCalc[" << calcid_ << "/" << nbcalc_ << "]::run() NPaqProc=" << totnbpaq_ 
               << " TTSec=" << getCurTimeTagSeconds() << " TotMegaFLOP=" << (uint_8)TotNbMegaFLOP() 
               << " VisibMtx_FileNum=" << numfile_ << endl;
        }
        SaveVisibilityMatrix();
      }

      for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) {   // boucle sur les paquets d'une zone  
        fgallfibok[jp]=fgokallfibers_=true;
        for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
          vpaq_[fib].Set(fbuff_[fib]+jp*paqsz);
          vfgok_[fib] = vpchk_[fib].Check(vpaq_[fib],curfc_[fib]);
          if (!vfgok_[fib])  fgallfibok[jp]=fgokallfibers_=false;
          //Pas utile     if (procpaqsz>0)    vprocpaq_[fib] = fprocbuff_[fib]+jp*procpaqsz;
        }
        if (fgokallfibers_)  {
          if (totprocnpaq_==0) {
            for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  {
              fcfirst_[fib]=curfc_[fib];
              ttfirst_[fib]=vpaq_[fib].TimeTag();
            }
            threshold_tt_=vpaq_[0].TimeTag();
            first_fc_=curfc_[0];
            first_tt_=vpaq_[0].TimeTag();
          }
          totprocnpaq_++;
          moyfc_ += curfc_[0];
          moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
          if ((fgcktt_)&&(calcid_==0))  CheckTimeTag();
          if (npaqcumul_==0) {
            first_fc_=curfc_[0];
            first_tt_= vpaq_[0].TimeTag();
          }
          npaqcumul_++;
          totnbpaq_++;
        }
      }  // Fin de boucle sur les paquets 

      // Execution parallele  pour calcul des visibilites par bandes de frequence
      int rcpex=paralex_.execute();
      if (rcpex!=0)  cout << " BRVisibilityCalculator[" << calcid_ << "]::run() / Error Rc[paralex_.execute()]=" << rcpex << endl;

      memgr_.FreeMemZone(mid, mmsta_);  // (MemZS_ProcA);
    } // Fin de boucle sur les zones a traiter 
    //------------------------------------
    cout << " --------- END BRVisibilityCalculator[" << calcid_ << "]::run() , TotNbProcPaq=" << totprocnpaq_ << endl;
    /*
    for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  vpchk_[fib].Print();
    cout << " ------------------------------------ " << endl;
    */
    delete[] fgallfibok;
  }
  catch (std::exception& exc) {
    cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched std::exception " << exc.what() << endl;
    setRC(98);  
    return; 
  }
  catch(...) {
    cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched unknown ... exception " << endl;
    setRC(99);  
    return; 
  }
  
}


/* --Methode-- */
int BRVisibilityCalculator::Process()
{
  if (totnbpaq_==0) {    // qques initialisations a partir du premier paquet 
    UpdateChanIds();  // Appele ici pour etre sur que le thread de remplissage a mis l'info a jour.
    threshold_tt_=vpaq_[0].TimeTag();
    first_fc_=curfc_[0];
    first_tt_= vpaq_[0].TimeTag();
  }
  if (fgdataraw_) {  // Donnees firmware RAW apres TF soft 
    for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
      vpdatar_[2*fib] = reinterpret_cast< complex<r_4>* > (vprocpaq_[fib]);
      vpdatar_[2*fib+1] = reinterpret_cast< complex<r_4>* >(vprocpaq_[fib]+memgr_.ProcPaqSize()/2) ;
    }
  } 
  else {   // donnees firmware FFT 
    for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
      vpdata_[2*fib] = vpaq_[fib].Data1C();
      vpdata_[2*fib+1] = vpaq_[fib].Data2C();
    }
  }
  
  if ( CheckInterval4VisMtxSave() ) {
    if ((prtlev_>0)&&(numfile_%prtmodulo_==0)) {
      cout << numfile_ << "-BRVisCalc[" << calcid_ << "/" << nbcalc_ << "]::Process() NPaqProc="<< totnbpaq_ 
           << " TTSec=" << getCurTimeTagSeconds() << " TotMegaFLOP=" << (uint_8)TotNbMegaFLOP() 
           << " VisibMtx_FileNum=" << numfile_ << endl;
    }
    SaveVisibilityMatrix();
  }

// kpair=numero sequentiel de la paire: 0->(0,0), 1->(0,1), 2->(0,2), 3->(0,3), 4->(1,1), 5->(1,2) ...
  sa_size_t kpair=0;
  sa_size_t k=0;     // numero de ligne dans la matrice des visibilites 
  for(size_t i=0; i<vpdata_.size(); i++) {
    for(size_t j=i; j<vpdata_.size(); j++) {
      kpair++;
      if (kpair<(pairst_+1))  continue;
      if (kpair>=(pairst_+nbpairs_+1))  break;
      if (fgpimp_&&(i!=j)&&((i+j)%2==0))  continue;  // calcul des visib avec numero pair-impair + autocorrel
      TVector< complex<r_4> > vis = vismtx_.Row(k);  k++;   

      if (fgdataraw_) {  // Donnees firmware RAW apres TF soft 
        for(sa_size_t f=1; f<vis.Size(); f++) {
          vis(f) += vpdatar_[i][f] * conj(vpdatar_[j][f]);
        }
      }
      else {   // donnees firmware FFT 
        for(sa_size_t f=1; f<vis.Size(); f++) {
          vis(f) += complex<r_4>((r_4)vpdata_[i][f].realB(), (r_4)vpdata_[i][f].imagB()) * 
            complex<r_4>((r_4)vpdata_[j][f].realB(), -(r_4)vpdata_[j][f].imagB());
        }
      }
      nb_flop_ += (8.*(r_8)(vis.Size()-1));
    }
  }

  npaqcumul_++;
  moyfc_ += curfc_[0];
  moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
  if ((fgcktt_)&&(calcid_==0))  CheckTimeTag();
  totnbpaq_++;
  return 0;
}
/* --Methode-- */
void BRVisibilityCalculator::SaveVisibilityMatrix()
{
  if (totnbpaq_ > 0) {
    moyfc_/=npaqcumul_;
    moytt_/=npaqcumul_;
    UpdateVisMtxInfo();    // add/update keywords in the Info DVList 
    // ATTENTION : Matrice visibilites non moyennee 
    char nfile[48];
    if (nbcalc_==1) 
      sprintf(nfile,"vismtx%d.%s",numfile_,OutFileExtension());
    else 
      sprintf(nfile,"vismtx_%d_%d.%s",(int)calcid_,numfile_,OutFileExtension());
    string flnm=outpath_+nfile;
    if (fgfitsout_) {  // Ecriture au format FITS
      FitsInOutFile fo(flnm, FitsInOutFile::Fits_Create);
      TArray<r_4> arvismtx = ArrCastC2R(vismtx_);
      arvismtx.Info()=vismtx_.Info();
      fo << arvismtx;
    }
    else {  // Format PPF 
      POutPersist po(flnm);
      po << vismtx_;
    }
    if (fgvisdt_) FillVisibTable(moyfc_, moytt_);
    numfile_++;
  }
  vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
  npaqcumul_=0;
  threshold_tt_=vpaq_[0].TimeTag();
  moyfc_=moytt_=0.;
  first_fc_=curfc_[0];
  first_tt_=vpaq_[0].TimeTag();
  // first_tmstamp_.SetNow();  // Current date and time 
  first_tmstamp_=cts_;  // Current date and time 
  
  return;
}

/* --Methode-- */
void BRVisibilityCalculator::UpdateVisMtxInfo()
{
  string ikey,ikdesc;
  ikey="DATEOBS";  ikdesc=" Date,Time corresponding to TimeTagFirst";
  vismtx_.Info().SetS(ikey,first_tmstamp_.ToString());
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="FirstFC";  ikdesc="First FrameCounter";
  vismtx_.Info().SetI(ikey,first_fc_);
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="FirstTT";  ikdesc="First TimeTag";
  vismtx_.Info().SetI(ikey,(first_tt_-ttfirst_[0]));
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="LastFC";  ikdesc="Last FrameCounter";
  vismtx_.Info().SetI(ikey,curfc_[0]);
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="LastTT";  ikdesc="Last TimeTag";
  vismtx_.Info().SetI(ikey,vpaq_[0].TimeTag()-ttfirst_[0]);
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="MeanFC";  ikdesc="Mean FrameCounter";
  vismtx_.Info().SetD(ikey,moyfc_); 
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="MeanTT";  ikdesc="Mean TimeTag";
  vismtx_.Info().SetD(ikey,moytt_);
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="NPAQSUM";  ikdesc="Number of paquets summed";
  vismtx_.Info().SetI(ikey,npaqcumul_);
  vismtx_.Info().SetComment(ikey,ikdesc);
  ikey="DELTIME";  ikdesc="visib cumul time (seconds)";
  vismtx_.Info().SetD(ikey,(vpaq_[0].TimeTag()-first_tt_)/1.25e8);
  vismtx_.Info().SetComment(ikey,ikdesc);
}


/* --Methode-- */
int BRVisibilityCalculator::execute(int tid)
{
  vector<TwoByteComplex*>  pvpdata(2*memgr_.NbFibres());
  vector< complex<r_4>* >  pvpdatar(2*memgr_.NbFibres());
  size_t paqsz=memgr_.PaqSize();
  size_t procpaqsz=memgr_.ProcPaqSize();
  BRPaquet ppaq(paqsz);


  sa_size_t fdelt = vismtx_.NCols()/nparthr_;
  sa_size_t fdeb = tid*fdelt;
  sa_size_t ffin = (tid+1)*fdelt;

  if (fdeb<1) fdeb=1;
  if ((ffin>vismtx_.NCols())||(tid==(nparthr_-1))) ffin=vismtx_.NCols();

  for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) {   // boucle sur les paquets d'une zone  
    if (!fgallfibok[jp])  continue;
    if (fgdataraw_) {  // Donnees firmware RAW apres TF soft 
      for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) { 
        Byte* procdatap=fprocbuff_[fib]+jp*procpaqsz;
        pvpdatar[2*fib] = reinterpret_cast< complex<r_4>* > (procdatap); 
        pvpdatar[2*fib+1] = reinterpret_cast< complex<r_4>* >(procdatap+procpaqsz/2) ;
      }
    } 
    else {  // donnees firmware FFT 
      for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
        ppaq.Set(fbuff_[fib]+jp*paqsz);
        pvpdata[2*fib] = ppaq.Data1C();
        pvpdata[2*fib+1] = ppaq.Data2C();
      }
    }


// kpair=numero sequentiel de la paire: 0->(0,0), 1->(0,1), 2->(0,2), 3->(0,3), 4->(1,1), 5->(1,2) ...
    sa_size_t kpair=0;
    sa_size_t k=0;     // numero de ligne dans la matrice des visibilites 
    for(size_t i=0; i<vpdata_.size(); i++) {
      for(size_t j=i; j<vpdata_.size(); j++) {
        kpair++;
        if (kpair<(pairst_+1))  continue;
        if (kpair>=(pairst_+nbpairs_+1))  break;
        if (fgpimp_&&(i!=j)&&((i+j)%2==0))  continue;  // calcul des visib avec numero pair-impair + autocorrel
        TVector< complex<r_4> > vis = vismtx_.Row(k);  k++;  
        if (fgdataraw_) {  // Donnees firmware RAW apres TF soft 
          for(sa_size_t f=fdeb; f<ffin; f++) {
            vis(f) += pvpdatar[i][f] * conj(pvpdatar[j][f]);
          }
        }
        else {  // donnees firmware FFT 
          for(sa_size_t f=fdeb; f<ffin; f++) {
            vis(f) += complex<r_4>((r_4)pvpdata[i][f].realB(), (r_4)pvpdata[i][f].imagB()) * 
              complex<r_4>((r_4)pvpdata[j][f].realB(), -(r_4)pvpdata[j][f].imagB());
          }
        }
        nb_flop_ += (8.*(r_8)(ffin-fdeb));
      }
    }

  }  // Fin de boucle sur les paquets 
  
  return 0;
}

/* --Methode-- */
int BRVisibilityCalculator::FillVisibTable(double fcm, double ttm)
{
  double xnt[10];
  xnt[0]=fcm;  xnt[1]=ttm/1.25e8;

  if (djf_<2) {
    for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
      for(sa_size_t jf=jf1_; jf<jf2_; jf++) {
        xnt[2]=jf;
        xnt[3]=chanpairnumall_(rv+pairst_);
        xnt[4]=vismtx_(rv,jf).real()/(r_4)(nmean_);
        xnt[5]=vismtx_(rv,jf).imag()/(r_4)(nmean_);
        visdt_.AddRow(xnt);
      }
    }
  }
  else {
    for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
      for(sa_size_t jf=jf1_; jf<jf2_; jf+=djf_) {
        r_4 moyreal=0.;
        r_4 moyimag=0.;
        sa_size_t jjfmx=jf+djf_;
        if (jjfmx > vismtx_.NCols()) jjfmx=vismtx_.NCols();
        for(sa_size_t jjf=jf; jjf<jjfmx; jjf++) {
          moyreal+=vismtx_(rv,jjf).real();
          moyimag+=vismtx_(rv,jjf).imag();
        }
        xnt[2]=jf+djf_/2;
        xnt[3]=chanpairnumall_(rv+pairst_);
        xnt[4]=moyreal/(r_4)(nmean_*djf_);
        xnt[5]=moyimag/(r_4)(nmean_*djf_);
        visdt_.AddRow(xnt);
      }
    }
  }
  return 0;
}

/* --Methode-- */
int BRVisibilityCalculator::CheckTimeTag()
{
  if (totnbpaq_==0) {
    for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  {
      vlasttt_[fib]=ttfirst_[fib];
      if (ttmtx_.NCols()>0) {
        fcmtx_(fib,totnbpaq_) = curfc_[fib]; 
        ttmtx_(fib,totnbpaq_) = vlasttt_[fib]; 
      }
    }
    return 0;
  }
  for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    int_8 ld = (int_8)vpaq_[fib].TimeTag()-(int_8)vlasttt_[fib];
    int_8 fd = (int_8)vpaq_[fib].TimeTag()-(int_8)ttfirst_[fib]-(int_8)vpaq_[0].TimeTag()+(int_8)ttfirst_[0];
    /*    if ( (ld < mindeltatt_) || (fd<-5) || (fd>5)) { vbadtt_[fib]++; vnsamett_[fib]++; }
    else {
      if (fd!=0)  vnsamett_[fib]++;
    }
    */
    if  (ld < mindeltatt_) vbadtt_[fib]++;
    else {
      if (fd != 0)  vnsamett_[fib]++;
      if ((fd<-5)||(fd>5))  vndiff5tt_[fib]++;
    }
    vlasttt_[fib]=vpaq_[fib].TimeTag();
    if (totnbpaq_<ttmtx_.NCols()) {
      fcmtx_(fib,totnbpaq_) = curfc_[fib]; 
      ttmtx_(fib,totnbpaq_) = vlasttt_[fib]; 
    }
  }
  return 0;
}

//-------------------------------------------------------------------------------
// Classe Groupe (ensemble) de Calculateur de Visibilites, tournant en parallele
//-------------------------------------------------------------------------------
/*!
   \class BRVisCalcGroup
   \ingroup TAcq

   \brief Set of visibility calculators (parallel computation)

  Each BRVisibilityCalculator object computes visibilities for a subset of pairs.
*/

/* --Methode-- */
BRVisCalcGroup::BRVisCalcGroup(size_t nbcalc, RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean,
                               uint_4 pair1, uint_4 nbpairs, bool fgpimp, size_t nthr)
  : tm_(false)
{
  if ((nbcalc<1)||(nbcalc>6)) 
    throw ParmError("BRVisCalcGroup::BRVisCalcGroup NbCalc > 6 !");
  for(size_t i=0; i<nbcalc; i++) {
    BRVisibilityCalculator * viscp=new BRVisibilityCalculator(memgr, outpath, nmean, nthr);
    viscp->DefineRank(nbcalc, i, pair1, nbpairs, fgpimp);
    viscalcp_.push_back(viscp); 
  }
}
/* --Methode-- */
BRVisCalcGroup::~BRVisCalcGroup()
{
  for(size_t i=0; i<viscalcp_.size(); i++) 
    delete viscalcp_[i];
}
/* --Methode-- */
MemZStatus BRVisCalcGroup::SetMemZAction(MemZaction mmact)
{
  MemZaction mmzas[10]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF,
                        MemZA_ProcG,MemZA_ProcH,MemZA_ProcI,MemZA_ProcJ};
  MemZStatus mmzss[10]={MemZS_ProcA,MemZS_ProcB,MemZS_ProcC,MemZS_ProcD,MemZS_ProcE,MemZS_ProcF,
                        MemZS_ProcG,MemZS_ProcH,MemZS_ProcI,MemZS_ProcJ};
  size_t ia=9999;
  for(int i=0; i<10; i++) {
    if (mmact==mmzas[i]) { ia=i; break; }
  }
  if (ia>=10)  
    throw ParmError("BRVisCalcGroup::SetMemZAction() Bad MemZaction mmact !");
  if ((ia+viscalcp_.size())>10)  
    throw ParmError("BRVisCalcGroup::SetMemZAction() MemZaction mmact too high !");
  for(size_t i=0; i<viscalcp_.size(); i++) {
    viscalcp_[i]->SetMemZAction(mmzas[ia]);   ia++;
  }
  return mmzss[ia-1];
}
/* --Methode-- */
int BRVisCalcGroup::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
{
  int rc=0;
  for(size_t i=0; i<viscalcp_.size(); i++) 
    rc=viscalcp_[i]->SelectFreqBinning(freq1, freq2, nbfreq);
  return rc;
}
/* --Methode-- */
void BRVisCalcGroup::ActivateVisDTable(bool fgfdt)
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->ActivateVisDTable(fgfdt);
}
/* --Methode-- */
void BRVisCalcGroup::SetFitsOutput()
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetFitsOutput();
} 
/* --Methode-- */
void BRVisCalcGroup::SetPPFOutput()
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetPPFOutput();
}
/* --Methode-- */
void BRVisCalcGroup::SetFFTData()
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetFFTData();
}
/* --Methode-- */
void BRVisCalcGroup::SetRawData()
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetRawData();
}  
/* --Methode-- */
void BRVisCalcGroup::SetNPaqIntervalMode(uint_4 nmean)
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetNPaqIntervalMode(nmean);
}
/* --Methode-- */
void BRVisCalcGroup::SetTimeIntervalMode(double dtime)
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetTimeIntervalMode(dtime);
}  
/* --Methode-- */
void BRVisCalcGroup::SetPrintLevel(int lev, uint_8 prtmodulo)
{
  for(size_t i=0; i<viscalcp_.size(); i++)
    viscalcp_[i]->SetPrintLevel(lev,prtmodulo);
} 
/* --Methode-- */
void BRVisCalcGroup::start()
{
  for(size_t i=0; i<viscalcp_.size(); i++) 
    viscalcp_[i]->start();
  tm_.SplitQ();
}
/* --Methode-- */
void BRVisCalcGroup::join()
{
  r_8 totflop=0.;
  for(size_t i=0; i<viscalcp_.size(); i++) {
    viscalcp_[i]->join();
    //    cout << " BRVisCalcGroup::join()/ VisibCalc[" << i << "]->TotNbMegaFLOP()=" 
    //   << viscalcp_[i]->TotNbFLOP()/1024e3 << endl;
    totflop += viscalcp_[i]->TotNbFLOP();
  }
  tm_.SplitQ();
  cout << "      ----------------------------------------------------------" << endl;
  cout << " BRVisCalcGroup::join() : Finished " << viscalcp_.size() << " VisibilityCalculator(s)" << endl;
  cout << " ... Elaspsed time: " << tm_.PartialElapsedTimems() 
       << " ms (total:" << tm_.TotalElapsedTimems() << ")"  << endl;
  double mflopsrate=totflop/(r_8)tm_.PartialElapsedTimems()/(1024.);
  cout << " ... TotalMegaFLOP= " << totflop/(1024.e3) << " @ " 
       << mflopsrate << " MFLOP/s" << " (=" << mflopsrate/(r_8)viscalcp_.size() << "/VisibCalcObject)" << endl; 
  cout << "      ----------------------------------------------------------" << endl;
  return;
}