Changeset 3872 in Sophya for trunk/AddOn/TAcq/brproc.cc


Ignore:
Timestamp:
Aug 28, 2010, 12:43:04 AM (15 years ago)
Author:
ansari
Message:

Ajout classe BRFFTCalculator et programme specmfib.cc, Reza 28/08/2010

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/AddOn/TAcq/brproc.cc

    r3780 r3872  
    1 #include "racquproc.h"
     1//----------------------------------------------------------------
     2// Projet BAORadio - (C) LAL/IRFU  2008-2010
     3// Classes de threads de traitement pour BAORadio
     4//----------------------------------------------------------------
    25
    36#include <stdlib.h>
     
    2629#include "brproc.h"
    2730
    28 //---------------------------------------------------------------------
    29 // Classe de traitement - calcul de visibilite pour n fibres
    30 //---------------------------------------------------------------------
    31 
    32 /* --Methode-- */
    33 BRVisibilityCalculator::BRVisibilityCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean, size_t nthr)
    34   : BRBaseProcessor(memgr), paralex_(*this, nthr), nparthr_(nthr),
    35     outpath_(outpath), nmean_(nmean), nbcalc_(1), calcid_(0), vpdata_(2*memgr.NbFibres())
    36     // , dtfos_(outpath+"visdt.fits", Fits_Create), visdt_(dtfos_, 1024, true);
    37 
    38   DefineRank(1,0);
    39 
    40   uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
    41   chanum_.SetSize(maxnpairs);
    42   sa_size_t k=0;
    43   for(size_t i=0; i<2*memgr_.NbFibres(); i++)  vpdata_[i]=NULL;
    44   for(size_t i=0; i<2*memgr_.NbFibres(); i++) {
    45     for(size_t j=i; j<2*memgr_.NbFibres(); j++) {
    46       chanum_(k) = (i+1)*100+(j+1);  k++;
    47     }
    48   }
    49   SelectPairs();
    50 
    51   //  visdt_.AddFloatColumn("mfc");
    52   visdt_.AddDoubleColumn("mfc");
    53   visdt_.AddDoubleColumn("mtt");
    54   visdt_.AddIntegerColumn("jfreq");
    55   visdt_.AddIntegerColumn("numch");
    56   visdt_.AddFloatColumn("vre");
    57   visdt_.AddFloatColumn("vim");
    58 
    59   if (nmean_ < 1)  nmean_=memgr_.NbPaquets();
    60   if (nmean_ < 1)  nmean_=1;
    61 
    62   cout << " BRVisibilityCalculator::/Info  nmean=" << nmean_ << endl;
    63 
    64   totnbpaq_=0;
    65   numfile_=0;
    66   nb_flop_=0.;
    67   moyfc_=moytt_=0.;
    68 
    69   fgallfibok=NULL;
    70   fgcktt_=false;
    71   setNameId("viscalc", 0);
    72 }
    73 
    74 /* --Methode-- */
    75 BRVisibilityCalculator::~BRVisibilityCalculator()
    76 {
    77   cout << " BRVisibilityCalculator - Visibility Datatable : " << endl;
    78   cout << visdt_;
    79   string filename;
    80   filename = outpath_+"visdt.ppf";
    81   if (nbcalc_>1) {
    82     char sbuff[32];
    83     sprintf(sbuff,"visdt_%d",(int)calcid_);
    84     filename = outpath_+sbuff;
    85   }
    86   POutPersist po(filename);
    87   po << visdt_;
    88   if (calcid_ == 0) {
    89     POutPersist poc(outpath_+"chanum.ppf");
    90     poc << chanum_;
    91  
    92     if (fgcktt_) {
    93       cout << " BRVisibilityCalculator -  Check TimeTag done: TotNPaqProc= " << totnbpaq_ << endl;
    94       for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  {
    95         cout << " BRTTCheck-Fiber[" << fib << "] NBadTT=" << vbadtt_[fib] << " NDiffTT>5="
    96              << vndiff5tt_[fib] << "  NotSameTT=" << vnsamett_[fib] << endl;
    97       }
    98       POutPersist pott(outpath_+"ttfcmtx.ppf");
    99       pott << PPFNameTag("FC") << fcmtx_;
    100       pott << PPFNameTag("TT") << ttmtx_;     
    101     }
    102   }
    103 }
    104 
    105 /* --Methode-- */
    106 void BRVisibilityCalculator::DefineRank(uint_4 nbc, uint_4 cid)
    107 {
    108   if ((nbc>6)||(cid>=nbc))
    109     throw ParmError("BRVisibilityCalculator::DefineRank() NbCalc > 6 !");
    110   nbcalc_=nbc;
    111   calcid_=cid;
    112   if (nbcalc_>1) {
    113     uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
    114     uint_4 npairs=maxnpairs/nbcalc_;
    115     if (calcid_==(nbcalc_-1))
    116       SelectPairs(calcid_*npairs, maxnpairs-calcid_*npairs);
    117     else
    118       SelectPairs(calcid_*npairs, npairs);
    119     MemZaction mmzas[6]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF};
    120     SetMemZAction(mmzas[calcid_]);
    121     setNameId("viscalc_grp", calcid_);
    122   }
    123   return ;
    124 }
    125 
    126 /* --Methode-- */
    127 uint_4 BRVisibilityCalculator::SelectPairs(uint_4 pair1, uint_4 nbpairs)
    128 {
    129   BRPaquet paq(memgr_.PaqSize());
    130   uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
    131 
    132   if (pair1 >= maxnpairs)  pair1=maxnpairs-1;
    133   if (nbpairs > maxnpairs-pair1)  nbpairs=maxnpairs-pair1;
    134   pairst_=pair1;
    135   nbpairs_=nbpairs;
    136   vismtx_.SetSize(nbpairs_, paq.DataSize()/4);
    137   return nbpairs_;
    138 }
    139 
    140 /* --Methode-- */
    141 int BRVisibilityCalculator::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
    142 {
    143   jf1_=freq1; jf2_=freq2;
    144   if ((jf1_<1)||(jf1_>=vismtx_.NCols()))  jf1_=1;
    145   if ((jf2_<1)||(jf2_>=vismtx_.NCols())||(jf2_<jf1_))  jf2_=vismtx_.NCols()-1;
    146   if (nbfreq<1) nbfreq=1;
    147   djf_=(jf2_-jf1_)/nbfreq;
    148   if (djf_<1) djf_=0;
    149   cout << " BRVisibilityCalculator::SelectFreqBinning/Info  JF1=" << jf1_
    150        << " JF2=" << jf2_ << " DJF=" << djf_ << endl;
    151 
    152 }
    153 
    154 
    155 /* --Methode-- */
    156 int BRVisibilityCalculator::ActivateTimeTagCheck(uint_8 maxnpaq)
    157 {
    158   mindeltatt_=memgr_.PaqSize()/2;
    159   if (mindeltatt_<1) mindeltatt_=1;
    160   fcmtx_.SetSize(memgr_.NbFibres(), maxnpaq);
    161   ttmtx_.SetSize(memgr_.NbFibres(), maxnpaq);
    162   vlasttt_.resize(memgr_.NbFibres(), 0);
    163   vbadtt_.resize(memgr_.NbFibres(), 0);
    164   vnsamett_.resize(memgr_.NbFibres(), 0);
    165   vndiff5tt_.resize(memgr_.NbFibres(), 0);
    166 
    167   fgcktt_=true;
    168   cout << " BRVisibilityCalculator::ActivateTimeTagCheck() - TT/Fc matrix NCols=" << maxnpaq
    169        << " MinDeltaTT=" << mindeltatt_ << endl;
    170 
    171   return 0;
    172 }
    173 
    174 
    175 /* --Methode-- */
    176 void BRVisibilityCalculator::run()
    177 {
    178   cout << " BRVisibilityCalculator[" << calcid_ << "/" << nbcalc_
    179        << "]::run() - Starting " << " NFibers=" << memgr_.NbFibres()
    180        << " NChan=" << 2*memgr_.NbFibres() << " NPairs=" << nbpairs_ << " First:" << pairst_ << endl;   
    181  
    182   if (nparthr_ < 2)  return BRBaseProcessor::run();
    183   // Execution multithread parallele
    184   setRC(1);     
    185   int rc=0;
    186   try {
    187     if ((nmean_%memgr_.NbPaquets())!=0) {
    188       uint_4 mnmean = (nmean_/memgr_.NbPaquets()+1)*memgr_.NbPaquets();
    189       cout << " BRVisibilityCalculator::run()/Info changing nmean=" << nmean_  << " to multiple of"
    190            << " memgr_.NbPaquets() -> " << mnmean << endl;
    191       nmean_=mnmean;
    192     }
    193     paralex_.SetParallelTask(*this);
    194     cout << " BRVisibilityCalculator::run()/Info : starting ParallelExecutor with nThreads="
    195          << paralex_.nThreads() << " ... " << endl;
    196     paralex_.start();
    197 
    198     fgallfibok = new bool[memgr_.NbPaquets()];
    199 
    200     size_t paqsz=memgr_.PaqSize();
    201     bool fgrun=true;
    202     while (fgrun) {
    203       if (stop_) break;
    204       if (memgr_.GetRunState() == MemZR_Stopped) break;
    205       int mid = memgr_.FindMemZoneId(mmact_);  // (MemZA_ProcA);
    206       //      Byte* buffg = memgr_.GetMemZone(mid);
    207       //      if (buffg == NULL) {
    208       if (mid < 0) {
    209         cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR FindMemZoneId("
    210              << (int)mmact_ << ") ->" << mid << ") -> NULL" << endl;
    211         setRC(7);      fgrun=false;             
    212         break; 
    213       }
    214       for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    215         fbuff_[fib] = memgr_.GetMemZone(mid,fib);
    216         if (fbuff_[fib] == NULL) { // cela ne devrait pas arriver
    217           cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetMemZone(" << mid << "," << fib << ") -> NULL" << endl;
    218           setRC(9);       fgrun=false;         
    219           break;
    220         }
    221       }
    222 
    223       if (totnbpaq_%nmean_ == 0) {
    224         if (totnbpaq_ > 0) {
    225           moyfc_/=nmean_;
    226           moytt_/=nmean_;
    227           vismtx_.Info()["MeanFC"] = moyfc_;
    228           vismtx_.Info()["MeanTT"] = moytt_;
    229           vismtx_.Info()["NPAQSUM"] = nmean_;
    230          
    231           // ATTENTION : Matrice visibilites non moyennee
    232           char nfile[48];
    233           if (nbcalc_==1)
    234             sprintf(nfile,"vismtx%d.ppf",numfile_);
    235           else
    236             sprintf(nfile,"vismtx_%d_%d.ppf",(int)calcid_,numfile_);
    237           string flnm=outpath_+nfile;
    238           POutPersist po(flnm);
    239           po << vismtx_;
    240           cout << numfile_ << "-BRVisibilityCalculator[" << calcid_ << "]::run() NPaqProc="
    241                << totnbpaq_ << "  -> Visibility Matrix in " << flnm << endl;
    242           FillVisibTable(moyfc_, moytt_);
    243           numfile_++;
    244         }
    245         vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
    246         moyfc_=moytt_=0.;
    247       }
    248 
    249       for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) {   // boucle sur les paquets d'une zone 
    250         fgallfibok[jp]=fgokallfibers_=true;
    251         for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    252           vpaq_[fib].Set(fbuff_[fib]+jp*paqsz);
    253           vfgok_[fib] = vpchk_[fib].Check(vpaq_[fib],curfc_[fib]);
    254           if (!vfgok_[fib])  fgallfibok[jp]=fgokallfibers_=false;
    255         }
    256         if (fgokallfibers_)  {
    257           if (totprocnpaq_==0) {
    258             for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  {
    259               fcfirst_[fib]=curfc_[fib];
    260               ttfirst_[fib]=vpaq_[fib].TimeTag();
    261             }
    262           }
    263           totprocnpaq_++;
    264           moyfc_ += curfc_[0];
    265           moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
    266           if ((fgcktt_)&&(calcid_==0))  CheckTimeTag();
    267           totnbpaq_++;
    268         }
    269       }  // Fin de boucle sur les paquets
    270      
    271       // Execution parallele  pour calcul des visibilites par bandes de frequence
    272       int rcpex=paralex_.execute();
    273       if (rcpex!=0)  cout << " BRVisibilityCalculator[" << calcid_ << "]::run() / Error Rc[paralex_.execute()]=" << rcpex << endl;
    274 
    275       memgr_.FreeMemZone(mid, mmsta_);  // (MemZS_ProcA);
    276     } // Fin de boucle sur les zones a traiter
    277     //------------------------------------
    278     cout << " --------- END BRVisibilityCalculator[" << calcid_ << "]::run() , TotNbProcPaq=" << totprocnpaq_ << endl;
    279     /*
    280     for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  vpchk_[fib].Print();
    281     cout << " ------------------------------------ " << endl;
    282     */
    283     delete[] fgallfibok;
    284   }
    285   catch (std::exception& exc) {
    286     cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched std::exception " << exc.what() << endl;
    287     setRC(98); 
    288     return;
    289   }
    290   catch(...) {
    291     cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched unknown ... exception " << endl;
    292     setRC(99); 
    293     return;
    294   }
    295  
    296 }
    297 
    298 /* --Methode-- */
    299 int BRVisibilityCalculator::Process()
    300 {
    301    
    302   for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    303     vpdata_[2*fib] = vpaq_[fib].Data1C();
    304     vpdata_[2*fib+1] = vpaq_[fib].Data2C();
    305   }
    306  
    307   if (totnbpaq_%nmean_ == 0) {
    308     if (totnbpaq_ > 0) {
    309       moyfc_/=nmean_;
    310       moytt_/=nmean_;
    311       vismtx_.Info()["MeanFC"] = moyfc_;
    312       vismtx_.Info()["MeanTT"] = moytt_;
    313       vismtx_.Info()["NPAQSUM"] = nmean_;
    314 
    315       // ATTENTION : Matrice visibilites non moyennee
    316       char nfile[48];
    317       if (nbcalc_==1)
    318         sprintf(nfile,"vismtx%d.ppf",numfile_);
    319       else
    320         sprintf(nfile,"vismtx_%d_%d.ppf",(int)calcid_,numfile_);
    321       string flnm=outpath_+nfile;
    322       POutPersist po(flnm);
    323       po << vismtx_;
    324       cout << numfile_ << "-BRVisibilityCalculator::Process() NPaqProc="
    325            << totnbpaq_ << "  -> Visibility Matrix in " << flnm << endl;
    326       FillVisibTable(moyfc_, moytt_);
    327       numfile_++;
    328     }
    329     vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
    330     moyfc_=moytt_=0.;
    331   }
    332 
    333   sa_size_t k=0;
    334   for(size_t i=0; i<vpdata_.size(); i++) {
    335     for(size_t j=i; j<vpdata_.size(); j++) {
    336       size_t kpair=i*vpdata_.size()+j;
    337       if (kpair<pairst_)  continue;
    338       if (kpair>=(pairst_+nbpairs_))  break;
    339       TVector< complex<r_4> > vis = vismtx_.Row(k);   k++;
    340       for(sa_size_t f=1; f<vis.Size(); f++) {
    341         vis(f) += complex<r_4>((r_4)vpdata_[i][f].realB(), (r_4)vpdata_[i][f].imagB()) *
    342           complex<r_4>((r_4)vpdata_[j][f].realB(), -(r_4)vpdata_[j][f].imagB());
    343       }
    344       nb_flop_ += (7.*(r_8)vis.Size());
    345     }
    346   }
    347 
    348   moyfc_ += curfc_[0];
    349   moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
    350   if ((fgcktt_)&&(calcid_==0))  CheckTimeTag();
    351   totnbpaq_++;
    352   return 0;
    353 }
    354 
    355 /* --Methode-- */
    356 int BRVisibilityCalculator::execute(int tid)
    357 {
    358   vector<TwoByteComplex*>  pvpdata(2*memgr_.NbFibres());
    359   size_t paqsz=memgr_.PaqSize();
    360   BRPaquet ppaq(paqsz);
    361 
    362   sa_size_t fdelt = vismtx_.NCols()/nparthr_;
    363   sa_size_t fdeb = tid*fdelt;
    364   sa_size_t ffin = (tid+1)*fdelt;
    365 
    366   if (fdeb<1) fdeb=1;
    367   if ((ffin>vismtx_.NCols())||(tid==(nparthr_-1))) ffin=vismtx_.NCols();
    368 
    369   for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) {   // boucle sur les paquets d'une zone 
    370     if (!fgallfibok[jp])  continue;
    371     for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    372       ppaq.Set(fbuff_[fib]+jp*paqsz);
    373       pvpdata[2*fib] = ppaq.Data1C();
    374       pvpdata[2*fib+1] = ppaq.Data2C();
    375     }
    376     sa_size_t k=0;
    377     for(size_t i=0; i<vpdata_.size(); i++) {
    378       for(size_t j=i; j<vpdata_.size(); j++) {
    379         size_t kpair=i*vpdata_.size()+j;
    380         if (kpair<pairst_)  continue;
    381         if (kpair>=(pairst_+nbpairs_))  break;
    382         TVector< complex<r_4> > vis = vismtx_.Row(k);   k++;
    383         for(sa_size_t f=fdeb; f<ffin; f++) {
    384           vis(f) += complex<r_4>((r_4)pvpdata[i][f].realB(), (r_4)pvpdata[i][f].imagB()) *
    385             complex<r_4>((r_4)pvpdata[j][f].realB(), -(r_4)pvpdata[j][f].imagB());
    386         }
    387         nb_flop_ += (7.*(r_8)(ffin-fdeb));
    388       }
    389     }
    390 
    391   }  // Fin de boucle sur les paquets
    392  
    393   return 0;
    394 }
    395 
    396 /* --Methode-- */
    397 int BRVisibilityCalculator::FillVisibTable(double fcm, double ttm)
    398 {
    399   double xnt[10];
    400   xnt[0]=fcm;  xnt[1]=ttm/1.25e8;
    401 
    402   if (djf_<2) {
    403     for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
    404       for(sa_size_t jf=jf1_; jf<jf2_; jf++) {
    405         xnt[2]=jf;
    406         xnt[3]=chanum_(rv+pairst_);
    407         xnt[4]=vismtx_(rv,jf).real()/(r_4)(nmean_);
    408         xnt[5]=vismtx_(rv,jf).imag()/(r_4)(nmean_);
    409         visdt_.AddRow(xnt);
    410       }
    411     }
    412   }
    413   else {
    414     for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
    415       for(sa_size_t jf=jf1_; jf<jf2_; jf+=djf_) {
    416         r_4 moyreal=0.;
    417         r_4 moyimag=0.;
    418         sa_size_t jjfmx=jf+djf_;
    419         if (jjfmx > vismtx_.NCols()) jjfmx=vismtx_.NCols();
    420         for(sa_size_t jjf=jf; jjf<jjfmx; jjf++) {
    421           moyreal+=vismtx_(rv,jjf).real();
    422           moyimag+=vismtx_(rv,jjf).imag();
    423         }
    424         xnt[2]=jf+djf_/2;
    425         xnt[3]=chanum_(rv+pairst_);
    426         xnt[4]=moyreal/(r_4)(nmean_*djf_);
    427         xnt[5]=moyimag/(r_4)(nmean_*djf_);
    428         visdt_.AddRow(xnt);
    429       }
    430     }
    431   }
    432   return 0;
    433 }
    434 
    435 /* --Methode-- */
    436 int BRVisibilityCalculator::CheckTimeTag()
    437 {
    438   if (totnbpaq_==0) {
    439     for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++)  {
    440       vlasttt_[fib]=ttfirst_[fib];
    441       if (ttmtx_.NCols()>0) {
    442         fcmtx_(fib,totnbpaq_) = curfc_[fib];
    443         ttmtx_(fib,totnbpaq_) = vlasttt_[fib];
    444       }
    445     }
    446     return 0;
    447   }
    448   for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    449     int_8 ld = (int_8)vpaq_[fib].TimeTag()-(int_8)vlasttt_[fib];
    450     int_8 fd = (int_8)vpaq_[fib].TimeTag()-(int_8)ttfirst_[fib]-(int_8)vpaq_[0].TimeTag()+(int_8)ttfirst_[0];
    451     /*    if ( (ld < mindeltatt_) || (fd<-5) || (fd>5)) { vbadtt_[fib]++; vnsamett_[fib]++; }
    452     else {
    453       if (fd!=0)  vnsamett_[fib]++;
    454     }
    455     */
    456     if  (ld < mindeltatt_) vbadtt_[fib]++;
    457     else {
    458       if (fd != 0)  vnsamett_[fib]++;
    459       if ((fd<-5)||(fd>5))  vndiff5tt_[fib]++;
    460     }
    461     vlasttt_[fib]=vpaq_[fib].TimeTag();
    462     if (totnbpaq_<ttmtx_.NCols()) {
    463       fcmtx_(fib,totnbpaq_) = curfc_[fib];
    464       ttmtx_(fib,totnbpaq_) = vlasttt_[fib];
    465     }
    466   }
    467   return 0;
    468 }
    469 
    470 //-------------------------------------------------------------------------------
    471 // Classe Groupe (ensemble) de Calculateur de Visibilites, tourant en parallele
    472 //-------------------------------------------------------------------------------
    473 
    474 /* --Methode-- */
    475 BRVisCalcGroup::BRVisCalcGroup(size_t nbcalc, RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean, size_t nthr)
    476   : tm_(false)
    477 {
    478   if ((nbcalc<1)||(nbcalc>6))
    479     throw ParmError("BRVisCalcGroup::BRVisCalcGroup NbCalc > 6 !");
    480   for(size_t i=0; i<nbcalc; i++) {
    481     BRVisibilityCalculator * viscp=new BRVisibilityCalculator(memgr, outpath, nmean, nthr);
    482     viscp->DefineRank(nbcalc, i);
    483     viscalcp_.push_back(viscp);
    484   }
    485 }
    486 /* --Methode-- */
    487 BRVisCalcGroup::~BRVisCalcGroup()
    488 {
    489   for(size_t i=0; i<viscalcp_.size(); i++)
    490     delete viscalcp_[i];
    491 }
    492 /* --Methode-- */
    493 int BRVisCalcGroup::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
    494 {
    495   int rc=0;
    496   for(size_t i=0; i<viscalcp_.size(); i++)
    497     rc=viscalcp_[i]->SelectFreqBinning(freq1, freq2, nbfreq);
    498   return rc;
    499 }
    500 /* --Methode-- */
    501 void BRVisCalcGroup::start()
    502 {
    503   for(size_t i=0; i<viscalcp_.size(); i++)
    504     viscalcp_[i]->start();
    505   tm_.SplitQ();
    506 }
    507 /* --Methode-- */
    508 void BRVisCalcGroup::join()
    509 {
    510   r_8 totflop=0.;
    511   for(size_t i=0; i<viscalcp_.size(); i++) {
    512     viscalcp_[i]->join();
    513     totflop += viscalcp_[i]->TotNbFLOP();
    514   }
    515   tm_.SplitQ();
    516   cout << "      ----------------------------------------------------------" << endl;
    517   cout << " BRVisCalcGroup::join() : Finished , Elaspsed time: " << tm_.PartialElapsedTimems()
    518        << " ms (total:" << tm_.TotalElapsedTimems() << ")" << endl;
    519   cout << " ... TotalMegaFLOP= " << totflop/(1024.e3) << " @ "
    520        << totflop/(r_8)tm_.PartialElapsedTimems()/(1024) << " MFLOP/s" << endl;
    521   cout << "      ----------------------------------------------------------" << endl;
    522   return;
    523 }
    52431
    52532
     
    52835//---------------------------------------------------------------------
    52936/* --Methode-- */
    530 BRMeanSpecCalculator::BRMeanSpecCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean)
    531   : BRBaseProcessor(memgr), outpath_(outpath), nmean_(nmean)
    532 {
     37BRMeanSpecCalculator::BRMeanSpecCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean,
     38                                           bool fgdatafft, bool fgsinglechan)
     39  : BRBaseProcessor(memgr), outpath_(outpath), nmean_(nmean),
     40    fgdatafft_(fgdatafft), fgsinglechannel_(fgsinglechan)
     41{
     42  setNameId("meanSpecCalc",1);
    53343  BRPaquet paq(memgr_.PaqSize());
    534   mspecmtx_.SetSize(2*memgr_.NbFibres(), paq.DataSize()/4);
     44  if (fgsinglechannel_)
     45    mspecmtx_.SetSize(memgr_.NbFibres(), paq.DataSize()/2);
     46  else
     47    mspecmtx_.SetSize(2*memgr_.NbFibres(), paq.DataSize()/4);
     48  mspecmtx_=(r_8)(0.);
    53549  numfile_=0;
     50  nbpaq4mean_=0;
    53651  totnbpaq_=0;
    53752}
     
    54055BRMeanSpecCalculator::~BRMeanSpecCalculator()
    54156{
    542 }
    543 
     57  if (nbpaq4mean_>1)  SaveSpectra();
     58}
     59
     60
     61static inline r_4 Zmod2(complex<r_4> z)
     62{ return (z.real()*z.real()+z.imag()*z.imag()); }
    54463
    54564/* --Methode-- */
     
    54766{
    54867     
    549   if (totnbpaq_%nmean_ == 0) {
    550     if (totnbpaq_ > 0) {
    551       mspecmtx_.Info()["NPAQSUM"] = nmean_;
    552       mspecmtx_ /= (double)nmean_;
    553       char nfile[32];
    554       sprintf(nfile,"mspecmtx%d.ppf",numfile_);
    555       string flnm=outpath_+nfile;
    556       POutPersist po(flnm);
    557       po << mspecmtx_;
    558       cout << numfile_ << "-BRMeanSpecCalculator::Process() NPaqProc="
    559            << totnbpaq_ << "  -> Mean spectra Matrix in " << flnm << endl;
    560       numfile_++;
     68  if (nbpaq4mean_==nmean_)  SaveSpectra();
     69   
     70  if (fgdatafft_) {  // Donnees firmware FFT
     71    for(sa_size_t i=0; i<(size_t)mspecmtx_.NRows(); i++) {
     72      TwoByteComplex* zp=NULL;
     73      if (fgsinglechannel_) {
     74        zp=vpaq_[i].Data1C();
     75      }
     76      else {
     77        zp=vpaq_[i/2].Data1C();
     78        if (i%2==1)  zp=vpaq_[i/2].Data2C();
     79      }
     80      TVector< r_4 > spec = mspecmtx_.Row(i);
     81      for(sa_size_t f=1; f<spec.Size(); f++)  spec(f) += zp[f].module2F();
    56182    }
    562     mspecmtx_ = (r_8)(0.);
    563   }
    564 
    565   sa_size_t k=0;
    566   for(size_t i=0; i<(size_t)2*memgr_.NbFibres(); i++) {
    567     TwoByteComplex* zp=vpaq_[i/2].Data1C();
    568     if (i%2==1)  zp=vpaq_[i/2].Data2C();
    569     TVector< r_4 > spec = mspecmtx_.Row(k);   k++;
    570     for(sa_size_t f=1; f<spec.Size(); f++) {
    571       spec(f) += zp[f].module2F();
    572       }
     83  }
     84  else {  // Donnees RAW qui ont du etre processe par BRFFTCalculator
     85    for(sa_size_t i=0; i<(size_t)mspecmtx_.NRows(); i++) {
     86      complex<ODT>* zp=NULL;
     87      if (fgsinglechannel_) {
     88        zp=reinterpret_cast< complex<ODT>* > (vprocpaq_[i]);
     89      }
     90      else {
     91        zp=reinterpret_cast< complex<ODT>* > (vprocpaq_[i/2]);
     92        if (i%2==1)  zp= reinterpret_cast< complex<ODT>* >(vprocpaq_[i/2]+memgr_.ProcPaqSize()/2) ;
     93      }
     94      TVector< r_4 > spec = mspecmtx_.Row(i);
     95      for(sa_size_t f=1; f<spec.Size(); f++)  spec(f) += Zmod2(zp[f]);
     96     
    57397    }
    574 
    575   totnbpaq_++;
     98  }
     99  nbpaq4mean_++;   totnbpaq_++;
    576100  return 0;
     101}
     102
     103/* --Methode-- */
     104void BRMeanSpecCalculator::SaveSpectra()
     105{
     106  mspecmtx_.Info()["NPAQSUM"] = nbpaq4mean_;
     107  mspecmtx_ /= (double)nbpaq4mean_;
     108  char nfile[32];
     109  sprintf(nfile,"mspecmtx%d.ppf",numfile_);
     110  string flnm=outpath_+nfile;
     111  POutPersist po(flnm);
     112  po << mspecmtx_;
     113  cout << numfile_ << "-BRMeanSpecCalculator::SaveSpectra() NPaqProc="
     114       << totnbpaq_ << "  -> Mean spectra Matrix in " << flnm << endl;
     115  numfile_++;  nbpaq4mean_=0;
     116  mspecmtx_ = (r_8)(0.);
     117  return;
     118}
     119
     120//---------------------------------------------------------------------
     121// Classe de thread de calcul de FFT sur donnees RAW
     122//---------------------------------------------------------------------
     123/* --Methode-- */
     124BRFFTCalculator::BRFFTCalculator(RAcqMemZoneMgr& memgr, bool fgsinglechannel)
     125  : BRBaseProcessor(memgr), fgsinglechannel_(fgsinglechannel), totnbfftpaq_(0)
     126{
     127  BRPaquet paq(memgr_.PaqSize());
     128  setNameId("FFTCalc",2);
     129  ffts_.SetInDataSize((fgsinglechannel_)?paq.DataSize():paq.DataSize()/2);
     130}
     131
     132/* --Methode-- */
     133BRFFTCalculator::~BRFFTCalculator()
     134{
     135}
     136
     137
     138/* --Methode-- */
     139int BRFFTCalculator::Process()
     140{
     141  for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
     142    ffts_.DoFFT( reinterpret_cast<IDT *>(vpaq_[fib].Data1() ), 
     143                 reinterpret_cast< complex<ODT>* > (vprocpaq_[fib]) );
     144    totnbfftpaq_++;
     145    if ( fgsinglechannel_ ) continue;
     146    ffts_.DoFFT( reinterpret_cast<IDT *>(vpaq_[fib].Data2() ), 
     147                 reinterpret_cast< complex<ODT>* > (vprocpaq_[fib]+memgr_.ProcPaqSize()/2) );
     148    totnbfftpaq_++;
     149  }   
     150  return 0;
     151}
     152
     153
     154//-------------------------------------------------------------------------
     155// Classe WBRFFT : Calcul de TF sur donnees brutes (firmware RAW)
     156//-------------------------------------------------------------------------
     157ZMutex* WBRFFT::mtx_fftwp_ = NULL;
     158
     159/* --Methode-- */
     160WBRFFT::WBRFFT(uint_4 sz)
     161  : sz_(sz)
     162{
     163  if (mtx_fftwp_ == NULL) mtx_fftwp_ = new ZMutex;
     164  if (sz>0)  SetInDataSize(sz);
     165}
     166
     167/* --Methode-- */
     168WBRFFT::~WBRFFT()
     169{
     170}
     171
     172/* --Methode-- */
     173void WBRFFT::SetInDataSize(uint_4 sz)
     174{
     175  sz_ = sz;
     176  if (sz_<1) return;
     177  inp.SetSize(sz);
     178  outfc.SetSize(sz/2+1);
     179  mtx_fftwp_->lock();
     180  myplan_ = fftwf_plan_dft_r2c_1d(inp.Size(), inp.Data(),
     181                       (fftwf_complex*)outfc.Data(), FFTW_ESTIMATE);
     182  mtx_fftwp_->unlock();
     183}
     184
     185/* --Methode-- */
     186void WBRFFT::DoFFT( IDT *indata, complex<ODT> * ofc)
     187{
     188  if (sz_<1) return;
     189  for(uint_4 k=0; k<inp.Size(); k++)   inp(k)=(ODT)indata[k];
     190  fftwf_execute(myplan_);
     191  for(uint_4 k=0; k<outfc.Size(); k++)   ofc[k]=outfc(k);
     192  return;
     193}
     194
     195/* --Methode-- */
     196void WBRFFT::PrintData(IDT *indata,  complex<ODT> * ofc, uint_4 sz)
     197{
     198    cout << " --- WBRFFT::PrintData() size=" << sz << endl;
     199  for(uint_4 k=0; k<sz; k+=8) {
     200    IDT* in = indata+k;
     201    cout << " Indata[" << k << "..." << k+8 << "]= ";
     202    for(uint_4 i=0; i<8; i++)   cout << (IIDT)in[i] << "  ";
     203    cout << endl; 
     204  }
     205  cout << endl;
     206  for(uint_4 k=0; k<sz/2; k+=4) {
     207    complex< ODT>* out = ofc+k;
     208    cout << " OutFC[" << k << "..." << k+4 << "]= ";
     209    for(uint_4 i=0; i<4; i++)   cout << out[i] << "  ";
     210    cout << endl; 
     211  }
     212  return;
     213
    577214}
    578215
     
    580217//---------------------------------------------------------------
    581218// Classe thread de traitement donnees ADC avec 2 voies par frame
     219//         !!!! OBSOLETE  !!!!   
    582220//---------------------------------------------------------------
    583221
     
    612250
    613251
    614 static inline r_4 Zmod2(complex<r_4> z)
    615 { return (z.real()*z.real()+z.imag()*z.imag()); }
     252
    616253
    617254static inline string card2name_(int card)
     
    1004641}   
    1005642
     643
    1006644//---------------------------------------------------------------------
    1007645// Classe thread de traitement 2 x 2 voies/frames (Apres BRProcA2C)
     646//         !!!! OBSOLETE  !!!!   
    1008647//---------------------------------------------------------------------
    1009648
Note: See TracChangeset for help on using the changeset viewer.