Changeset 3956 in Sophya for trunk/AddOn/TAcq/brviscalc.cc


Ignore:
Timestamp:
Mar 2, 2011, 3:31:39 PM (15 years ago)
Author:
ansari
Message:

Amelioration du processeur de calcul de visibilite (BRVisibilityCalculator) et du programme vismfib.cc pour permettre la prise en charge des donnees raw-2c pour le calcul des visibilites et ajout de la possibilite d ecrire les fichiers de sortie (matrices de visibilites) au format FITS, Reza 02/03/2011

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/AddOn/TAcq/brviscalc.cc

    r3923 r3956  
    1212#include "pexceptions.h"
    1313#include "fioarr.h"
     14#include "fitsarrhand.h"
    1415#include "matharr.h"
     16#include "arrctcast.h"
    1517#include "timestamp.h"
    1618#include "ctimer.h"
    … …  
    2628BRVisibilityCalculator::BRVisibilityCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean, size_t nthr)
    2729  : BRBaseProcessor(memgr), paralex_(*this, nthr), nparthr_(nthr),
    28     outpath_(outpath), nmean_(nmean), nbcalc_(1), calcid_(0), vpdata_(2*memgr.NbFibres())
     30    outpath_(outpath), nmean_(nmean), nbcalc_(1), calcid_(0),
     31    vpdata_(2*memgr.NbFibres()), vpdatar_(2*memgr.NbFibres())
    2932    // , dtfos_(outpath+"visdt.fits", Fits_Create), visdt_(dtfos_, 1024, true);
    3033{ 
     34  SetFFTData();
    3135  DefineRank(1,0);
     36  SetPPFOutput();
    3237
    3338  uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
    … …  
    3540  chanpairnumall_.SetSize(maxnpairs);
    3641  chanpairsall_.SetSize(maxnpairs,2);
    37   for(size_t i=0; i<2*memgr_.NbFibres(); i++)   vpdata_[i]=NULL;
     42  for(size_t i=0; i<2*memgr_.NbFibres(); i++)   {
     43    vpdata_[i]=NULL;   vpdatar_[i]=NULL; 
     44  }
    3845  SelectPairs();
    3946
    … …  
    113120    else
    114121      SelectPairs(calcid_*npairspth+pair1, npairspth, fgpimp);
    115     MemZaction mmzas[6]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF};
     122    MemZaction mmzas[10]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF,
     123                          MemZA_ProcG,MemZA_ProcH,MemZA_ProcI,MemZA_ProcJ};
    116124    SetMemZAction(mmzas[calcid_]);
    117125    setNameId("viscalc_grp", calcid_);
    … …  
    197205
    198206  string filename;
    199   filename = outpath_+"chanum.ppf";
     207  filename = (outpath_+"chanum.")+OutFileExtension();
    200208  if (nbcalc_>1) {
    201209    char sbuff[32];
    202     sprintf(sbuff,"chanum_%d.ppf",(int)calcid_);
     210    sprintf(sbuff,"chanum_%d.%s",(int)calcid_,OutFileExtension());
    203211    filename = outpath_+sbuff;
    204212  }
    205   POutPersist poc(filename);
    206   poc << PPFNameTag("chanids") << chanids_;
    207   poc << PPFNameTag("chanpairs") << chanpairs_;
    208   poc << PPFNameTag("chanpairnum") << chanpairnum_;
    209   poc << PPFNameTag("chanpairsall") << chanpairsall_;
    210   poc << PPFNameTag("chanpairnumall") << chanpairnumall_;
    211   cout << "BRVisibilityCalculator[" << calcid_ << "]::UpdateChanIds() Channel Ids/Pairs saved to PPF file "
     213  if (fgfitsout_) {  // Ecriture au format FITS
     214    FitsInOutFile foc(filename, FitsInOutFile::Fits_Create);
     215    foc.SetNextExtensionName("chanids");
     216    foc << chanids_;
     217    foc.SetNextExtensionName("chanpairs");
     218    foc << chanpairs_;
     219    foc.SetNextExtensionName("chanpairnum");
     220    foc << chanpairnum_;
     221    foc.SetNextExtensionName("chanpairsall");
     222    foc << chanpairsall_;
     223    foc.SetNextExtensionName("chanpairnumall");
     224    foc << chanpairnumall_;
     225  }
     226  else {  // Format PPF
     227    POutPersist poc(filename);
     228    poc << PPFNameTag("chanids") << chanids_;
     229    poc << PPFNameTag("chanpairs") << chanpairs_;
     230    poc << PPFNameTag("chanpairnum") << chanpairnum_;
     231    poc << PPFNameTag("chanpairsall") << chanpairsall_;
     232    poc << PPFNameTag("chanpairnumall") << chanpairnumall_;
     233  }
     234  cout << "BRVisibilityCalculator[" << calcid_ << "]::UpdateChanIds() Channel Ids/Pairs saved to file "
    212235       << filename << endl;
    213236  cout << " ... NbVisib=NbChanPairs=" << chanpairs_.NRows() << " ChannelPairs= " ;
    … …  
    266289
    267290    size_t paqsz=memgr_.PaqSize();
     291    size_t procpaqsz=memgr_.ProcPaqSize();
    268292    bool fgrun=true;
    269293    while (fgrun) {
    … …  
    289313          break;
    290314        }
     315        if ((procpaqsz>0)&&((fprocbuff_[fib]=memgr_.GetProcMemZone(mid,fib))==NULL)) {   // cela ne devrait pas arriver non plus
     316          cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetProcMemZone("
     317               << mid << "," << fib << ") -> NULL" << endl;
     318          setRC(9);       fgrun=false;         
     319          break;
     320        }
    291321      }
    292322
    … …  
    299329          char nfile[48];
    300330          if (nbcalc_==1)
    301             sprintf(nfile,"vismtx%d.ppf",numfile_);
     331            sprintf(nfile,"vismtx%d.%s",numfile_,OutFileExtension());
    302332          else
    303             sprintf(nfile,"vismtx_%d_%d.ppf",(int)calcid_,numfile_);
     333            sprintf(nfile,"vismtx_%d_%d.%s",(int)calcid_,numfile_,OutFileExtension());
    304334          string flnm=outpath_+nfile;
    305           POutPersist po(flnm);
    306           po << vismtx_;
     335          if (fgfitsout_) {  // Ecriture au format FITS
     336            FitsInOutFile fo(flnm, FitsInOutFile::Fits_Create);
     337            TArray<r_4> arvismtx = ArrCastC2R(vismtx_);
     338            arvismtx.Info()=vismtx_.Info();
     339            fo << arvismtx;
     340          }
     341          else {  // Format PPF
     342            POutPersist po(flnm);
     343            po << vismtx_;
     344          }
    307345          if ((prtlev_>0)&&(numfile_%prtmodulo_==0)) {
    308346            cout << numfile_ << "-BRVisCalc[" << calcid_ << "/" << nbcalc_ << "]::run() NPaqProc="
    … …  
    324362          vfgok_[fib] = vpchk_[fib].Check(vpaq_[fib],curfc_[fib]);
    325363          if (!vfgok_[fib])  fgallfibok[jp]=fgokallfibers_=false;
     364          if (procpaqsz>0)    vprocpaq_[fib] = fprocbuff_[fib]+jp*procpaqsz;
    326365        }
    327366        if (fgokallfibers_)  {
    … …  
    375414{
    376415  if (totnbpaq_==0) UpdateChanIds();  // Appele ici pour etre sur que le thread de remplissage a mis l'info a jour.
    377   for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    378     vpdata_[2*fib] = vpaq_[fib].Data1C();
    379     vpdata_[2*fib+1] = vpaq_[fib].Data2C();
     416  if (fgdataraw_) {  // Donnees firmware RAW apres TF soft
     417    for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
     418      vpdatar_[2*fib] = reinterpret_cast< complex<r_4>* > (vprocpaq_[fib]);
     419      vpdatar_[2*fib+1] = reinterpret_cast< complex<r_4>* >(vprocpaq_[fib]+memgr_.ProcPaqSize()/2) ;
     420    }
     421  }
     422  else {   // donnees firmware FFT
     423    for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
     424      vpdata_[2*fib] = vpaq_[fib].Data1C();
     425      vpdata_[2*fib+1] = vpaq_[fib].Data2C();
     426    }
    380427  }
    381428 
    … …  
    388435      char nfile[48];
    389436      if (nbcalc_==1)
    390         sprintf(nfile,"vismtx%d.ppf",numfile_);
     437        sprintf(nfile,"vismtx%d.%s",numfile_,OutFileExtension());
    391438      else
    392         sprintf(nfile,"vismtx_%d_%d.ppf",(int)calcid_,numfile_);
     439        sprintf(nfile,"vismtx_%d_%d.%s",(int)calcid_,numfile_,OutFileExtension());
    393440      string flnm=outpath_+nfile;
    394       POutPersist po(flnm);
    395       po << vismtx_;
     441      if (fgfitsout_) {  // Ecriture au format FITS
     442        FitsInOutFile fo(flnm, FitsInOutFile::Fits_Create);
     443        TArray<r_4> arvismtx = ArrCastC2R(vismtx_);
     444        arvismtx.Info()=vismtx_.Info();
     445        fo << arvismtx;
     446      }
     447      else {  // Format PPF
     448        POutPersist po(flnm);
     449        po << vismtx_;
     450      }
    396451      if ((prtlev_>0)&&(numfile_%prtmodulo_==0)) {
    397452        cout << numfile_ << "-BRVisCalc[" << calcid_ << "/" << nbcalc_ << "]::Process() NPaqProc="
    … …  
    419474      if (fgpimp_&&(i!=j)&&((i+j)%2==0))  continue;  // calcul des visib avec numero pair-impair + autocorrel
    420475      TVector< complex<r_4> > vis = vismtx_.Row(k);  k++;   
    421       for(sa_size_t f=1; f<vis.Size(); f++) {
    422         vis(f) += complex<r_4>((r_4)vpdata_[i][f].realB(), (r_4)vpdata_[i][f].imagB()) *
    423           complex<r_4>((r_4)vpdata_[j][f].realB(), -(r_4)vpdata_[j][f].imagB());
     476
     477      if (fgdataraw_) {  // Donnees firmware RAW apres TF soft
     478        for(sa_size_t f=1; f<vis.Size(); f++) {
     479          vis(f) += vpdatar_[i][f] * vpdatar_[j][f];
     480        }
     481      }
     482      else {   // donnees firmware FFT
     483        for(sa_size_t f=1; f<vis.Size(); f++) {
     484          vis(f) += complex<r_4>((r_4)vpdata_[i][f].realB(), (r_4)vpdata_[i][f].imagB()) *
     485            complex<r_4>((r_4)vpdata_[j][f].realB(), -(r_4)vpdata_[j][f].imagB());
     486        }
    424487      }
    425488      nb_flop_ += (8.*(r_8)(vis.Size()-1));
    … …  
    468531{
    469532  vector<TwoByteComplex*>  pvpdata(2*memgr_.NbFibres());
     533  vector< complex<r_4>* >  pvpdatar(2*memgr_.NbFibres());
    470534  size_t paqsz=memgr_.PaqSize();
    471535  BRPaquet ppaq(paqsz);
    … …  
    481545  for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) {   // boucle sur les paquets d'une zone 
    482546    if (!fgallfibok[jp])  continue;
    483     for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
    484       ppaq.Set(fbuff_[fib]+jp*paqsz);
    485       pvpdata[2*fib] = ppaq.Data1C();
    486       pvpdata[2*fib+1] = ppaq.Data2C();
     547    if (fgdataraw_) {  // Donnees firmware RAW apres TF soft
     548      size_t procpaqsz=memgr_.ProcPaqSize();
     549      for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
     550        pvpdatar[2*fib] = reinterpret_cast< complex<r_4>* > (fprocbuff_[fib]+jp*procpaqsz);
     551        pvpdatar[2*fib+1] = reinterpret_cast< complex<r_4>* >(vprocpaq_[fib]+jp*procpaqsz+procpaqsz/2) ;
     552      }
     553    }
     554    else {  // donnees firmware FFT
     555      for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
     556        ppaq.Set(fbuff_[fib]+jp*paqsz);
     557        pvpdata[2*fib] = ppaq.Data1C();
     558        pvpdata[2*fib+1] = ppaq.Data2C();
     559      }
    487560    }
    488561
    … …  
    498571        if (fgpimp_&&(i!=j)&&((i+j)%2==0))  continue;  // calcul des visib avec numero pair-impair + autocorrel
    499572        TVector< complex<r_4> > vis = vismtx_.Row(k);  k++; 
    500         for(sa_size_t f=fdeb; f<ffin; f++) {
    501           vis(f) += complex<r_4>((r_4)pvpdata[i][f].realB(), (r_4)pvpdata[i][f].imagB()) *
    502             complex<r_4>((r_4)pvpdata[j][f].realB(), -(r_4)pvpdata[j][f].imagB());
     573        if (fgdataraw_) {  // Donnees firmware RAW apres TF soft
     574          for(sa_size_t f=fdeb; f<ffin; f++) {
     575            vis(f) += pvpdatar[i][f] * pvpdatar[j][f];
     576          }
     577        }
     578        else {  // donnees firmware FFT
     579          for(sa_size_t f=fdeb; f<ffin; f++) {
     580            vis(f) += complex<r_4>((r_4)pvpdata[i][f].realB(), (r_4)pvpdata[i][f].imagB()) *
     581              complex<r_4>((r_4)pvpdata[j][f].realB(), -(r_4)pvpdata[j][f].imagB());
     582          }
    503583        }
    504584        nb_flop_ += (8.*(r_8)(ffin-fdeb));
    … …  
    609689}
    610690/* --Methode-- */
     691MemZStatus BRVisCalcGroup::SetMemZAction(MemZaction mmact)
     692{
     693  MemZaction mmzas[10]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF,
     694                        MemZA_ProcG,MemZA_ProcH,MemZA_ProcI,MemZA_ProcJ};
     695  MemZStatus mmzss[10]={MemZS_ProcA,MemZS_ProcB,MemZS_ProcC,MemZS_ProcD,MemZS_ProcE,MemZS_ProcF,
     696                        MemZS_ProcG,MemZS_ProcH,MemZS_ProcI,MemZS_ProcJ};
     697  size_t ia=9999;
     698  for(int i=0; i<10; i++) {
     699    if (mmact==mmzas[i]) { ia=i; break; }
     700  }
     701  if (ia>=10) 
     702    throw ParmError("BRVisCalcGroup::SetMemZAction() Bad MemZaction mmact !");
     703  if ((ia+viscalcp_.size())>10) 
     704    throw ParmError("BRVisCalcGroup::SetMemZAction() MemZaction mmact too high !");
     705  for(size_t i=0; i<viscalcp_.size(); i++) {
     706    viscalcp_[i]->SetMemZAction(mmzas[ia]);   ia++;
     707  }
     708  return mmzss[ia-1];
     709}
     710/* --Methode-- */
    611711int BRVisCalcGroup::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
    612712{
    … …  
    622722    viscalcp_[i]->ActivateVisDTable(fgfdt);
    623723}
     724/* --Methode-- */
     725void BRVisCalcGroup::SetFitsOutput()
     726{
     727  for(size_t i=0; i<viscalcp_.size(); i++)
     728    viscalcp_[i]->SetFitsOutput();
     729}
     730/* --Methode-- */
     731void BRVisCalcGroup::SetPPFOutput()
     732{
     733  for(size_t i=0; i<viscalcp_.size(); i++)
     734    viscalcp_[i]->SetPPFOutput();
     735}
    624736/* --Methode-- */
    625737void BRVisCalcGroup::SetPrintLevel(int lev, uint_8 prtmodulo)
Note: See TracChangeset for help on using the changeset viewer.