source: Sophya/trunk/AddOn/TAcq/brproc.cc@ 3778

Last change on this file since 3778 was 3778, checked in by ansari, 15 years ago

Ajout numero d'identification a BRBaseProcesseur, Reza 18/05/2010

File size: 41.8 KB
RevLine 
[3635]1#include "racquproc.h"
2
3#include <stdlib.h>
[3642]4#include <string.h>
[3635]5#include <unistd.h>
6#include <fstream>
7#include <signal.h>
8
9#include "pexceptions.h"
10#include "tvector.h"
[3646]11#include "ntuple.h"
[3647]12#include "datatable.h"
[3652]13#include "histos.h"
[3635]14#include "fioarr.h"
[3655]15#include "matharr.h"
[3635]16#include "timestamp.h"
17#include "ctimer.h"
18#include "fftpserver.h"
19#include "fftwserver.h"
20
21#include "FFTW/fftw3.h"
22
23
24#include "pciewrap.h"
25#include "brpaqu.h"
26#include "brproc.h"
27
[3683]28//---------------------------------------------------------------------
29// Classe de traitement - calcul de visibilite pour n fibres
30//---------------------------------------------------------------------
[3635]31
[3683]32/* --Methode-- */
[3776]33BRVisibilityCalculator::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())
[3692]36 // , dtfos_(outpath+"visdt.fits", Fits_Create), visdt_(dtfos_, 1024, true);
37{
[3776]38 DefineRank(1,0);
[3689]39
[3776]40 uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
41 chanum_.SetSize(maxnpairs);
[3692]42 sa_size_t k=0;
[3694]43 for(size_t i=0; i<2*memgr_.NbFibres(); i++) vpdata_[i]=NULL;
[3692]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 }
[3776]49 SelectPairs();
50
[3692]51 // visdt_.AddFloatColumn("mfc");
[3696]52 visdt_.AddDoubleColumn("mfc");
53 visdt_.AddDoubleColumn("mtt");
[3692]54 visdt_.AddIntegerColumn("jfreq");
55 visdt_.AddIntegerColumn("numch");
56 visdt_.AddFloatColumn("vre");
57 visdt_.AddFloatColumn("vim");
[3689]58
[3692]59 if (nmean_ < 1) nmean_=memgr_.NbPaquets();
60 if (nmean_ < 1) nmean_=1;
61
[3776]62 cout << " BRVisibilityCalculator::/Info nmean=" << nmean_ << endl;
63
[3692]64 totnbpaq_=0;
65 numfile_=0;
[3777]66 nb_flop_=0.;
[3692]67 moyfc_=moytt_=0.;
68
[3726]69 fgallfibok=NULL;
[3705]70 fgcktt_=false;
[3778]71 setNameId("viscalc", 0);
[3683]72}
73
74/* --Methode-- */
75BRVisibilityCalculator::~BRVisibilityCalculator()
76{
[3692]77 cout << " BRVisibilityCalculator - Visibility Datatable : " << endl;
78 cout << visdt_;
[3776]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);
[3692]87 po << visdt_;
[3776]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_;
[3705]101 }
102 }
[3683]103}
104
105/* --Methode-- */
[3776]106void BRVisibilityCalculator::DefineRank(uint_4 nbc, uint_4 cid)
107{
108 if ((nbc>4)||(cid>=nbc))
109 throw ParmError("BRVisibilityCalculator::DefineRank() NbCalc > 4 !");
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[4]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD};
120 SetMemZAction(mmzas[calcid_]);
[3778]121 setNameId("viscalc_grp", calcid_);
[3776]122 }
123 return ;
124}
125
126/* --Methode-- */
127uint_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-- */
141int 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-- */
156int 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-- */
[3724]176void BRVisibilityCalculator::run()
177{
[3776]178 cout << " BRVisibilityCalculator[" << calcid_ << "/" << nbcalc_
179 << "]::run() - Starting " << " NFibers=" << memgr_.NbFibres()
180 << " NChan=" << 2*memgr_.NbFibres() << " NPairs=" << nbpairs_ << " First:" << pairst_ << endl;
181
[3726]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 }
[3724]193 paralex_.SetParallelTask(*this);
194 cout << " BRVisibilityCalculator::run()/Info : starting ParallelExecutor with nThreads="
195 << paralex_.nThreads() << " ... " << endl;
196 paralex_.start();
[3726]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;
[3774]205 int mid = memgr_.FindMemZoneId(mmact_); // (MemZA_ProcA);
[3726]206 Byte* buffg = memgr_.GetMemZone(mid);
207 if (buffg == NULL) {
[3778]208 cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetMemZone(" << mid << ") -> NULL" << endl;
[3726]209 setRC(7); fgrun=false;
210 break;
211 }
212 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
213 fbuff_[fib] = memgr_.GetMemZone(mid,fib);
214 if (fbuff_[fib] == NULL) { // cela ne devrait pas arriver
[3778]215 cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetMemZone(" << mid << "," << fib << ") -> NULL" << endl;
[3726]216 setRC(9); fgrun=false;
217 break;
218 }
219 }
220
221 if (totnbpaq_%nmean_ == 0) {
222 if (totnbpaq_ > 0) {
223 moyfc_/=nmean_;
224 moytt_/=nmean_;
225 vismtx_.Info()["MeanFC"] = moyfc_;
226 vismtx_.Info()["MeanTT"] = moytt_;
227 vismtx_.Info()["NPAQSUM"] = nmean_;
228
229 // ATTENTION : Matrice visibilites non moyennee
[3776]230 char nfile[48];
231 if (nbcalc_==1)
232 sprintf(nfile,"vismtx%d.ppf",numfile_);
233 else
234 sprintf(nfile,"vismtx_%d_%d.ppf",(int)calcid_,numfile_);
[3726]235 string flnm=outpath_+nfile;
236 POutPersist po(flnm);
237 po << vismtx_;
[3778]238 cout << numfile_ << "-BRVisibilityCalculator[" << calcid_ << "]::run() NPaqProc="
[3726]239 << totnbpaq_ << " -> Visibility Matrix in " << flnm << endl;
240 FillVisibTable(moyfc_, moytt_);
241 numfile_++;
242 }
243 vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
244 moyfc_=moytt_=0.;
245 }
246
247 for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) { // boucle sur les paquets d'une zone
248 fgallfibok[jp]=fgokallfibers_=true;
249 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
250 vpaq_[fib].Set(fbuff_[fib]+jp*paqsz);
251 vfgok_[fib] = vpchk_[fib].Check(vpaq_[fib],curfc_[fib]);
252 if (!vfgok_[fib]) fgallfibok[jp]=fgokallfibers_=false;
253 }
254 if (fgokallfibers_) {
255 if (totprocnpaq_==0) {
256 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
257 fcfirst_[fib]=curfc_[fib];
258 ttfirst_[fib]=vpaq_[fib].TimeTag();
259 }
260 }
261 totprocnpaq_++;
262 moyfc_ += curfc_[0];
263 moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
[3776]264 if ((fgcktt_)&&(calcid_==0)) CheckTimeTag();
[3726]265 totnbpaq_++;
266 }
267 } // Fin de boucle sur les paquets
268
269 // Execution parallele pour calcul des visibilites par bandes de frequence
270 int rcpex=paralex_.execute();
[3778]271 if (rcpex!=0) cout << " BRVisibilityCalculator[" << calcid_ << "]::run() / Error Rc[paralex_.execute()]=" << rcpex << endl;
[3726]272
[3774]273 memgr_.FreeMemZone(mid, mmsta_); // (MemZS_ProcA);
[3726]274 } // Fin de boucle sur les zones a traiter
275 //------------------------------------
[3778]276 cout << " --------- END BRVisibilityCalculator[" << calcid_ << "]::run() , TotNbProcPaq=" << totprocnpaq_ << endl;
[3726]277 /*
278 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) vpchk_[fib].Print();
279 cout << " ------------------------------------ " << endl;
280 */
281 delete[] fgallfibok;
[3724]282 }
[3726]283 catch (std::exception& exc) {
[3778]284 cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched std::exception " << exc.what() << endl;
[3726]285 setRC(98);
286 return;
287 }
288 catch(...) {
[3778]289 cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched unknown ... exception " << endl;
[3726]290 setRC(99);
291 return;
292 }
293
[3724]294}
295
296/* --Methode-- */
[3683]297int BRVisibilityCalculator::Process()
298{
[3696]299
[3689]300 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
[3694]301 vpdata_[2*fib] = vpaq_[fib].Data1C();
302 vpdata_[2*fib+1] = vpaq_[fib].Data2C();
303 }
[3696]304
[3692]305 if (totnbpaq_%nmean_ == 0) {
306 if (totnbpaq_ > 0) {
[3701]307 moyfc_/=nmean_;
308 moytt_/=nmean_;
309 vismtx_.Info()["MeanFC"] = moyfc_;
310 vismtx_.Info()["MeanTT"] = moytt_;
311 vismtx_.Info()["NPAQSUM"] = nmean_;
312
[3692]313 // ATTENTION : Matrice visibilites non moyennee
[3776]314 char nfile[48];
315 if (nbcalc_==1)
316 sprintf(nfile,"vismtx%d.ppf",numfile_);
317 else
318 sprintf(nfile,"vismtx_%d_%d.ppf",(int)calcid_,numfile_);
[3692]319 string flnm=outpath_+nfile;
320 POutPersist po(flnm);
321 po << vismtx_;
322 cout << numfile_ << "-BRVisibilityCalculator::Process() NPaqProc="
323 << totnbpaq_ << " -> Visibility Matrix in " << flnm << endl;
324 FillVisibTable(moyfc_, moytt_);
325 numfile_++;
[3689]326 }
[3692]327 vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
328 moyfc_=moytt_=0.;
[3689]329 }
[3724]330
[3726]331 sa_size_t k=0;
332 for(size_t i=0; i<vpdata_.size(); i++) {
333 for(size_t j=i; j<vpdata_.size(); j++) {
[3776]334 size_t kpair=i*vpdata_.size()+j;
335 if (kpair<pairst_) continue;
336 if (kpair>=(pairst_+nbpairs_)) break;
[3726]337 TVector< complex<r_4> > vis = vismtx_.Row(k); k++;
338 for(sa_size_t f=1; f<vis.Size(); f++) {
339 vis(f) += complex<r_4>((r_4)vpdata_[i][f].realB(), (r_4)vpdata_[i][f].imagB()) *
340 complex<r_4>((r_4)vpdata_[j][f].realB(), -(r_4)vpdata_[j][f].imagB());
[3724]341 }
[3777]342 nb_flop_ += (7.*(r_8)vis.Size());
[3724]343 }
344 }
[3726]345
[3724]346 moyfc_ += curfc_[0];
347 moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
[3776]348 if ((fgcktt_)&&(calcid_==0)) CheckTimeTag();
[3724]349 totnbpaq_++;
350 return 0;
351}
352
353/* --Methode-- */
354int BRVisibilityCalculator::execute(int tid)
355{
[3726]356 vector<TwoByteComplex*> pvpdata(2*memgr_.NbFibres());
357 size_t paqsz=memgr_.PaqSize();
358 BRPaquet ppaq(paqsz);
359
[3724]360 sa_size_t fdelt = vismtx_.NCols()/nparthr_;
361 sa_size_t fdeb = tid*fdelt;
362 sa_size_t ffin = (tid+1)*fdelt;
[3726]363
[3724]364 if (fdeb<1) fdeb=1;
365 if ((ffin>vismtx_.NCols())||(tid==(nparthr_-1))) ffin=vismtx_.NCols();
[3726]366
367 for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) { // boucle sur les paquets d'une zone
368 if (!fgallfibok[jp]) continue;
369 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
370 ppaq.Set(fbuff_[fib]+jp*paqsz);
371 pvpdata[2*fib] = ppaq.Data1C();
372 pvpdata[2*fib+1] = ppaq.Data2C();
373 }
374 sa_size_t k=0;
375 for(size_t i=0; i<vpdata_.size(); i++) {
376 for(size_t j=i; j<vpdata_.size(); j++) {
[3776]377 size_t kpair=i*vpdata_.size()+j;
378 if (kpair<pairst_) continue;
379 if (kpair>=(pairst_+nbpairs_)) break;
[3726]380 TVector< complex<r_4> > vis = vismtx_.Row(k); k++;
381 for(sa_size_t f=fdeb; f<ffin; f++) {
382 vis(f) += complex<r_4>((r_4)pvpdata[i][f].realB(), (r_4)pvpdata[i][f].imagB()) *
383 complex<r_4>((r_4)pvpdata[j][f].realB(), -(r_4)pvpdata[j][f].imagB());
384 }
[3777]385 nb_flop_ += (7.*(r_8)(ffin-fdeb));
[3692]386 }
387 }
[3726]388
389 } // Fin de boucle sur les paquets
390
[3686]391 return 0;
[3683]392}
393
[3692]394/* --Methode-- */
395int BRVisibilityCalculator::FillVisibTable(double fcm, double ttm)
396{
[3776]397 double xnt[10];
398 xnt[0]=fcm; xnt[1]=ttm/1.25e8;
[3692]399
400 if (djf_<2) {
401 for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
402 for(sa_size_t jf=jf1_; jf<jf2_; jf++) {
[3776]403 xnt[2]=jf;
404 xnt[3]=chanum_(rv+pairst_);
405 xnt[4]=vismtx_(rv,jf).real()/(r_4)(nmean_);
406 xnt[5]=vismtx_(rv,jf).imag()/(r_4)(nmean_);
407 visdt_.AddRow(xnt);
[3692]408 }
409 }
410 }
411 else {
412 for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
413 for(sa_size_t jf=jf1_; jf<jf2_; jf+=djf_) {
414 r_4 moyreal=0.;
415 r_4 moyimag=0.;
[3705]416 sa_size_t jjfmx=jf+djf_;
417 if (jjfmx > vismtx_.NCols()) jjfmx=vismtx_.NCols();
418 for(sa_size_t jjf=jf; jjf<jjfmx; jjf++) {
[3696]419 moyreal+=vismtx_(rv,jjf).real();
420 moyimag+=vismtx_(rv,jjf).imag();
[3692]421 }
[3776]422 xnt[2]=jf+djf_/2;
423 xnt[3]=chanum_(rv+pairst_);
424 xnt[4]=moyreal/(r_4)(nmean_*djf_);
425 xnt[5]=moyimag/(r_4)(nmean_*djf_);
426 visdt_.AddRow(xnt);
[3692]427 }
428 }
429 }
430 return 0;
431}
432
[3705]433/* --Methode-- */
434int BRVisibilityCalculator::CheckTimeTag()
435{
436 if (totnbpaq_==0) {
437 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
438 vlasttt_[fib]=ttfirst_[fib];
439 if (ttmtx_.NCols()>0) {
440 fcmtx_(fib,totnbpaq_) = curfc_[fib];
441 ttmtx_(fib,totnbpaq_) = vlasttt_[fib];
442 }
443 }
444 return 0;
445 }
446 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
[3709]447 int_8 ld = (int_8)vpaq_[fib].TimeTag()-(int_8)vlasttt_[fib];
448 int_8 fd = (int_8)vpaq_[fib].TimeTag()-(int_8)ttfirst_[fib]-(int_8)vpaq_[0].TimeTag()+(int_8)ttfirst_[0];
[3708]449 /* if ( (ld < mindeltatt_) || (fd<-5) || (fd>5)) { vbadtt_[fib]++; vnsamett_[fib]++; }
[3705]450 else {
451 if (fd!=0) vnsamett_[fib]++;
452 }
[3708]453 */
454 if (ld < mindeltatt_) vbadtt_[fib]++;
[3709]455 else {
456 if (fd != 0) vnsamett_[fib]++;
[3710]457 if ((fd<-5)||(fd>5)) vndiff5tt_[fib]++;
[3709]458 }
[3705]459 vlasttt_[fib]=vpaq_[fib].TimeTag();
460 if (totnbpaq_<ttmtx_.NCols()) {
461 fcmtx_(fib,totnbpaq_) = curfc_[fib];
462 ttmtx_(fib,totnbpaq_) = vlasttt_[fib];
463 }
464 }
465 return 0;
466}
467
[3776]468//-------------------------------------------------------------------------------
469// Classe Groupe (ensemble) de Calculateur de Visibilites, tourant en parallele
470//-------------------------------------------------------------------------------
[3774]471
[3776]472/* --Methode-- */
473BRVisCalcGroup::BRVisCalcGroup(size_t nbcalc, RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean, size_t nthr)
[3777]474 : tm_(false)
[3776]475{
476 if ((nbcalc<1)||(nbcalc>4))
477 throw ParmError("BRVisCalcGroup::BRVisCalcGroup NbCalc > 4 !");
478 for(size_t i=0; i<nbcalc; i++) {
479 BRVisibilityCalculator * viscp=new BRVisibilityCalculator(memgr, outpath, nmean, nthr);
480 viscp->DefineRank(nbcalc, i);
481 viscalcp_.push_back(viscp);
482 }
483}
484/* --Methode-- */
485BRVisCalcGroup::~BRVisCalcGroup()
486{
487 for(size_t i=0; i<viscalcp_.size(); i++)
488 delete viscalcp_[i];
489}
490/* --Methode-- */
491int BRVisCalcGroup::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
492{
493 int rc=0;
494 for(size_t i=0; i<viscalcp_.size(); i++)
495 rc=viscalcp_[i]->SelectFreqBinning(freq1, freq2, nbfreq);
496 return rc;
497}
498/* --Methode-- */
499void BRVisCalcGroup::start()
500{
501 for(size_t i=0; i<viscalcp_.size(); i++)
502 viscalcp_[i]->start();
[3777]503 tm_.SplitQ();
[3776]504}
505/* --Methode-- */
506void BRVisCalcGroup::join()
507{
[3777]508 r_8 totflop=0.;
509 for(size_t i=0; i<viscalcp_.size(); i++) {
[3776]510 viscalcp_[i]->join();
[3777]511 totflop += viscalcp_[i]->TotNbFLOP();
512 }
513 tm_.SplitQ();
514 cout << " ----------------------------------------------------------" << endl;
515 cout << " BRVisCalcGroup::join() : Finished , Elaspsed time: " << tm_.PartialElapsedTimems()
516 << " ms (total:" << tm_.TotalElapsedTimems() << ")" << endl;
517 cout << " ... TotalMegaFLOP= " << totflop/(1024.e3) << " @ "
518 << totflop/(r_8)tm_.PartialElapsedTimems()/(1024) << " MFLOP/s" << endl;
519 cout << " ----------------------------------------------------------" << endl;
520 return;
[3776]521}
522
523
[3774]524//---------------------------------------------------------------------
525// Classe de traitement simple - calcul de spectres moyennes / voie
526//---------------------------------------------------------------------
527/* --Methode-- */
528BRMeanSpecCalculator::BRMeanSpecCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean)
529 : BRBaseProcessor(memgr), outpath_(outpath), nmean_(nmean)
530{
531 BRPaquet paq(memgr_.PaqSize());
532 mspecmtx_.SetSize(2*memgr_.NbFibres(), paq.DataSize()/4);
533 numfile_=0;
534 totnbpaq_=0;
535}
536
537/* --Methode-- */
538BRMeanSpecCalculator::~BRMeanSpecCalculator()
539{
540}
541
542
543/* --Methode-- */
544int BRMeanSpecCalculator::Process()
545{
546
547 if (totnbpaq_%nmean_ == 0) {
548 if (totnbpaq_ > 0) {
549 mspecmtx_.Info()["NPAQSUM"] = nmean_;
550 mspecmtx_ /= (double)nmean_;
551 char nfile[32];
552 sprintf(nfile,"mspecmtx%d.ppf",numfile_);
553 string flnm=outpath_+nfile;
554 POutPersist po(flnm);
555 po << mspecmtx_;
556 cout << numfile_ << "-BRMeanSpecCalculator::Process() NPaqProc="
557 << totnbpaq_ << " -> Mean spectra Matrix in " << flnm << endl;
558 numfile_++;
559 }
560 mspecmtx_ = (r_8)(0.);
561 }
562
563 sa_size_t k=0;
564 for(size_t i=0; i<(size_t)2*memgr_.NbFibres(); i++) {
565 TwoByteComplex* zp=vpaq_[i/2].Data1C();
566 if (i%2==1) zp=vpaq_[i/2].Data2C();
567 TVector< r_4 > spec = mspecmtx_.Row(k); k++;
568 for(sa_size_t f=1; f<spec.Size(); f++) {
569 spec(f) += zp[f].module2F();
570 }
571 }
572
573 totnbpaq_++;
574 return 0;
575}
576
577
[3635]578//---------------------------------------------------------------
579// Classe thread de traitement donnees ADC avec 2 voies par frame
580//---------------------------------------------------------------
581
[3649]582// Mutex pour eviter le plantage du a FFTW qui ne semble pas thread-safe
583static ZMutex* pmutfftw=NULL;
584
[3645]585/* --Methode-- */
[3683]586BRProcA2C::BRProcA2C(RAcqMemZoneMgr& mem, string& path, bool fgraw, uint_4 nmean,
[3656]587 uint_4 nmax, bool fghist, uint_4 nfsmap, bool fgnotrl, int card)
[3635]588 : memgr(mem)
589{
[3683]590 fgraw_ = fgraw;
[3635]591 nmax_ = nmax;
592 nmean_ = nmean;
[3683]593 if (fgraw_) cout << " BRProcA2C::BRProcA2C() - constructeur RAW data - NMean=" << nmean_ << endl;
594 else cout << " BRProcA2C::BRProcA2C() - constructeur FFT data - NMean=" << nmean_ << endl;
[3656]595 nfsmap_ = nfsmap;
[3635]596 stop_ = false;
597 path_ = path;
[3640]598 fgnotrl_ = fgnotrl;
[3652]599 fghist_ = fghist;
[3645]600 card_ = card;
[3649]601 if (pmutfftw==NULL) pmutfftw=new ZMutex;
[3635]602}
603
[3645]604/* --Methode-- */
[3683]605void BRProcA2C::Stop()
[3635]606{
607 stop_=true;
[3683]608 // cout <<" BRProcA2C::Stop ... > STOP " << endl;
[3635]609}
610
611
612static inline r_4 Zmod2(complex<r_4> z)
613{ return (z.real()*z.real()+z.imag()*z.imag()); }
614
[3645]615static inline string card2name_(int card)
616{
617 if (card==2) return " (Chan3,4) ";
618 else return " (Chan1,2) ";
619}
620/* --Methode-- */
[3683]621void BRProcA2C::run()
[3635]622{
623 setRC(1);
624 try {
[3683]625 Timer tm("BRProcA2C", false);
[3635]626 TimeStamp ts;
[3646]627 BRPaqChecker pcheck(!fgnotrl_); // Verification/comptage des paquets
[3640]628
629 size_t totnbytesout = 0;
630 size_t totnbytesproc = 0;
631
[3683]632 cout << " BRProcA2C::run() - Starting " << ts << " NMaxMemZones=" << nmax_
[3645]633 << " NMean=" << nmean_ << card2name_(card_) << endl;
[3683]634 cout << " BRProcA2C::run()... - Output Data Path: " << path_ << endl;
[3635]635 char fname[512];
636// sprintf(fname,"%s/proc.log",path_.c_str());
637// ofstream filog(fname);
[3683]638// filog << " BRProcA2C::run() - starting log file " << ts << endl;
[3635]639// filog << " ... NMaxMemZones=" << nmax_ << " NMean=" << nmean_ << " Step=" << step_ << endl;
640
[3647]641/*----DELETE NTuple
[3646]642 const char* nnames[8] = {"fcs","tts","s1","s2","s12","s12re","s12im","s12phi"};
643 NTuple nt(8, nnames);
644 double xnt[10];
645 uint_4 nmnt = 0;
646 double ms1,ms2,ms12,ms12re,ms12im,ms12phi;
[3647]647----*/
[3683]648// Time sample (raw data) /FFT coeff histograms
649 Histo* ph1=NULL;
650 Histo* ph2=NULL;
651 if (fghist_) {
652 if (fgraw_) {
653 ph1 = new Histo(-0.5, 255.5, 256);
654 ph2 = new Histo(-0.5, 255.5, 256);
655 }
656 else {
657 ph1 = new Histo(-128.5, 128.5, 257);
658 ph2 = new Histo(-128.5, 128.5, 257);
659 }
660 }
661
[3635]662// Initialisation pour calcul FFT
[3640]663 TVector< complex<r_4> > cfour1; // composant TF
[3635]664 uint_4 paqsz = memgr.PaqSize();
665 uint_4 procpaqsz = memgr.ProcPaqSize();
[3646]666
667
[3635]668 BRPaquet pq(NULL, NULL, paqsz);
669 TVector<r_4> vx(pq.DataSize()/2);
[3648]670 int szfour = pq.DataSize()/2/2+1;
671 cfour1.SetSize(szfour);
672/*
[3635]673 vx = (r_4)(0.);
674 FFTPackServer ffts;
[3640]675 ffts.FFTForward(vx, cfour1);
[3648]676 szfour = cfour1.Size();
677*/
678
[3656]679 bool fgtimfreq = false; // true->cartes temps<>frequences
680 if (nfsmap_>0) fgtimfreq=true;
681
[3640]682 TVector< complex<r_4> > cfour2(cfour1.Size());
[3635]683
[3640]684 TVector<r_4> spectreV1(cfour1.Size());
685 TVector<r_4> spectreV2(cfour1.Size());
[3655]686 TVector<r_4> moyspecV1(cfour1.Size()); // Moyenne des Spectres
687 TVector<r_4> moyspecV2(cfour1.Size());
688 TVector<r_4> sigspecV1(cfour1.Size()); // Sigma des Spectres
689 TVector<r_4> sigspecV2(cfour1.Size());
[3640]690 TVector< complex<r_4> > visiV12( cfour1.Size() );
[3635]691
[3656]692 TMatrix<r_4> timfreqV1, timfreqV2; // Cartes temps<>frequences
693 if (fgtimfreq) {
694 timfreqV1.SetSize(nmean_, spectreV1.Size()/nfsmap_);
695 timfreqV2.SetSize(nmean_, spectreV2.Size()/nfsmap_);
696 }
[3683]697 cout << " *DBG*BRProcA2C PaqSz=" << paqsz << " ProcPaqSize=" << procpaqsz
[3646]698 << " procpaqsz/2=" << procpaqsz/2 << " cfour1.Size()=" << cfour1.Size()
699 << " *8=" << cfour1.Size()*8 << endl;
[3635]700
[3649]701 pmutfftw->lock();
[3640]702 fftwf_plan plan1 = fftwf_plan_dft_r2c_1d(vx.Size(), vx.Data(),
703 (fftwf_complex*)cfour1.Data(), FFTW_ESTIMATE);
704 fftwf_plan plan2 = fftwf_plan_dft_r2c_1d(vx.Size(), vx.Data(),
705 (fftwf_complex*)cfour2.Data(), FFTW_ESTIMATE);
[3649]706 pmutfftw->unlock();
[3640]707
[3635]708 uint_4 ifile = 0;
[3655]709 uint_4 nzm = 0; // Nb de paquets moyennes pour le calcul de chaque spectre
710 uint_4 nmoyspec = 0; // Nb de spectres moyennes
[3647]711
712 uint_4 curfc=0;
713 uint_8 curtt=0;
714 uint_8 firsttt=0;
715 bool fgfirst=true;
[3658]716 double moysig[2]={0.,0.};
717 double sigsig[2]={0.,0.};
718 uint_8 nbsig[2]={0,0};
719
[3635]720 for (uint_4 kmz=0; kmz<nmax_; kmz++) {
721 if (stop_) break;
722 int mid = memgr.FindMemZoneId(MemZA_ProcA);
723 Byte* buff = memgr.GetMemZone(mid);
724 if (buff == NULL) {
[3683]725 cout << " BRProcA2C::run()/ERROR memgr.GetMemZone(" << mid << ") -> NULL" << endl;
[3640]726 break;
[3635]727 }
728 Byte* procbuff = memgr.GetProcMemZone(mid);
729 if (procbuff == NULL) {
[3683]730 cout << " BRProcA2C::run()/ERROR memgr.GetProcMemZone(" << mid << ") -> NULL" << endl;
[3640]731 break;
[3635]732 }
[3647]733//---- DELETE nmnt=0; ms1=ms2=ms12=ms12re=ms12im=ms12phi=0.;
[3645]734 for(uint_4 i=0; i<memgr.NbPaquets(); i++) {
[3640]735 BRPaquet paq(NULL, buff+i*paqsz, paqsz);
736 if (!pcheck.Check(paq)) continue; // on ne traite que les paquets OK
[3647]737 if (fgfirst) { firsttt=paq.TimeTag(); fgfirst=false; }
738 curfc=paq.FrameCounter();
739 curtt=paq.TimeTag()-firsttt;
[3635]740// Traitement voie 1
[3652]741 if (fghist_) {
742 for(sa_size_t j=0; j<vx.Size(); j++) {
[3683]743 r_4 vts=(fgraw_)?((r_4)(*(paq.Data1()+j))):((r_4)(*(paq.Data1S()+j)));
744 ph1->Add((r_8)vts);
[3658]745 moysig[0] += (double)vts;
746 sigsig[0] += ((double)vts)*((double)vts);
747 nbsig[0]++;
[3652]748 }
[3683]749 for(sa_size_t j=0; j<vx.Size(); j++) {
750 r_4 vts=(fgraw_)?((r_4)(*(paq.Data2()+j))):((r_4)(*(paq.Data2S()+j)));
751 ph2->Add((r_8)vts);
[3658]752 moysig[1] += (double)vts;
753 sigsig[1] += ((double)vts)*((double)vts);
754 nbsig[1]++;
[3652]755 }
756 }
[3683]757 if (fgraw_) {
758 for(sa_size_t j=0; j<vx.Size(); j++)
759 vx(j) = (r_4)(*(paq.Data1()+j))-127.5;
760 // fftwf_complex* coeff1 = (fftwf_complex*)(procbuff+i*procpaqsz);
761 fftwf_execute(plan1);
762 // Traitement voie 2
763 for(sa_size_t j=0; j<vx.Size(); j++)
764 vx(j) = (r_4)(*(paq.Data2()+j))-127.5;
765 fftwf_execute(plan2);
766 }
767 else {
768 for(sa_size_t j=1; j<cfour1.Size()-1; j++) {
769 cfour1(j) = complex<r_4>((r_4)paq.Data1C()[j].realB(), (r_4)paq.Data1C()[j].imagB());
770 cfour2(j) = complex<r_4>((r_4)paq.Data2C()[j].realB(), (r_4)paq.Data2C()[j].imagB());
771 }
772 cfour1(0) = complex<r_4>((r_4)paq.Data1C()[0].realB(), (r_4)0.);
773 cfour1(cfour1.Size()-1) = complex<r_4>((r_4)paq.Data1C()[0].imagB(), (r_4)0.);
774 cfour2(0) = complex<r_4>((r_4)paq.Data2C()[0].realB(), (r_4)0.);
775 cfour2(cfour2.Size()-1) = complex<r_4>((r_4)paq.Data2C()[0].imagB(), (r_4)0.);
776 }
777 for(sa_size_t j=0; j<spectreV1.Size(); j++)
778 spectreV1(j) += Zmod2(cfour1(j));
779 memcpy(procbuff+i*procpaqsz, cfour1.Data(), sizeof(complex<r_4>)*cfour1.Size());
780 if (fgtimfreq) { // Remplissage tableau temps-frequence
781 for(sa_size_t c=1; c<timfreqV1.NCols(); c++) {
782 for(sa_size_t j=c*nfsmap_; j<(c+1)*nfsmap_; j++)
783 timfreqV1(nzm, c) += Zmod2(cfour1(j));
784 }
785 }
[3635]786 for(sa_size_t j=0; j<spectreV2.Size(); j++)
[3640]787 spectreV2(j) += Zmod2(cfour2(j)); // Zmod2(zp2[j]);
788 memcpy(procbuff+i*procpaqsz+procpaqsz/2, cfour2.Data(), sizeof(complex<r_4>)*cfour2.Size());
[3656]789 if (fgtimfreq) { // Remplissage tableau temps-frequence
790 for(sa_size_t c=1; c<timfreqV2.NCols(); c++) {
791 for(sa_size_t j=c*nfsmap_; j<(c+1)*nfsmap_; j++)
792 timfreqV2(nzm,c) += Zmod2(cfour2(j));
793 }
794 }
[3635]795
796// Calcul correlation (visibilite V1 * V2)
[3640]797 for(sa_size_t j=0; j<visiV12.Size(); j++)
798 visiV12(j)+=cfour1(j)*conj(cfour2(j));
799// for(sa_size_t j=0; j<visiV12.Size(); j++) visiV12(j)+=zp1[j]*zp2[j];
[3647]800 if (nzm==0) {
801 spectreV1.Info()["StartFC"] = curfc;
802 spectreV2.Info()["StartFC"] = curfc;
803 visiV12.Info()["StartFC"] = curfc;
804 spectreV1.Info()["StartTT"] = curtt;
805 spectreV2.Info()["StartTT"] = curtt;
806 visiV12.Info()["StartTT"] = curtt;
807 }
[3640]808 nzm++;
[3647]809/*----DELETE
[3646]810 if (nmnt==0) { xnt[0]=paq.FrameCounter(); xnt[1]=paq.TimeTag(); }
811 for(sa_size_t j=2700; j<2800; j++) {
812 ms1 += Zmod2(cfour1(j)); ms2 += Zmod2(cfour2(j));
813 complex<r_4> zvis = cfour1(j)*conj(cfour2(j));
814 ms12 += Zmod2(zvis); ms12re += zvis.real(); ms12im += zvis.imag();
815 ms12phi+= atan2(zvis.imag(),zvis.real());
816 }
817 nmnt++;
[3647]818----*/
[3640]819 totnbytesproc += paq.DataSize();
820 totnbytesout += (2*sizeof(complex<r_4>)*cfour1.Size());
[3635]821
[3640]822 } // Fin de boucle sur les paquets d'une zone
[3647]823
824/*---- DELETE
[3646]825 if (nmnt>0) {
826 double fnorm = (double)nmnt*(2800-2700);
827 xnt[2] = ms1 /= fnorm;
828 xnt[3] = ms2 /= fnorm;
829 xnt[4] = ms12 /= fnorm;
830 xnt[5] = ms12re /= fnorm;
831 xnt[6] = ms12im /= fnorm;
832 xnt[7] = ms12phi /= fnorm;
833 nt.Fill(xnt);
834 }
[3647]835----*/
[3640]836 if ((nzm >= nmean_) || ((kmz==(nmax_-1))&&(nzm>1))) {
[3635]837 spectreV1 /= (r_4)(nzm);
838 spectreV2 /= (r_4)(nzm);
839
[3655]840 // pour le calcul des moyennes et sigmas de ces spectres
841 moyspecV1 += spectreV1;
842 moyspecV2 += spectreV2;
843 sigspecV1 += (spectreV1 && spectreV1);
844 sigspecV2 += (spectreV2 && spectreV2);
845 nmoyspec++;
846
[3635]847 visiV12 /= complex<r_4>((r_4)nzm, 0.);
848
849 spectreV1.Info()["NPaqMoy"] = nzm;
850 spectreV2.Info()["NPaqMoy"] = nzm;
851 visiV12.Info()["NPaqMoy"] = nzm;
[3647]852 spectreV1.Info()["EndFC"] = curfc;
853 spectreV2.Info()["EndFC"] = curfc;
854 visiV12.Info()["EndFC"] = curfc;
855 spectreV1.Info()["EndTT"] = curtt;
856 spectreV2.Info()["EndTT"] = curtt;
857 visiV12.Info()["EndTT"] = curtt;
[3652]858 {
[3635]859 sprintf(fname,"%s_%d.ppf",path_.c_str(),(int)ifile);
860 POutPersist po(fname);
[3645]861 string tag1="specV1";
862 string tag2="specV2";
863 string tag12="visiV12";
[3652]864 string tagh1="tshV1";
865 string tagh2="tshV2";
[3656]866 string tagtf1="timfreqV1";
867 string tagtf2="timfreqV2";
[3645]868 if (card_==2) {
869 tag1 = "specV3";
870 tag2 = "specV4";
[3652]871 tagh1 = "tshV1";
872 tagh2 = "tshV2";
[3645]873 tag12="visiV34";
[3656]874 tagtf1="timfreqV3";
875 tagtf2="timfreqV4";
[3645]876 }
877 po << PPFNameTag(tag1) << spectreV1;
878 po << PPFNameTag(tag2) << spectreV2;
879 po << PPFNameTag(tag12) << visiV12;
[3652]880 if (fghist_) {
[3683]881 po << PPFNameTag(tagh1) << (*ph1);
882 po << PPFNameTag(tagh2) << (*ph2);
[3658]883
884 double sspvmax[3] = {0.,0.,0.};
885 int_4 sspvmaxidx[3] = {-1,-1,-1};
886 for(int jji=1;jji<visiV12.Size()-1;jji++) {
887 r_4 zmv2 = Zmod2(visiV12(jji));
888 if (zmv2>sspvmax[2]) { sspvmax[2]=zmv2; sspvmaxidx[2]=jji; }
889 }
890 TVector<r_4>& sspv = spectreV1;
891 for(int ic=0; ic<2; ic++) {
892 if (ic==1) sspv = spectreV2;
893 for(int jji=1;jji<sspv.Size()-1;jji++)
894 if (sspv(jji)>sspvmax[ic]) { sspvmax[ic]=sspv(jji); sspvmaxidx[ic]=jji; }
895 if (nbsig[ic] < 1) { moysig[ic]=sigsig[ic]=-1.; }
896 else {
897 moysig[ic] /= (double)nbsig[ic];
898 sigsig[ic] /= (double)nbsig[ic];
899 sigsig[ic] -= (moysig[ic]*moysig[ic]);
900 sigsig[ic] = sqrt(sigsig[ic]);
901 cout << "===Voie " << ic << " Moy=" << moysig[ic] << " Sig=" << sigsig[ic]
902 << " MaxSpec Amp= " << sqrt(sspvmax[ic])/double(pq.DataSize()/2/2)
903 << " Pos=" << sspvmaxidx[ic] << " (NPts=" << nbsig[ic] << ")" << endl;
904 }
905 }
906 cout << "=== Voie1x2 MaxSpec Amp= " << sqrt(sqrt(sspvmax[2])/double(pq.DataSize()/2/2))
907 << " Pos=" << sspvmaxidx[2] << endl;
908 } // fin if (fghist_)
909
[3656]910 if (fgtimfreq) {
911 timfreqV1 /= (r_4)nzm;
912 timfreqV2 /= (r_4)nzm;
913 po << PPFNameTag(tagtf1) << timfreqV1;
914 po << PPFNameTag(tagtf2) << timfreqV2;
915 }
[3652]916 }
[3635]917 spectreV1 = (r_4)(0.);
918 spectreV2 = (r_4)(0.);
919 visiV12 = complex<r_4>(0., 0.);
[3652]920 if (fghist_) {
[3683]921 ph1->Zero();
922 ph2->Zero();
[3658]923 moysig[0]=moysig[1]=0.;
924 sigsig[0]=sigsig[1]=0.;
925 nbsig[0]=nbsig[1]=0;
[3652]926 }
[3656]927 if (fgtimfreq) {
928 timfreqV1 = (r_4)(0.);
929 timfreqV2 = (r_4)(0.);
930 }
[3635]931 nzm = 0; ifile++;
932// ts.SetNow();
933// filog << ts << " : proc file " << fname << endl;
[3683]934 cout << " BRProcA2C::run() created file " << fname << card2name_(card_) << endl;
[3640]935 }
[3635]936
937 memgr.FreeMemZone(mid, MemZS_ProcA);
[3640]938 } // Fin de boucle sur les zones a traiter
[3683]939 cout << " ------------ BRProcA2C::run() END " << card2name_(card_)
[3645]940 << " ------------ " << endl;
[3647]941/*---- DELETE
942 {
943 nt.Info()["FirstTT"]=firsttt;
[3646]944 cout << nt;
945 sprintf(fname,"%s_nt.ppf",path_.c_str());
946 POutPersist po(fname);
947 po << PPFNameTag("ntv12") << nt;
[3683]948 cout << " BRProcA2C::run() created NTuple file " << fname << card2name_(card_) << endl;
[3647]949 }
950---- */
[3655]951 if (nmoyspec>0) { // Calcul des moyennes et sigmas des spectres
952 r_4 fnms = nmoyspec;
953 moyspecV1 /= fnms;
954 moyspecV2 /= fnms;
955 sigspecV1 /= fnms;
956 sigspecV2 /= fnms;
957 sigspecV1 -= (moyspecV1 && moyspecV1);
958 sigspecV2 -= (moyspecV2 && moyspecV2);
959 sigspecV1 = Sqrt(sigspecV1);
960 sigspecV2 = Sqrt(sigspecV2);
961 TVector<r_4> rsbV1, rsbV2; // Rapport signal/bruit
962 moyspecV1.DivElt(sigspecV1, rsbV1, false, true);
963 moyspecV2.DivElt(sigspecV2, rsbV2, false, true);
964 sprintf(fname,"%s_ms.ppf",path_.c_str());
965 POutPersist po(fname);
966 po << PPFNameTag("moyspecV1") << moyspecV1;
967 po << PPFNameTag("moyspecV2") << moyspecV2;
968 po << PPFNameTag("sigspecV1") << sigspecV1;
969 po << PPFNameTag("sigspecV2") << sigspecV2;
970 po << PPFNameTag("rsbV1") << rsbV1;
971 po << PPFNameTag("rsbV2") << rsbV2;
[3683]972 cout << " BRProcA2C::run() created moysigspec file " << fname << card2name_(card_) << endl;
[3655]973 }
974
[3683]975 if (fghist_) {
976 delete ph1;
977 delete ph2;
978 }
[3640]979 ts.SetNow();
980 tm.SplitQ();
981 cout << " TotalProc= " << totnbytesproc/(1024*1024) << " MBytes, rate= "
982 << (double)(totnbytesproc)/1024./tm.PartialElapsedTimems() << " MB/s"
983 << " ProcDataOut=" << totnbytesout/(1024*1024) << " MB" << endl;
984 cout << pcheck;
[3683]985 cout << " BRProcA2C::run()/Timing: " << card2name_(card_) << endl;
[3640]986 tm.Print();
987 cout << " ---------------------------------------------------------- " << endl;
988
[3635]989 }
990 catch (PException& exc) {
[3683]991 cout << " BRProcA2C::run()/catched PException " << exc.Msg() << endl;
[3635]992 setRC(3);
993 return;
994 }
995 catch(...) {
[3683]996 cout << " BRProcA2C::run()/catched unknown ... exception " << endl;
[3635]997 setRC(4);
998 return;
999 }
1000 setRC(0);
1001 return;
1002}
1003
[3645]1004//---------------------------------------------------------------------
[3683]1005// Classe thread de traitement 2 x 2 voies/frames (Apres BRProcA2C)
[3645]1006//---------------------------------------------------------------------
[3635]1007
[3645]1008/* --Methode-- */
[3683]1009BRProcB4C::BRProcB4C(RAcqMemZoneMgr& mem1, RAcqMemZoneMgr& mem2, string& path,
1010 bool fgraw, uint_4 nmean, uint_4 nmax, bool fgnotrl)
[3645]1011 : memgr1(mem1), memgr2(mem2)
1012{
[3683]1013 fgraw_ = fgraw;
[3645]1014 nmax_ = nmax;
1015 nmean_ = nmean;
[3683]1016 if (fgraw_) cout << " BRProcB4C::BRProcB4C() - constructeur RAW data - NMean= " << nmean_ << endl;
1017 else cout << " BRProcB4C::BRProcB4C() - constructeur FFT data - NMean= " << nmean_ << endl;
[3645]1018 stop_ = false;
1019 path_ = path;
1020 fgnotrl_ = fgnotrl;
1021}
[3635]1022
[3645]1023/* --Methode-- */
[3683]1024void BRProcB4C::Stop()
[3645]1025{
1026 stop_=true;
[3683]1027 // cout <<" BRProcB4C::Stop ... > STOP " << endl;
[3645]1028}
[3635]1029
[3645]1030
1031/* --Methode-- */
[3683]1032void BRProcB4C::run()
[3645]1033{
1034 setRC(1);
1035 try {
[3683]1036 Timer tm("BRProcB4C", false);
[3645]1037 TimeStamp ts;
[3646]1038 BRPaqChecker pcheck1(!fgnotrl_); // Verification/comptage des paquets
1039 BRPaqChecker pcheck2(!fgnotrl_); // Verification/comptage des paquets
[3645]1040
1041 size_t totnbytesout = 0;
1042 size_t totnbytesproc = 0;
1043
[3683]1044 cout << " BRProcB4C::run() - Starting " << ts << " NMaxMemZones=" << nmax_
[3645]1045 << " NMean=" << nmean_ << endl;
[3683]1046 cout << " BRProcB4C::run()... - Output Data Path: " << path_ << endl;
[3645]1047
1048 uint_4 paqsz = memgr1.PaqSize();
1049 uint_4 procpaqsz = memgr1.ProcPaqSize();
1050 if ((paqsz != memgr2.PaqSize())||(procpaqsz!= memgr2.ProcPaqSize())) {
[3683]1051 cout << "BRProcB4C::run()/ERROR : different paquet size -> stop \n ...(PaqSz1="
[3645]1052 << paqsz << " Sz2=" << memgr2.PaqSize() << " ProcPaqSz1="
1053 << procpaqsz << " Sz2=" << memgr2.ProcPaqSize() << " )" << endl;
1054 setRC(9);
1055 return;
1056 }
1057
1058 TVector< complex<r_4> > cfour; // composant TF
[3648]1059 BRPaquet pq(NULL, NULL, paqsz);
[3646]1060/*
[3645]1061 TVector<r_4> vx(pq.DataSize()/2);
1062 vx = (r_4)(0.);
1063 FFTPackServer ffts;
1064 ffts.FFTForward(vx, cfour);
1065
1066 TVector< complex<r_4> > visiV13( cfour.Size() );
1067 TVector< complex<r_4> > visiV14( cfour.Size() );
1068 TVector< complex<r_4> > visiV23( cfour.Size() );
1069 TVector< complex<r_4> > visiV24( cfour.Size() );
[3646]1070*/
[3648]1071 int szfour = pq.DataSize()/2/2+1;
1072// int szfour = (paqsz-40)/2+1;
[3646]1073 TVector< complex<r_4> > visiV13( szfour );
1074 TVector< complex<r_4> > visiV14( szfour );
1075 TVector< complex<r_4> > visiV23( szfour );
1076 TVector< complex<r_4> > visiV24( szfour );
1077 // cout << " *DBG*AAAAA ---- Vectors OK" << endl;
[3683]1078 cout << " *DBG*BRProcB4C PaqSz=" << paqsz << " ProcPaqSize=" << procpaqsz
[3646]1079 << " procpaqsz/2=" << procpaqsz/2 << " cfour.Size()=" << szfour
1080 << " *8=" << szfour*8 << endl;
[3645]1081
[3647]1082 DataTable dt;
1083 dt.AddLongColumn("fc1");
[3651]1084 dt.AddLongColumn("tt1");
[3647]1085 dt.AddLongColumn("fc2");
1086 dt.AddLongColumn("tt2");
1087 DataTableRow dtr = dt.EmptyRow();
1088
[3645]1089 uint_4 nzm = 0;
1090 uint_4 totnoksfc = 0;
1091 uint_4 totnokpaq = 0;
1092 uint_4 totnpaq = 0;
1093 uint_4 ifile = 0;
[3647]1094
1095 uint_4 curfc=0;
1096 uint_8 curtt=0;
1097 uint_4 curfc2=0;
1098 uint_8 curtt2=0;
1099 uint_8 firsttt=0;
1100 uint_8 firsttt2=0;
1101 bool fgfirst=true;
[3645]1102 for (uint_4 kmz=0; kmz<nmax_; kmz++) {
1103 uint_4 noksfc = 0;
1104 uint_4 nokpaq = 0;
1105 if (stop_) break;
[3646]1106 // cout << " *DBG*BBBBB" << kmz << endl;
1107
[3645]1108 int mid1 = memgr1.FindMemZoneId(MemZA_ProcB);
1109 Byte* buff1 = memgr1.GetMemZone(mid1);
1110 if (buff1 == NULL) {
[3683]1111 cout << " BRProcB4C::run()/ERROR memgr.GetMemZone(" << mid1 << ") -> NULL" << endl;
[3645]1112 break;
1113 }
1114 Byte* procbuff1 = memgr1.GetProcMemZone(mid1);
1115 if (procbuff1 == NULL) {
[3683]1116 cout << " BRProcB4C::run()/ERROR memgr.GetProcMemZone(" << mid1 << ") -> NULL" << endl;
[3645]1117 break;
1118 }
1119 int mid2 = memgr2.FindMemZoneId(MemZA_ProcB);
1120 Byte* buff2 = memgr2.GetMemZone(mid2);
1121 if (buff1 == NULL) {
[3683]1122 cout << " BRProcB4C::run()/ERROR memgr.GetMemZone(" << mid2 << ") -> NULL" << endl;
[3645]1123 break;
1124 }
1125 Byte* procbuff2 = memgr2.GetProcMemZone(mid2);
1126 if (procbuff2 == NULL) {
[3683]1127 cout << " BRProcB4C::run()/ERROR memgr.GetProcMemZone(" << mid2 << ") -> NULL" << endl;
[3645]1128 break;
1129 }
1130 uint_4 i1,i2;
1131 i1=i2=0;
[3646]1132// cout << " *DBG*CCCCCC " << kmz << " memgr1.NbPaquets() =" << memgr1.NbPaquets() << endl;
[3645]1133 while((i1<memgr1.NbPaquets())&&(i2<memgr2.NbPaquets())) {
1134 BRPaquet paq1(NULL, buff1+i1*paqsz, paqsz);
1135 BRPaquet paq2(NULL, buff2+i2*paqsz, paqsz);
1136 totnpaq++;
1137// cout << " DBG["<<kmz<<"] i1,i2=" << i1 <<","<<i2<<" FC1,FC2=" <<paq1.FrameCounter()
1138//<<","<<paq2.FrameCounter()<<endl;
1139 // on ne traite que les paquets OK
1140 if (!pcheck1.Check(paq1)) { i1++; continue; }
1141 if (!pcheck2.Check(paq2)) { i2++; continue; }
1142 nokpaq++;
1143 if (paq1.FrameCounter()<paq2.FrameCounter()) { i1++; continue; }
1144 if (paq2.FrameCounter()<paq1.FrameCounter()) { i2++; continue; }
1145// cout << " DBG["<<kmz<<"]OKOK i1,i2=" << i1 <<","<<i2<<" FC1,FC2=" <<paq1.FrameCounter()
1146// <<","<<paq2.FrameCounter()<<endl;
1147
[3646]1148 if ((i1>=memgr1.NbPaquets())||(i2>=memgr1.NbPaquets())) {
1149 cout << " *BUG*BUG i1=" << i1 << " i2=" << i2 << endl;
1150 break;
1151 }
[3645]1152 // Les deux framecounters sont identiques ...
1153 noksfc++;
[3647]1154 curfc=paq1.FrameCounter();
1155 curfc2=paq2.FrameCounter();
1156 if (fgfirst) {
[3651]1157 firsttt=paq1.TimeTag(); firsttt2=paq2.TimeTag();
[3683]1158 cout << " BRProcB4C()/Info First FC="<<curfc<<" , "<<curfc2<<" -> TT="
[3647]1159 << firsttt<<" , "<<firsttt2 <<endl;
1160 fgfirst=false;
1161 }
1162 curtt=paq1.TimeTag()-firsttt;
1163 curtt2=paq2.TimeTag()-firsttt2;
1164 dtr[0]=curfc; dtr[1]=curtt;
1165 dtr[2]=curfc2; dtr[3]=curtt2;
1166 dt.AddRow(dtr);
1167
[3645]1168 complex<r_4>* zp1 = (complex<r_4>*)(procbuff1+i1*procpaqsz);
1169 complex<r_4>* zp2 = (complex<r_4>*)(procbuff1+i1*procpaqsz+procpaqsz/2);
1170 complex<r_4>* zp3 = (complex<r_4>*)(procbuff2+i2*procpaqsz);
1171 complex<r_4>* zp4 = (complex<r_4>*)(procbuff2+i2*procpaqsz+procpaqsz/2);
1172 for(sa_size_t j=0; j<visiV13.Size(); j++) {
1173 visiV13(j)+=zp1[j]*conj(zp3[j]);
1174 visiV14(j)+=zp1[j]*conj(zp4[j]);
1175 visiV23(j)+=zp2[j]*conj(zp3[j]);
1176 visiV24(j)+=zp2[j]*conj(zp4[j]);
1177 }
[3647]1178 if (nzm==0) {
1179 visiV13.Info()["StartFC"] = curfc;
1180 visiV14.Info()["StartFC"] = curfc;
1181 visiV23.Info()["StartFC"] = curfc;
1182 visiV24.Info()["StartFC"] = curfc;
1183 visiV13.Info()["StartTT"] = curtt;
1184 visiV14.Info()["StartTT"] = curtt;
1185 visiV23.Info()["StartTT"] = curtt;
1186 visiV24.Info()["StartTT"] = curtt;
1187 }
[3645]1188 nzm++; i1++; i2++;
1189 totnbytesproc += 2*paq1.DataSize();
1190 } // Fin de boucle sur les paquets d'une zone
[3646]1191 memgr1.FreeMemZone(mid1, MemZS_ProcB);
1192 memgr2.FreeMemZone(mid2, MemZS_ProcB);
1193
[3645]1194 if ((nzm >= nmean_) || ((kmz==(nmax_-1))&&(nzm>1))) {
1195 visiV13 /= complex<r_4>((r_4)nzm, 0.);
1196 visiV14 /= complex<r_4>((r_4)nzm, 0.);
1197 visiV23 /= complex<r_4>((r_4)nzm, 0.);
1198 visiV24 /= complex<r_4>((r_4)nzm, 0.);
1199 visiV13.Info()["NPaqMoy"] = nzm;
1200 visiV14.Info()["NPaqMoy"] = nzm;
1201 visiV23.Info()["NPaqMoy"] = nzm;
1202 visiV24.Info()["NPaqMoy"] = nzm;
[3647]1203 visiV13.Info()["EndFC"] = curfc;
1204 visiV14.Info()["EndFC"] = curfc;
1205 visiV23.Info()["EndFC"] = curfc;
1206 visiV24.Info()["EndFC"] = curfc;
1207 visiV13.Info()["EndTT"] = curtt;
1208 visiV14.Info()["EndTT"] = curtt;
1209 visiV23.Info()["EndTT"] = curtt;
1210 visiV24.Info()["EndTT"] = curtt;
[3645]1211 char fname[512];
1212 {
1213 sprintf(fname,"%s_%d.ppf",path_.c_str(),(int)ifile);
1214 POutPersist po(fname);
1215 po << PPFNameTag("visiV13") << visiV13;
1216 po << PPFNameTag("visiV14") << visiV14;
1217 po << PPFNameTag("visiV23") << visiV23;
1218 po << PPFNameTag("visiV24") << visiV24;
1219 }
1220 visiV13 = complex<r_4>(0., 0.);
1221 visiV14 = complex<r_4>(0., 0.);
1222 visiV23 = complex<r_4>(0., 0.);
1223 visiV24 = complex<r_4>(0., 0.);
[3646]1224 nzm = 0; ifile++;
[3645]1225// ts.SetNow();
1226// filog << ts << " : proc file " << fname << endl;
[3683]1227 cout << " BRProcB4C::run() created file " << fname << endl;
[3645]1228 }
1229 double okfrac = (nokpaq>1)?((double)noksfc/(double)nokpaq*100.):0.;
[3683]1230 cout << "BRProcB4C ["<<kmz<<"] NOKPaq=" << nokpaq << " NSameFC=" << noksfc
[3645]1231 << " (" << okfrac << " %)" << endl;
1232 totnokpaq += nokpaq;
1233 totnoksfc += noksfc;
1234 } // Fin de boucle sur les zones a traiter
[3683]1235 cout << " ------------------ BRProcB4C::run() END ----------------- " << endl;
[3647]1236 {
1237 dt.Info()["FirstTT1"]=firsttt;
1238 dt.Info()["FirstTT2"]=firsttt2;
1239 cout << dt;
1240 char fname[512];
1241 sprintf(fname,"%s_fctt.ppf",path_.c_str());
1242 POutPersist po(fname);
1243 po << PPFNameTag("ttfc") << dt;
[3683]1244 cout << " BRProcB4C::run() created TimeTag/FrameCounter file " << fname << endl;
[3647]1245 }
[3645]1246 ts.SetNow();
1247 tm.SplitQ();
1248 cout << " TotalProc= " << totnbytesproc/(1024*1024) << " MBytes, rate= "
1249 << (double)(totnbytesproc)/1024./tm.PartialElapsedTimems() << " MB/s" << endl;
1250 double totokfrac = (totnokpaq>1)?((double)totnoksfc/(double)totnokpaq*100.):0.;
1251 cout << " NOkPaq1,2=" << totnokpaq << " /TotNPaq=" << totnpaq << " TotNSameFC="
1252 << totnoksfc << " (" << totokfrac << " %)" << endl;
1253// cout << pcheck1;
1254// cout << pcheck2;
[3683]1255 cout << " BRProcB4C::run()/Timing: \n";
[3645]1256 tm.Print();
1257 cout << " ---------------------------------------------------------- " << endl;
1258}
1259 catch (PException& exc) {
[3683]1260 cout << " BRProcB4C::run()/catched PException " << exc.Msg() << endl;
[3645]1261 setRC(3);
1262 return;
1263 }
1264 catch(...) {
[3683]1265 cout << " BRProcB4C::run()/catched unknown ... exception " << endl;
[3645]1266 setRC(4);
1267 return;
1268 }
1269 setRC(0);
1270 return;
1271}
1272
1273
Note: See TracBrowser for help on using the repository browser.