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

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

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