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

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

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

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