Changeset 1484 in Sophya


Ignore:
Timestamp:
Apr 30, 2001, 5:20:03 PM (24 years ago)
Author:
ansari
Message:

Ameliorations diverses - Reza 30/4/2001

Location:
trunk/ArchTOIPipe
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/ArchTOIPipe/Kernel/toiseqbuff.cc

    r1464 r1484  
    1111
    1212TOISeqBuffered::TOISeqBuffered(int wsz) {
     13  data = NULL;
     14  flags = NULL;
    1315  AllocBuffer(wsz);
     16  setName("toiseqbuff");
    1417}
    1518
    1619TOISeqBuffered::TOISeqBuffered(string nm, int wsz) {
     20  data = NULL;
     21  flags = NULL;
    1722  AllocBuffer(wsz);
    1823  setName(nm);
    … …  
    2934  wsize = wsz;
    3035  buffsize = 2*wsz;
     36  if (data)   delete[] data;
     37  if (flags)  delete[] flags;
    3138  data = new double[buffsize];
    3239  flags = new int_8[buffsize];
    … …  
    4350void TOISeqBuffered::PrintStatus(ostream & os) const
    4451{
    45   os << "---TOISeqBuffered::PrintStatus() - Name=" << getName() << endl;
     52  os << "---TOISeqBuffered::PrintStatus() - Name=" << getName()
     53     << "\n  WindowSize= " << wsize << " BufferSize= " << buffsize << endl;
    4654  os << "Index: FirstIn= " << getFirstIn() << " LastIn= " << getLastIn()
    4755     << "  Total= " << getLastIn()-getFirstIn()+1 << endl;
    … …  
    7280
    7381void TOISeqBuffered::doWontNeedBefore(int i) {
    74   next_out = i;
     82  // $CHECK$  Reza 30/4/2001 - Je ne sais pas a quoi ca sert !
     83  //  next_out = i; $CHECK$  Reza 30/4/2001
    7584}
    7685

  • trunk/ArchTOIPipe/Kernel/toiseqbuff.h

    r1464 r1484  
    55#include "toi.h"
    66
     7// -------------------  TOISeqBuffered ---------------------------
    78// Classe de TOI avec buffer, et echantillonnage regulier.
    89// Pour le moment au moins,
    910// il faut que les providers fassent arriver les donnees par samplenum croissant.
     11// ---------------------------------------------------------------
     12
    1013class TOISeqBuffered : public TOIRegular {
    1114public:
    … …  
    1316  TOISeqBuffered(string nm, int wsz=8192);
    1417  virtual ~TOISeqBuffered();
     18
     19 
     20  inline  void SetBufferSize(int wsz)  // ATTENTION - Ne doit pas etre appele
     21               { AllocBuffer(wsz); }   // apres le demarrage des threads
    1522
    1623  virtual void PrintStatus(ostream & os) const;

  • trunk/ArchTOIPipe/ProcWSophya/simtoipr.cc

    r1483 r1484  
    162162    uint_8 fgcur;
    163163    bool fgokcur=false;
     164
     165    int sx_refresh_count = 0;
     166    int sx_refresh_count_max = 16*wsize;
     167
    164168    // Boucle sur les sampleNum
    165 
    166169    int knext;
    167170    int kfin = sne-snb;
    … …  
    199202        }
    200203        else {
    201           s += (valadd-valsub);
    202           s2 += (valadd*valadd-valsub*valsub);
     204          if (sx_refresh_count >= sx_refresh_count_max) {
     205            //  On recalcule la somme
     206            s = vas.Sum();
     207            s2 = vas.SumX2();
     208            sx_refresh_count = 0;
     209          }
     210          else {
     211            s += (valadd-valsub);
     212            s2 += (valadd*valadd-valsub*valsub);
     213            sx_refresh_count++;
     214          }
    203215          mean = s/wsize;
    204216          sigma = sqrt(s2/wsize-mean*mean);
    … …  
    616628    throw ParmError("SimpleFourierFilter::SimpleFourierFilter() WSize<16 !");
    617629  KeepSpectra("spectra.ppf", 0);
     630  ComputeMeanSpectra(false);
     631  totnscount = 0;
     632  totnbblock = 0;
    618633}
    619634
    … …  
    629644     << WSize() << endl;
    630645  TOIProcessor::PrintStatus(os);
    631   os << " Coeff= " << ffcoef << endl;
    632   os << " ProcessedSampleCount=" << ProcessedSampleCount() << endl;
     646  os << " Coeff (Size= " << ffcoef.Size() << "): " << endl; 
     647  for(int i=0; i<16; i++) {
     648    os << ffcoef(i) << "  " ;
     649    if (i == 7) os << endl;
     650  }
     651  os << "  .... " << endl;
     652  os << " ProcessedSampleCount=" << ProcessedSampleCount()
     653     << " NbFFTBlocks= " << totnbblock << endl;
    633654  os << " ------------------------------------------------------ " << endl;
    634655}
    … …  
    657678    throw ParmError("SimpleFourierFilter::run() Input TOI (in) not connected!");
    658679  }
    659   if (!fgout) {
    660     cerr << " SimpleFourierFilter::run() - No Output TOI connected! "
    661          << endl;
    662     throw ParmError("SimpleFourierFilter::run() No output TOI connected!");
    663   }
     680
     681  // ---- On peut utiliser cette classe pour calculer un spectre de Fourier ----
     682  //  if (!fgout) {
     683  //    cerr << " SimpleFourierFilter::run() - No Output TOI connected! "
     684  //     << endl;
     685  //    throw ParmError("SimpleFourierFilter::run() No output TOI connected!");
     686  //  }
    664687
    665688  cout << " SimpleFourierFilter::run() SNRange=" << snb << " - " << sne << endl;
    … …  
    673696    TVector<int_8> vfg(wsize);
    674697    TVector< complex<r_8> > vfft, vfftmean;
     698    Vector meanpowerspectra;
    675699    TVector< complex<r_8> > zcoef(ffcoef.Size());
    676700    zcoef = ffcoef;
    … …  
    684708    int k,i,klast;
    685709    int nks = 0;
    686     int nblk = 0;
    687710    klast = snb-1;
     711    totnbblock = 0;
    688712    // Boucle sur les sampleNum
    689713    // 1er partie, on traite par paquets de wsize
    … …  
    692716        getData(0, k+i, vin(i), vfg(i));
    693717      ffts.FFTForward(vin, vfft);
    694       if (nblk == 0)  vfftmean = vfft;
    695       else vfftmean += vfft;
    696       nblk++;
     718      if (c_meanspectra) { // Compute mean-spectra
     719        if (totnbblock == 0)  {
     720          vfftmean = vfft;
     721          meanpowerspectra.ReSize(vfft.Size());
     722          for(i=0; i<meanpowerspectra.Size(); i++)
     723            meanpowerspectra(i) = sqrt(vfft(i).real()*vfft(i).real() +
     724                                    vfft(i).imag()*vfft(i).imag() );
     725        }
     726        else {
     727          vfftmean += vfft;
     728          for(i=0; i<meanpowerspectra.Size(); i++)
     729            meanpowerspectra(i) += sqrt(vfft(i).real()*vfft(i).real() +
     730                                        vfft(i).imag()*vfft(i).imag() );
     731        }
     732      }
     733      totnbblock++;
    697734      if (nks < nb_keep) {
    698735        TVector< complex<r_8> > vfftcopie;
    699736        vfftcopie = vfft;
    700         string nomvfft = "spectra" + MuTyV(nks);
     737        string nomvfft = "spectra" + (string)MuTyV(nks);
    701738        pout.PutObject(vfftcopie, nomvfft);
    702739        nks++;
    703740      }
    704       vfft.MulElt(zcoef);
    705       ffts.FFTBackward(vfft, vout);
    706       for(i=0; i<wsize; i++)
    707         putData(0,k+i,vout(i),vfg(i));
    708       if (fgincopie)
    709         for(i=0; i<wsize; i++)
     741      if (fgout) {
     742        vfft.MulElt(zcoef);
     743        ffts.FFTBackward(vfft, vout);
     744      }
     745      for(i=0; i<wsize; i++) {
     746        if (fgout)
     747          putData(0,k+i,vout(i),vfg(i));
     748        if (fgincopie) 
    710749          putData(1, k+i, vin(i), vfg(i));
     750      }
    711751      klast+=wsize;
    712752      totnscount+=wsize;
    … …  
    720760      for(k=klast+1; k<=sne; k++) {
    721761        getData(0, k, inval, inflg);
    722         putData(0, k, inval, inflg);
     762        if (fgout) putData(0, k, inval, inflg);
    723763        if (fgincopie)
    724764          putData(1, k, inval, inflg);
    … …  
    726766      }
    727767
    728    
    729     vfftmean /= complex<r_8>((r_8)nblk, 0.);   
    730     pout.PutObject(vfftmean, "meanspectra");
    731 
    732     cout << " SimpleFourierFilter::run() - End of processing " << endl;
     768    if (c_meanspectra) {
     769      vfftmean /= complex<r_8>((r_8)totnbblock, 0.);   
     770      pout.PutObject(vfftmean, "meanspectra");
     771      meanpowerspectra /= (r_8)totnbblock;
     772      pout.PutObject(vfftmean, "meanpowerspectra");
     773    }
     774    pout.PutObject(ffcoef, "fourierfilter");
     775
     776    cout << " SimpleFourierFilter::run() - End of processing "
     777         << " NbFFTBlocks= " << totnbblock << endl;
    733778  }  // Bloc try
    734779

  • trunk/ArchTOIPipe/ProcWSophya/simtoipr.h

    r1483 r1484  
    160160  inline  void  KeepSpectra(string outname, int nb)
    161161                { outppfname = outname; nb_keep = nb; }
     162  inline  void  ComputeMeanSpectra(bool fg)
     163                {  c_meanspectra = fg; }
    162164protected:
    163165  int_8 totnscount;   // Nombre total d'echantillon processe
     166  int_8 totnbblock;   // Nombre total de blocs pour FFT
    164167  int wsize;         // Taille de fenetre de travail
    165168  Vector ffcoef;     // Coefficients du filtre
     169  bool c_meanspectra;
    166170  int nb_keep;
    167171  string outppfname;

  • trunk/ArchTOIPipe/SMakefile

    r1464 r1484  
    1717            $(AOBJ)tstrztoi.o $(AOBJ)rztoi.o \
    1818            $(AOBJ)fits2ascii.o $(AOBJ)tstmap2toi.o $(AOBJ)tsttoi2map.o
    19 EXELIST := $(AOBJ)simtst $(AOBJ)tstrztoi $(AOBJ)fits2asc $(AOBJ)tsttoi \
     19EXELIST := $(AOBJ)sfftc $(AOBJ)simtst $(AOBJ)tstrztoi $(AOBJ)fits2asc $(AOBJ)tsttoi \
    2020           $(AOBJ)tsttoi2 $(AOBJ)fits2ascii $(AOBJ)tstmap2toi $(AOBJ)tsttoi2map
    2121
    … …  
    4141$(AOBJ)tstrztoi: $(AOBJ)tstrztoi.o $(AOBJ)rztoi.o $(LIBOLIST)
    4242        $(LINK.cc)  -o $@ $^ $(LIBS)
    43 $(AOBJ)simtst: $(AOBJ)simtst.o $(AOBJ)simtoipr.o $(LIBOLIST)
     43$(AOBJ)sfftc: $(AOBJ)sfftc.o  $(LIBOLIST)
     44        $(LINK.cc)  -o $@ $^ $(LIBS)
     45$(AOBJ)simtst: $(AOBJ)simtst.o  $(LIBOLIST)
    4446        $(LINK.cc)  -o $@ $^ $(LIBS)
    4547$(AOBJ)fits2asc: $(AOBJ)fits2asc.o $(LIBOLIST)
    … …  
    5658        $(LINK.cc)  -o $@ $^ $(LIBS)
    5759
     60$(AOBJ)simtst.o: simtst.cc $(INCLIST)
     61$(AOBJ)sfftc.o: sfftc.cc $(INCLIST)
    5862$(AOBJ)tstrztoi.o: tstrztoi.cc rztoi.h $(INCLIST)
    5963$(AOBJ)rztoi.o: rztoi.cc rztoi.h $(INCLIST)
    … …  
    6367$(AOBJ)fits2ascii.o: fits2ascii.cc $(INCLIST)
    6468$(AOBJ)tstmap2toi.o: tstmap2toi.cc $(INCLIST)
    65 $(AOBJ)tsttoi2toi.o: tsttoi2map.cc $(INCLIST)
     69$(AOBJ)tsttoi2map.o: tsttoi2map.cc $(INCLIST)
    6670
    6771$(AOBJ)asciitoiwtr.o: asciitoiwtr.cc asciitoiwtr.h toi.h config.h conf.h \

  • trunk/ArchTOIPipe/TestPipes/simtst.cc

    r1483 r1484  
    22
    33----------------   Exemple d'appel  ---------------------
    4 csh> simtst -start 104385384 -end 104399964 -range -500,500
    5             -intoi boloMuV_26 -wtoi 8192 -wdegli 1024 inputbolo.fits degli.fits xx.ppf
     4csh> simtst -start 104385384 -end 104399964 -range -500,500 \
     5            -intoi boloMuV_26 -wtoi 8192 -wdegli 1024 \
     6            inputbolo.fits filt.fits xx.ppf
     7# Avec Filtre de Fourier
     8csh> simtst -start 104385384 -end 104399964 -range -500,500 \
     9            -intoi boloMuV_26 -wtoi 8192 -wdegli 1024 \
     10            -wfft 4096 -keepfft 5 -killfreq 15,4,2.5 \
     11            inputbolo.fits filtfft.fits spectre.ppf
     12# Autre exemple b545k02
     13cool> ./simtst -start 104389122 -end 104649122 -range -1000.,1000. -intoi boloMu
     14V_23 -wtoi 8192 -wdegli 1024 -wfft 4096 -keepfft 5 -killfreq 15,4,2.5 Cmv/b545k02
     152_16.00_16.49.fits b545.fits b545.ppf
     16
    617*/
    718
    … …  
    2839  cout << endl;
    2940  if (fgerr) {
    30     cout << " simtst : Argument Error ! tstrztoi -h for usage " << endl;
     41    cout << " simtst : Argument Error ! simtst -h for usage " << endl;
    3142    exit(1);
    3243  }
    3344  else {
    34     cout << "\n Usage : simtst [-dbg] [-start snb] [-end sne] \n"
    35          << "         [-wtoi sz] [-wdegli sz] [-wfft sz] [-keepfft n] \n"
     45    cout << "\n Usage : simtst [-intoi toiname] [-start snb] [-end sne] \n"
     46         << "         [-dbg] [-wtoi sz] [-wdegli sz] [-range min,max] \n"
     47         <<"          [-wfft sz] [-keepfft n] [-killfreq bf,nharm,sigf \n"
    3648         << "         inFitsName outFitsName  outPPFName \n"
    3749         << "   -dbg : sets TOISeqBuffered debug level to 1 \n"
    … …  
    4557         << "   -wfft sz : Activate Fourier filter and sets its width \n"
    4658         << "   -keepfft n : Keeps n spectra (if Fourier filter activated) \n"
     59         << "   -killfreq bf,nharm,sigf : kills nharm frequency, basefreq=bf \n"
     60            "              with a gaussian with width=sigf \n"
    4761         << endl;
    4862    exit(0);
    … …  
    7286  string intoi = "boloMuV_27";
    7387  bool fg_f_filt = false;
     88  bool fg_killfreq = false;
     89  int bf_killfreq = 1;
     90  int nharm_killfreq = 1;
     91  double sigf_killfreq = 1.;
    7492
    7593  if (narg < 4) Usage(true);
    … …  
    103121      keepfft = atoi(arg[ia+1]); ia++;
    104122    }   
     123    else if (strcmp(arg[ia],"-killfreq") == 0) {
     124      if (ia == narg-1) Usage(true);
     125      sscanf(arg[ia+1],"%d,%d,%lf",&bf_killfreq, &nharm_killfreq, &sigf_killfreq);
     126      fg_killfreq = true;
     127      ia++;
     128    }   
    105129    else if (strcmp(arg[ia],"-range") == 0) {
    106130      if (ia == narg-1) Usage(true);
    … …  
    127151  InitTim();
    128152
    129   cout << ">>>> tstrztoi: Infile= " << infile << " outFile="
     153  cout << ">>>> simtst: Infile= " << infile << " outFile="
    130154       << outfile << endl;
    131155  cout << ">>>> Window Size WTOI= " << wtoi << " WDEGLI= " << wdegli
    … …  
    150174
    151175    int w1 = wtoi;
    152     int w2 = (fg_f_filt) ? 2*wfft+w1 : w1;
     176    int w2 = (fg_f_filt) ? wfft+w1 : w1;
    153177
    154178    //    char * colname[5] = {"MJD", "UTC","boloMuV_11","magnFlux","pivot"};
    … …  
    234258    filt.addOutput("incopie", toideglioffcopie);
    235259
    236     Vector fcoef(wfft/2);
    237     fcoef = RegularSequence(1., -2./wfft);
     260    Vector fcoef(wfft/2+1);
     261    //    fcoef = RegularSequence(1., -2./wfft);
     262    fcoef = 1.;
     263    for(int kfk=0; kfk<nharm_killfreq; kfk++) {
     264      double cfreq = bf_killfreq*(kfk+1);
     265      double xfreq = 0.;
     266      for(int kfj=cfreq-3.*sigf_killfreq; kfj<= cfreq+3.*sigf_killfreq; kfj++) {
     267        if ( (kfj < 0) || (kfj >= fcoef.Size() ) )  continue;
     268        xfreq = (cfreq-kfj)/sigf_killfreq;
     269        fcoef(kfj) -= exp(-xfreq*xfreq);
     270      }
     271    }
    238272    cout << " Creating Fourier Filter ... " << endl;
    239273    SimpleFourierFilter sfft(fcoef);
    240274    sfft.KeepSpectra(outppfname, keepfft);
     275    sfft.ComputeMeanSpectra(true);
    241276    TOISeqBuffered * toifft = NULL;
    242277    TOISeqBuffered * toifiltcopie = NULL;
    … …  
    259294      w.addInput("fftout", toifft);
    260295    }
    261     else   w.addInput("out", toiout);
     296    else   w.addInput("filtout", toiout);
    262297    w.addInput("degli", toidegli1);
    263298    w.addInput("deglioff", toideglioffcopie);
Note: See TracChangeset for help on using the changeset viewer.