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

Last change on this file since 4002 was 3993, checked in by ansari, 14 years ago

Ajout possibilite de faire un DataTable (NTuple) avec la puissance spectrale dans plusieurs bandes en fonction du temps, Reza 13/05/2011
File size: 26.1 KB
Line 
1//----------------------------------------------------------------
2// Projet BAORadio - (C) LAL/IRFU 2008-2010
3// Classes de threads de calcul de visibilites pour BAORadio
4//----------------------------------------------------------------
5
6#include <stdlib.h>
7#include <string.h>
8#include <unistd.h>
9#include <fstream>
10#include <signal.h>
11
12#include "pexceptions.h"
13#include "fioarr.h"
14#include "fitsarrhand.h"
15#include "matharr.h"
16#include "arrctcast.h"
17#include "timestamp.h"
18#include "ctimer.h"
19
20#include "brviscalc.h"
21
22
23//---------------------------------------------------------------------
24// Classe de traitement - calcul de visibilite pour n fibres
25//---------------------------------------------------------------------
26
27/* --Methode-- */
28BRVisibilityCalculator::BRVisibilityCalculator(RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean, size_t nthr)
29 : BRBaseProcessor(memgr), paralex_(*this, nthr), nparthr_(nthr),
30 outpath_(outpath), nmean_(nmean), nbcalc_(1), calcid_(0),
31 vpdata_(2*memgr.NbFibres()), vpdatar_(2*memgr.NbFibres())
32 // , dtfos_(outpath+"visdt.fits", Fits_Create), visdt_(dtfos_, 1024, true);
33{
34 SetFFTData();
35 DefineRank(1,0);
36 SetPPFOutput();
37
38 uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
39 chanids_.SetSize(2*memgr_.NbFibres());
40 chanpairnumall_.SetSize(maxnpairs);
41 chanpairsall_.SetSize(maxnpairs,2);
42 for(size_t i=0; i<2*memgr_.NbFibres(); i++) {
43 vpdata_[i]=NULL; vpdatar_[i]=NULL;
44 }
45 // Pas necessaire - appele par DefineRank- SelectPairs();
46
47 // visdt_.AddFloatColumn("mfc");
48 visdt_.AddDoubleColumn("mfc");
49 visdt_.AddDoubleColumn("mtt");
50 visdt_.AddIntegerColumn("jfreq");
51 visdt_.AddIntegerColumn("numch");
52 visdt_.AddFloatColumn("vre");
53 visdt_.AddFloatColumn("vim");
54 ActivateVisDTable();
55
56 if (nmean_ < 1) nmean_=memgr_.NbPaquets();
57 if (nmean_ < 1) nmean_=1;
58
59 cout << " BRVisibilityCalculator::/Info nmean=" << nmean_ << endl;
60
61 totnbpaq_=0;
62 numfile_=0;
63 nb_flop_=0.;
64 moyfc_=moytt_=0.;
65
66 fgallfibok=NULL;
67 fgcktt_=false;
68 setNameId("viscalc", 0);
69}
70
71/* --Methode-- */
72BRVisibilityCalculator::~BRVisibilityCalculator()
73{
74 if (totnbpaq_<1) return;
75 if (fgvisdt_) {
76 cout << " ~BRVisibilityCalculator - Visibility Datatable : " << endl;
77 cout << visdt_;
78
79 string filename;
80 filename = outpath_+"visdt.ppf";
81 if (nbcalc_>1) {
82 char sbuff[32];
83 sprintf(sbuff,"visdt_%d.ppf",(int)calcid_);
84 filename = outpath_+sbuff;
85 }
86 POutPersist po(filename);
87 po << visdt_;
88 }
89 if (calcid_ == 0) {
90 if (fgcktt_) {
91 cout << " BRVisibilityCalculator - Check TimeTag done: TotNPaqProc= " << totnbpaq_ << endl;
92 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
93 cout << " BRTTCheck-Fiber[" << fib << "] NBadTT=" << vbadtt_[fib] << " NDiffTT>5="
94 << vndiff5tt_[fib] << " NotSameTT=" << vnsamett_[fib] << endl;
95 }
96 POutPersist pott(outpath_+"ttfcmtx.ppf");
97 pott << PPFNameTag("FC") << fcmtx_;
98 pott << PPFNameTag("TT") << ttmtx_;
99 }
100 }
101}
102
103/* --Methode-- */
104void BRVisibilityCalculator::DefineRank(uint_4 nbc, uint_4 cid, uint_4 pair1, uint_4 nbpairs, bool fgpimp)
105{
106 if ((nbc>6)||(cid>=nbc))
107 throw ParmError("BRVisibilityCalculator::DefineRank() NbCalc > 6 !");
108 uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
109 if ((pair1>=maxnpairs)||(nbpairs<1))
110 throw ParmError("BRVisibilityCalculator::DefineRank() (pair1>=maxnpairs)||(nbpairs<1) !");
111 maxnpairs-=pair1;
112 if (maxnpairs>nbpairs) maxnpairs=nbpairs;
113
114 nbcalc_=nbc;
115 calcid_=cid;
116 if (nbcalc_>1) {
117 uint_4 npairspth=maxnpairs/nbcalc_; // nb de paires par thread
118 if (calcid_==(nbcalc_-1))
119 SelectPairs(calcid_*npairspth+pair1, maxnpairs-calcid_*npairspth, fgpimp);
120 else
121 SelectPairs(calcid_*npairspth+pair1, npairspth, fgpimp);
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};
124 SetMemZAction(mmzas[calcid_]);
125 setNameId("viscalc_grp", calcid_);
126 }
127 else SelectPairs(pair1,maxnpairs, fgpimp);
128 return ;
129}
130
131/* --Methode-- */
132uint_4 BRVisibilityCalculator::SelectPairs(uint_4 pair1, uint_4 nbpairs, bool fgpimp)
133{
134 BRPaquet paq(memgr_.PaqSize());
135 uint_4 maxnpairs = (2*memgr_.NbFibres()+1)*memgr_.NbFibres();
136
137 if (pair1 >= maxnpairs) pair1=maxnpairs-1;
138 if (nbpairs > maxnpairs-pair1) nbpairs=maxnpairs-pair1;
139
140 pairst_=pair1;
141 nbpairs_=nbpairs;
142 fgpimp_=fgpimp;
143
144 uint_4 nbvisicomp=0;
145 sa_size_t kpair=0;
146 for(size_t i=0; i<vpdata_.size(); i++) {
147 for(size_t j=i; j<vpdata_.size(); j++) {
148 kpair++;
149 if (kpair<(pairst_+1)) continue;
150 if (kpair>=(pairst_+nbpairs_+1)) break;
151 if (fgpimp_&&(i!=j)&&((i+j)%2==0)) continue; // calcul des visib avec numero pair-impair + autocorrel
152 nbvisicomp++;
153 }
154 }
155
156 vismtx_.SetSize(nbvisicomp, paq.DataSize()/4);
157
158 chanpairnum_.SetSize(nbvisicomp);
159 chanpairs_.SetSize(nbvisicomp,2);
160
161 return nbvisicomp;
162}
163
164/* --Methode-- */
165int BRVisibilityCalculator::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
166{
167 if (freq2<freq1) {
168 freq1=0; freq2=vismtx_.NCols()-1; nbfreq=1;
169 }
170 jf1_=freq1; jf2_=freq2;
171 if ((jf1_<1)||(jf1_>=vismtx_.NCols())) jf1_=1;
172 if ((jf2_<1)||(jf2_>=vismtx_.NCols())||(jf2_<jf1_)) jf2_=vismtx_.NCols()-1;
173 if (nbfreq<1) nbfreq=(jf2_-jf1_);
174 djf_=(jf2_-jf1_)/nbfreq;
175 if (djf_<1) djf_=1;
176 cout << " BRVisibilityCalculator::SelectFreqBinning/Info JF1=" << jf1_
177 << " JF2=" << jf2_ << " DJF=" << djf_ << endl;
178
179}
180
181/* --Methode-- */
182void BRVisibilityCalculator::UpdateChanIds()
183{
184
185 for(size_t i=0; i<memgr_.NbFibres(); i++) {
186 chanids_(2*i)=memgr_.FiberId(i)*2-1;
187 chanids_(2*i+1)=memgr_.FiberId(i)*2;
188 }
189 sa_size_t k=0; // numero de ligne dans la matrice des visibilites
190 for(size_t i=0; i<vpdata_.size(); i++) {
191 for(size_t j=i; j<vpdata_.size(); j++) {
192 chanpairnumall_(k)=chanids_(i)*CHANPAIRCONVFAC+chanids_(j);
193 chanpairsall_(k,0)=chanids_(i); chanpairsall_(k,1)=chanids_(j); k++;
194 }
195 }
196 sa_size_t kpair=0;
197 k=0; // numero de ligne dans la matrice des visibilites
198 for(size_t i=0; i<vpdata_.size(); i++) {
199 for(size_t j=i; j<vpdata_.size(); j++) {
200 kpair++;
201 if (kpair<(pairst_+1)) continue;
202 if (kpair>=(pairst_+nbpairs_+1)) break;
203 if (fgpimp_&&(i!=j)&&((i+j)%2==0)) continue; // calcul des visib avec numero pair-impair + autocorrel
204 chanpairnum_(k)=chanids_(i)*CHANPAIRCONVFAC+chanids_(j);
205 chanpairs_(k,0)=chanids_(i); chanpairs_(k,1)=chanids_(j); k++;
206 }
207 }
208
209 string filename;
210 filename = (outpath_+"chanum.")+OutFileExtension();
211 if (nbcalc_>1) {
212 char sbuff[32];
213 sprintf(sbuff,"chanum_%d.%s",(int)calcid_,OutFileExtension());
214 filename = outpath_+sbuff;
215 }
216 if (fgfitsout_) { // Ecriture au format FITS
217 FitsInOutFile foc(filename, FitsInOutFile::Fits_Create);
218 foc.SetNextExtensionName("chanids");
219 foc << chanids_;
220 foc.SetNextExtensionName("chanpairs");
221 foc << chanpairs_;
222 foc.SetNextExtensionName("chanpairnum");
223 foc << chanpairnum_;
224 foc.SetNextExtensionName("chanpairsall");
225 foc << chanpairsall_;
226 foc.SetNextExtensionName("chanpairnumall");
227 foc << chanpairnumall_;
228 }
229 else { // Format PPF
230 POutPersist poc(filename);
231 poc << PPFNameTag("chanids") << chanids_;
232 poc << PPFNameTag("chanpairs") << chanpairs_;
233 poc << PPFNameTag("chanpairnum") << chanpairnum_;
234 poc << PPFNameTag("chanpairsall") << chanpairsall_;
235 poc << PPFNameTag("chanpairnumall") << chanpairnumall_;
236 }
237 cout << "BRVisibilityCalculator[" << calcid_ << "]::UpdateChanIds() Channel Ids/Pairs saved to file "
238 << filename << endl;
239 cout << " ... NbVisib=NbChanPairs=" << chanpairs_.NRows() << " ChannelPairs= " ;
240 for(sa_size_t ir=0; ir<chanpairs_.NRows(); ir++) {
241 if (ir%10==0) cout << endl;
242 cout << "(" << chanpairs_(ir,0) << "," << chanpairs_(ir,1) << ") ";
243 }
244 cout << endl;
245 return;
246}
247
248/* --Methode-- */
249int BRVisibilityCalculator::ActivateTimeTagCheck(uint_8 maxnpaq)
250{
251 mindeltatt_=memgr_.PaqSize()/2;
252 if (mindeltatt_<1) mindeltatt_=1;
253 fcmtx_.SetSize(memgr_.NbFibres(), maxnpaq);
254 ttmtx_.SetSize(memgr_.NbFibres(), maxnpaq);
255 vlasttt_.resize(memgr_.NbFibres(), 0);
256 vbadtt_.resize(memgr_.NbFibres(), 0);
257 vnsamett_.resize(memgr_.NbFibres(), 0);
258 vndiff5tt_.resize(memgr_.NbFibres(), 0);
259
260 fgcktt_=true;
261 cout << " BRVisibilityCalculator::ActivateTimeTagCheck() - TT/Fc matrix NCols=" << maxnpaq
262 << " MinDeltaTT=" << mindeltatt_ << endl;
263
264 return 0;
265}
266
267
268/* --Methode-- */
269void BRVisibilityCalculator::run()
270{
271 cout << " BRVisibilityCalculator[" << calcid_ << "/" << nbcalc_
272 << "]::run() - Starting " << " NFibers=" << memgr_.NbFibres()
273 << " NChan=" << 2*memgr_.NbFibres() << " NPairs=" << nbpairs_ << " First:" << pairst_ << endl;
274
275 if (nparthr_ < 2) return BRBaseProcessor::run();
276 // Execution multithread parallele
277 setRC(1);
278 int rc=0;
279 try {
280 if ((nmean_%memgr_.NbPaquets())!=0) {
281 uint_4 mnmean = (nmean_/memgr_.NbPaquets()+1)*memgr_.NbPaquets();
282 cout << " BRVisibilityCalculator::run()/Info changing nmean=" << nmean_ << " to multiple of"
283 << " memgr_.NbPaquets() -> " << mnmean << endl;
284 nmean_=mnmean;
285 }
286 paralex_.SetParallelTask(*this);
287 cout << " BRVisibilityCalculator::run()/Info : starting ParallelExecutor with nThreads="
288 << paralex_.nThreads() << " ... " << endl;
289 paralex_.start();
290
291 fgallfibok = new bool[memgr_.NbPaquets()];
292
293 size_t paqsz=memgr_.PaqSize();
294 size_t procpaqsz=memgr_.ProcPaqSize();
295 bool fgrun=true;
296 while (fgrun) {
297 if (stop_) break;
298 if (memgr_.GetRunState() == MemZR_Stopped) break;
299 int mid = memgr_.FindMemZoneId(mmact_); // (MemZA_ProcA);
300 // Byte* buffg = memgr_.GetMemZone(mid);
301 // if (buffg == NULL) {
302 if (mid < 0) {
303 cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR FindMemZoneId("
304 << (int)mmact_ << ") ->" << mid << ") -> NULL" << endl;
305 setRC(7); fgrun=false;
306 break;
307 }
308 cts_=memgr_.GetAuxData(mid)->FillTime(); // get associated date/time (DATEOBS)
309 if (totnbpaq_==0) UpdateChanIds(); // Appele ici pour etre sur que le thread de remplissage a mis l'info a jour.
310
311 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
312 fbuff_[fib] = memgr_.GetMemZone(mid,fib);
313 if (fbuff_[fib] == NULL) { // cela ne devrait pas arriver
314 cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetMemZone(" << mid << "," << fib << ") -> NULL" << endl;
315 setRC(9); fgrun=false;
316 break;
317 }
318 if ((procpaqsz>0)&&((fprocbuff_[fib]=memgr_.GetProcMemZone(mid,fib))==NULL)) { // cela ne devrait pas arriver non plus
319 cout << "BRVisibilityCalculator[" << calcid_ << "]::run()/ERROR memgr.GetProcMemZone("
320 << mid << "," << fib << ") -> NULL" << endl;
321 setRC(9); fgrun=false;
322 break;
323 }
324 }
325
326 if (totnbpaq_%nmean_ == 0) {
327 if (totnbpaq_ > 0) {
328 moyfc_/=nmean_;
329 moytt_/=nmean_;
330 UpdateVisMtxInfo(); // add/update keywords in the Info DVList
331 // ATTENTION : Matrice visibilites non moyennee
332 char nfile[48];
333 if (nbcalc_==1)
334 sprintf(nfile,"vismtx%d.%s",numfile_,OutFileExtension());
335 else
336 sprintf(nfile,"vismtx_%d_%d.%s",(int)calcid_,numfile_,OutFileExtension());
337 string flnm=outpath_+nfile;
338 if (fgfitsout_) { // Ecriture au format FITS
339 FitsInOutFile fo(flnm, FitsInOutFile::Fits_Create);
340 TArray<r_4> arvismtx = ArrCastC2R(vismtx_);
341 arvismtx.Info()=vismtx_.Info();
342 fo << arvismtx;
343 }
344 else { // Format PPF
345 POutPersist po(flnm);
346 po << vismtx_;
347 }
348 if ((prtlev_>0)&&(numfile_%prtmodulo_==0)) {
349 cout << numfile_ << "-BRVisCalc[" << calcid_ << "/" << nbcalc_ << "]::run() NPaqProc="
350 << totnbpaq_ << " TotMegaFLOP=" << (uint_8)TotNbMegaFLOP() << " -> VisibMtx in " << flnm << endl;
351 }
352 if (fgvisdt_) FillVisibTable(moyfc_, moytt_);
353 numfile_++;
354 }
355 vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
356 moyfc_=moytt_=0.;
357 // first_tmstamp_.SetNow(); // Current date and time
358 first_tmstamp_=cts_; // Current date and time
359 }
360
361 for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) { // boucle sur les paquets d'une zone
362 fgallfibok[jp]=fgokallfibers_=true;
363 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
364 vpaq_[fib].Set(fbuff_[fib]+jp*paqsz);
365 vfgok_[fib] = vpchk_[fib].Check(vpaq_[fib],curfc_[fib]);
366 if (!vfgok_[fib]) fgallfibok[jp]=fgokallfibers_=false;
367 //Pas utile if (procpaqsz>0) vprocpaq_[fib] = fprocbuff_[fib]+jp*procpaqsz;
368 }
369 if (fgokallfibers_) {
370 if (totprocnpaq_==0) {
371 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
372 fcfirst_[fib]=curfc_[fib];
373 ttfirst_[fib]=vpaq_[fib].TimeTag();
374 }
375 }
376 totprocnpaq_++;
377 moyfc_ += curfc_[0];
378 moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
379 if ((fgcktt_)&&(calcid_==0)) CheckTimeTag();
380 if (totnbpaq_%nmean_ == 0) {
381 first_fc_=curfc_[0];
382 first_tt_= (vpaq_[0].TimeTag()-ttfirst_[0]);
383 }
384 totnbpaq_++;
385 }
386 } // Fin de boucle sur les paquets
387
388 // Execution parallele pour calcul des visibilites par bandes de frequence
389 int rcpex=paralex_.execute();
390 if (rcpex!=0) cout << " BRVisibilityCalculator[" << calcid_ << "]::run() / Error Rc[paralex_.execute()]=" << rcpex << endl;
391
392 memgr_.FreeMemZone(mid, mmsta_); // (MemZS_ProcA);
393 } // Fin de boucle sur les zones a traiter
394 //------------------------------------
395 cout << " --------- END BRVisibilityCalculator[" << calcid_ << "]::run() , TotNbProcPaq=" << totprocnpaq_ << endl;
396 /*
397 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) vpchk_[fib].Print();
398 cout << " ------------------------------------ " << endl;
399 */
400 delete[] fgallfibok;
401 }
402 catch (std::exception& exc) {
403 cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched std::exception " << exc.what() << endl;
404 setRC(98);
405 return;
406 }
407 catch(...) {
408 cout << " BRVisibilityCalculator[" << calcid_ << "]::run()/catched unknown ... exception " << endl;
409 setRC(99);
410 return;
411 }
412
413}
414
415
416/* --Methode-- */
417int BRVisibilityCalculator::Process()
418{
419 if (totnbpaq_==0) UpdateChanIds(); // Appele ici pour etre sur que le thread de remplissage a mis l'info a jour.
420 if (fgdataraw_) { // Donnees firmware RAW apres TF soft
421 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
422 vpdatar_[2*fib] = reinterpret_cast< complex<r_4>* > (vprocpaq_[fib]);
423 vpdatar_[2*fib+1] = reinterpret_cast< complex<r_4>* >(vprocpaq_[fib]+memgr_.ProcPaqSize()/2) ;
424 }
425 }
426 else { // donnees firmware FFT
427 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
428 vpdata_[2*fib] = vpaq_[fib].Data1C();
429 vpdata_[2*fib+1] = vpaq_[fib].Data2C();
430 }
431 }
432
433 if (totnbpaq_%nmean_ == 0) {
434 if (totnbpaq_ > 0) {
435 moyfc_/=nmean_;
436 moytt_/=nmean_;
437 UpdateVisMtxInfo(); // add/update keywords in the Info DVList
438 // ATTENTION : Matrice visibilites non moyennee
439 char nfile[48];
440 if (nbcalc_==1)
441 sprintf(nfile,"vismtx%d.%s",numfile_,OutFileExtension());
442 else
443 sprintf(nfile,"vismtx_%d_%d.%s",(int)calcid_,numfile_,OutFileExtension());
444 string flnm=outpath_+nfile;
445 if (fgfitsout_) { // Ecriture au format FITS
446 FitsInOutFile fo(flnm, FitsInOutFile::Fits_Create);
447 TArray<r_4> arvismtx = ArrCastC2R(vismtx_);
448 arvismtx.Info()=vismtx_.Info();
449 fo << arvismtx;
450 }
451 else { // Format PPF
452 POutPersist po(flnm);
453 po << vismtx_;
454 }
455 if ((prtlev_>0)&&(numfile_%prtmodulo_==0)) {
456 cout << numfile_ << "-BRVisCalc[" << calcid_ << "/" << nbcalc_ << "]::Process() NPaqProc="
457 << totnbpaq_ << " TotMegaFLOP=" << (uint_8)TotNbMegaFLOP() << " -> VisibMtx in " << flnm << endl;
458 }
459 if (fgvisdt_) FillVisibTable(moyfc_, moytt_);
460 numfile_++;
461 }
462 vismtx_ = complex<r_4>((r_4)0.,(r_4)0.);
463 moyfc_=moytt_=0.;
464 first_fc_=curfc_[0];
465 first_tt_= (vpaq_[0].TimeTag()-ttfirst_[0]);
466 // first_tmstamp_.SetNow(); // Current date and time
467 first_tmstamp_=cts_; // Current date and time
468 }
469
470// kpair=numero sequentiel de la paire: 0->(0,0), 1->(0,1), 2->(0,2), 3->(0,3), 4->(1,1), 5->(1,2) ...
471 sa_size_t kpair=0;
472 sa_size_t k=0; // numero de ligne dans la matrice des visibilites
473 for(size_t i=0; i<vpdata_.size(); i++) {
474 for(size_t j=i; j<vpdata_.size(); j++) {
475 kpair++;
476 if (kpair<(pairst_+1)) continue;
477 if (kpair>=(pairst_+nbpairs_+1)) break;
478 if (fgpimp_&&(i!=j)&&((i+j)%2==0)) continue; // calcul des visib avec numero pair-impair + autocorrel
479 TVector< complex<r_4> > vis = vismtx_.Row(k); k++;
480
481 if (fgdataraw_) { // Donnees firmware RAW apres TF soft
482 for(sa_size_t f=1; f<vis.Size(); f++) {
483 vis(f) += vpdatar_[i][f] * conj(vpdatar_[j][f]);
484 }
485 }
486 else { // donnees firmware FFT
487 for(sa_size_t f=1; f<vis.Size(); f++) {
488 vis(f) += complex<r_4>((r_4)vpdata_[i][f].realB(), (r_4)vpdata_[i][f].imagB()) *
489 complex<r_4>((r_4)vpdata_[j][f].realB(), -(r_4)vpdata_[j][f].imagB());
490 }
491 }
492 nb_flop_ += (8.*(r_8)(vis.Size()-1));
493 }
494 }
495
496 moyfc_ += curfc_[0];
497 moytt_ += (vpaq_[0].TimeTag()-ttfirst_[0]);
498 if ((fgcktt_)&&(calcid_==0)) CheckTimeTag();
499 totnbpaq_++;
500 return 0;
501}
502
503/* --Methode-- */
504void BRVisibilityCalculator::UpdateVisMtxInfo()
505{
506 string ikey,ikdesc;
507 ikey="DATEOBS"; ikdesc=" Date/Time corresponding to TimeTagFirst";
508 vismtx_.Info().SetS(ikey,first_tmstamp_.ToString());
509 vismtx_.Info().SetComment(ikey,ikdesc);
510 ikey="FirstFC"; ikdesc="First FrameCounter";
511 vismtx_.Info().SetI(ikey,first_fc_);
512 vismtx_.Info().SetComment(ikey,ikdesc);
513 ikey="FirstTT"; ikdesc="First TimeTag";
514 vismtx_.Info().SetI(ikey,first_tt_);
515 vismtx_.Info().SetComment(ikey,ikdesc);
516 ikey="LastFC"; ikdesc="Last FrameCounter";
517 vismtx_.Info().SetI(ikey,curfc_[0]);
518 vismtx_.Info().SetComment(ikey,ikdesc);
519 ikey="LastTT"; ikdesc="Last TimeTag";
520 vismtx_.Info().SetI(ikey,vpaq_[0].TimeTag()-ttfirst_[0]);
521 vismtx_.Info().SetComment(ikey,ikdesc);
522 ikey="MeanFC"; ikdesc="Mean FrameCounter";
523 vismtx_.Info().SetD(ikey,moyfc_);
524 vismtx_.Info().SetComment(ikey,ikdesc);
525 ikey="MeanTT"; ikdesc="Mean TimeTag";
526 vismtx_.Info().SetD(ikey,moytt_);
527 vismtx_.Info().SetComment(ikey,ikdesc);
528 ikey="NPAQSUM"; ikdesc="Number of paquets summed";
529 vismtx_.Info().SetI(ikey,nmean_);
530}
531
532
533/* --Methode-- */
534int BRVisibilityCalculator::execute(int tid)
535{
536 vector<TwoByteComplex*> pvpdata(2*memgr_.NbFibres());
537 vector< complex<r_4>* > pvpdatar(2*memgr_.NbFibres());
538 size_t paqsz=memgr_.PaqSize();
539 size_t procpaqsz=memgr_.ProcPaqSize();
540 BRPaquet ppaq(paqsz);
541
542
543 sa_size_t fdelt = vismtx_.NCols()/nparthr_;
544 sa_size_t fdeb = tid*fdelt;
545 sa_size_t ffin = (tid+1)*fdelt;
546
547 if (fdeb<1) fdeb=1;
548 if ((ffin>vismtx_.NCols())||(tid==(nparthr_-1))) ffin=vismtx_.NCols();
549
550 for(size_t jp=0; jp<memgr_.NbPaquets(); jp++) { // boucle sur les paquets d'une zone
551 if (!fgallfibok[jp]) continue;
552 if (fgdataraw_) { // Donnees firmware RAW apres TF soft
553 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
554 Byte* procdatap=fprocbuff_[fib]+jp*procpaqsz;
555 pvpdatar[2*fib] = reinterpret_cast< complex<r_4>* > (procdatap);
556 pvpdatar[2*fib+1] = reinterpret_cast< complex<r_4>* >(procdatap+procpaqsz/2) ;
557 }
558 }
559 else { // donnees firmware FFT
560 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
561 ppaq.Set(fbuff_[fib]+jp*paqsz);
562 pvpdata[2*fib] = ppaq.Data1C();
563 pvpdata[2*fib+1] = ppaq.Data2C();
564 }
565 }
566
567
568// kpair=numero sequentiel de la paire: 0->(0,0), 1->(0,1), 2->(0,2), 3->(0,3), 4->(1,1), 5->(1,2) ...
569 sa_size_t kpair=0;
570 sa_size_t k=0; // numero de ligne dans la matrice des visibilites
571 for(size_t i=0; i<vpdata_.size(); i++) {
572 for(size_t j=i; j<vpdata_.size(); j++) {
573 kpair++;
574 if (kpair<(pairst_+1)) continue;
575 if (kpair>=(pairst_+nbpairs_+1)) break;
576 if (fgpimp_&&(i!=j)&&((i+j)%2==0)) continue; // calcul des visib avec numero pair-impair + autocorrel
577 TVector< complex<r_4> > vis = vismtx_.Row(k); k++;
578 if (fgdataraw_) { // Donnees firmware RAW apres TF soft
579 for(sa_size_t f=fdeb; f<ffin; f++) {
580 vis(f) += pvpdatar[i][f] * conj(pvpdatar[j][f]);
581 }
582 }
583 else { // donnees firmware FFT
584 for(sa_size_t f=fdeb; f<ffin; f++) {
585 vis(f) += complex<r_4>((r_4)pvpdata[i][f].realB(), (r_4)pvpdata[i][f].imagB()) *
586 complex<r_4>((r_4)pvpdata[j][f].realB(), -(r_4)pvpdata[j][f].imagB());
587 }
588 }
589 nb_flop_ += (8.*(r_8)(ffin-fdeb));
590 }
591 }
592
593 } // Fin de boucle sur les paquets
594
595 return 0;
596}
597
598/* --Methode-- */
599int BRVisibilityCalculator::FillVisibTable(double fcm, double ttm)
600{
601 double xnt[10];
602 xnt[0]=fcm; xnt[1]=ttm/1.25e8;
603
604 if (djf_<2) {
605 for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
606 for(sa_size_t jf=jf1_; jf<jf2_; jf++) {
607 xnt[2]=jf;
608 xnt[3]=chanpairnumall_(rv+pairst_);
609 xnt[4]=vismtx_(rv,jf).real()/(r_4)(nmean_);
610 xnt[5]=vismtx_(rv,jf).imag()/(r_4)(nmean_);
611 visdt_.AddRow(xnt);
612 }
613 }
614 }
615 else {
616 for(sa_size_t rv=0; rv<vismtx_.NRows(); rv++) {
617 for(sa_size_t jf=jf1_; jf<jf2_; jf+=djf_) {
618 r_4 moyreal=0.;
619 r_4 moyimag=0.;
620 sa_size_t jjfmx=jf+djf_;
621 if (jjfmx > vismtx_.NCols()) jjfmx=vismtx_.NCols();
622 for(sa_size_t jjf=jf; jjf<jjfmx; jjf++) {
623 moyreal+=vismtx_(rv,jjf).real();
624 moyimag+=vismtx_(rv,jjf).imag();
625 }
626 xnt[2]=jf+djf_/2;
627 xnt[3]=chanpairnumall_(rv+pairst_);
628 xnt[4]=moyreal/(r_4)(nmean_*djf_);
629 xnt[5]=moyimag/(r_4)(nmean_*djf_);
630 visdt_.AddRow(xnt);
631 }
632 }
633 }
634 return 0;
635}
636
637/* --Methode-- */
638int BRVisibilityCalculator::CheckTimeTag()
639{
640 if (totnbpaq_==0) {
641 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
642 vlasttt_[fib]=ttfirst_[fib];
643 if (ttmtx_.NCols()>0) {
644 fcmtx_(fib,totnbpaq_) = curfc_[fib];
645 ttmtx_(fib,totnbpaq_) = vlasttt_[fib];
646 }
647 }
648 return 0;
649 }
650 for(size_t fib=0; fib<(size_t)memgr_.NbFibres(); fib++) {
651 int_8 ld = (int_8)vpaq_[fib].TimeTag()-(int_8)vlasttt_[fib];
652 int_8 fd = (int_8)vpaq_[fib].TimeTag()-(int_8)ttfirst_[fib]-(int_8)vpaq_[0].TimeTag()+(int_8)ttfirst_[0];
653 /* if ( (ld < mindeltatt_) || (fd<-5) || (fd>5)) { vbadtt_[fib]++; vnsamett_[fib]++; }
654 else {
655 if (fd!=0) vnsamett_[fib]++;
656 }
657 */
658 if (ld < mindeltatt_) vbadtt_[fib]++;
659 else {
660 if (fd != 0) vnsamett_[fib]++;
661 if ((fd<-5)||(fd>5)) vndiff5tt_[fib]++;
662 }
663 vlasttt_[fib]=vpaq_[fib].TimeTag();
664 if (totnbpaq_<ttmtx_.NCols()) {
665 fcmtx_(fib,totnbpaq_) = curfc_[fib];
666 ttmtx_(fib,totnbpaq_) = vlasttt_[fib];
667 }
668 }
669 return 0;
670}
671
672//-------------------------------------------------------------------------------
673// Classe Groupe (ensemble) de Calculateur de Visibilites, tournant en parallele
674//-------------------------------------------------------------------------------
675
676/* --Methode-- */
677BRVisCalcGroup::BRVisCalcGroup(size_t nbcalc, RAcqMemZoneMgr& memgr, string outpath, uint_4 nmean,
678 uint_4 pair1, uint_4 nbpairs, bool fgpimp, size_t nthr)
679 : tm_(false)
680{
681 if ((nbcalc<1)||(nbcalc>6))
682 throw ParmError("BRVisCalcGroup::BRVisCalcGroup NbCalc > 6 !");
683 for(size_t i=0; i<nbcalc; i++) {
684 BRVisibilityCalculator * viscp=new BRVisibilityCalculator(memgr, outpath, nmean, nthr);
685 viscp->DefineRank(nbcalc, i, pair1, nbpairs, fgpimp);
686 viscalcp_.push_back(viscp);
687 }
688}
689/* --Methode-- */
690BRVisCalcGroup::~BRVisCalcGroup()
691{
692 for(size_t i=0; i<viscalcp_.size(); i++)
693 delete viscalcp_[i];
694}
695/* --Methode-- */
696MemZStatus BRVisCalcGroup::SetMemZAction(MemZaction mmact)
697{
698 MemZaction mmzas[10]={MemZA_ProcA,MemZA_ProcB,MemZA_ProcC,MemZA_ProcD,MemZA_ProcE,MemZA_ProcF,
699 MemZA_ProcG,MemZA_ProcH,MemZA_ProcI,MemZA_ProcJ};
700 MemZStatus mmzss[10]={MemZS_ProcA,MemZS_ProcB,MemZS_ProcC,MemZS_ProcD,MemZS_ProcE,MemZS_ProcF,
701 MemZS_ProcG,MemZS_ProcH,MemZS_ProcI,MemZS_ProcJ};
702 size_t ia=9999;
703 for(int i=0; i<10; i++) {
704 if (mmact==mmzas[i]) { ia=i; break; }
705 }
706 if (ia>=10)
707 throw ParmError("BRVisCalcGroup::SetMemZAction() Bad MemZaction mmact !");
708 if ((ia+viscalcp_.size())>10)
709 throw ParmError("BRVisCalcGroup::SetMemZAction() MemZaction mmact too high !");
710 for(size_t i=0; i<viscalcp_.size(); i++) {
711 viscalcp_[i]->SetMemZAction(mmzas[ia]); ia++;
712 }
713 return mmzss[ia-1];
714}
715/* --Methode-- */
716int BRVisCalcGroup::SelectFreqBinning(uint_4 freq1, uint_4 freq2, uint_4 nbfreq)
717{
718 int rc=0;
719 for(size_t i=0; i<viscalcp_.size(); i++)
720 rc=viscalcp_[i]->SelectFreqBinning(freq1, freq2, nbfreq);
721 return rc;
722}
723/* --Methode-- */
724void BRVisCalcGroup::ActivateVisDTable(bool fgfdt)
725{
726 for(size_t i=0; i<viscalcp_.size(); i++)
727 viscalcp_[i]->ActivateVisDTable(fgfdt);
728}
729/* --Methode-- */
730void BRVisCalcGroup::SetFitsOutput()
731{
732 for(size_t i=0; i<viscalcp_.size(); i++)
733 viscalcp_[i]->SetFitsOutput();
734}
735/* --Methode-- */
736void BRVisCalcGroup::SetPPFOutput()
737{
738 for(size_t i=0; i<viscalcp_.size(); i++)
739 viscalcp_[i]->SetPPFOutput();
740}
741/* --Methode-- */
742void BRVisCalcGroup::SetFFTData()
743{
744 for(size_t i=0; i<viscalcp_.size(); i++)
745 viscalcp_[i]->SetFFTData();
746}
747/* --Methode-- */
748void BRVisCalcGroup::SetRawData()
749{
750 for(size_t i=0; i<viscalcp_.size(); i++)
751 viscalcp_[i]->SetRawData();
752}
753/* --Methode-- */
754void BRVisCalcGroup::SetPrintLevel(int lev, uint_8 prtmodulo)
755{
756 for(size_t i=0; i<viscalcp_.size(); i++)
757 viscalcp_[i]->SetPrintLevel(lev,prtmodulo);
758}
759/* --Methode-- */
760void BRVisCalcGroup::start()
761{
762 for(size_t i=0; i<viscalcp_.size(); i++)
763 viscalcp_[i]->start();
764 tm_.SplitQ();
765}
766/* --Methode-- */
767void BRVisCalcGroup::join()
768{
769 r_8 totflop=0.;
770 for(size_t i=0; i<viscalcp_.size(); i++) {
771 viscalcp_[i]->join();
772 // cout << " BRVisCalcGroup::join()/ VisibCalc[" << i << "]->TotNbMegaFLOP()="
773 // << viscalcp_[i]->TotNbFLOP()/1024e3 << endl;
774 totflop += viscalcp_[i]->TotNbFLOP();
775 }
776 tm_.SplitQ();
777 cout << " ----------------------------------------------------------" << endl;
778 cout << " BRVisCalcGroup::join() : Finished " << viscalcp_.size() << " VisibilityCalculator(s)" << endl;
779 cout << " ... Elaspsed time: " << tm_.PartialElapsedTimems()
780 << " ms (total:" << tm_.TotalElapsedTimems() << ")" << endl;
781 double mflopsrate=totflop/(r_8)tm_.PartialElapsedTimems()/(1024.);
782 cout << " ... TotalMegaFLOP= " << totflop/(1024.e3) << " @ "
783 << mflopsrate << " MFLOP/s" << " (=" << mflopsrate/(r_8)viscalcp_.size() << "/VisibCalcObject)" << endl;
784 cout << " ----------------------------------------------------------" << endl;
785 return;
786}
Note: See TracBrowser for help on using the repository browser.