source: Sophya/trunk/Poubelle/archTOI.old/toiiter.cc@ 398

Last change on this file since 398 was 394, checked in by ansari, 26 years ago

Integration detecteur d'etoiles DY
File size: 17.9 KB
Line 
1// toiiter.cc
2// Eric Aubourg CEA/DAPNIA/SPP juillet 1999
3
4
5#include "toiiter.h"
6#include "toiinterpolator.h"
7#include "archparam.h"
8#include <dirent.h>
9#include <iostream.h>
10
11// Format bloc GPS
12// $GPGGA,hhmmss.ss,ddmm.mmmm,n,dddmm.mmmm,e,q,ss,y.y,a.a,z,
13
14TOIIter::TOIIter() {
15 // Toutes les valeurs initialisees par defaut. Le TOISvr les positionnera,
16 // puis lancera l'initialisation
17
18 file = NULL;
19 directory = "";
20
21 files.clear();
22
23 isOnBoardRecorder = false;
24
25 imes=0;
26
27 tStart = -9.e99;
28 tEnd = 9.e99;
29
30 //tBlock0 = -1;
31 // perEch = -1;
32 // Values for Trapani 99 = default values
33 //tBlock0 = 1376.8358818;
34 //perEch = 0.005836818076;
35
36 trigMask = 0;
37 rawIter = NULL;
38 interp = NULL;
39 lastSample = -1;
40 maxLookAhead = 10000;
41
42 auxGPS = NULL;
43
44 initDone = false;
45}
46
47TOIIter::TOIIter(TOIIter const& x) {
48 directory = x.directory;
49 files = x.files;
50 // curFile = x.curFile; // $CHECK$ DANGER !!
51 curFile = files.find(*(x.curFile));
52 isOnBoardRecorder = x.isOnBoardRecorder;
53 imes = x.imes;
54 tStart = x.tStart;
55 tEnd = x.tEnd;
56 //tBlock0 = x.tBlock0;
57 //perEch = x.perEch;
58 trigMask = x.trigMask;
59 infos = x.infos;
60
61 if (x.file)
62 file = new ArcheopsFile(*x.file);
63 else
64 file = NULL;
65
66 if (x.rawIter) {
67 rawIter = new TOIIter(*x.rawIter);
68 interp = new TOIInterpolator[infos.size()];
69 for (int i=0; i<infos.size(); i++)
70 interp[i] = x.interp[i];
71 } else {
72 rawIter = NULL;
73 interp = NULL;
74 }
75
76 lastSample = x.lastSample;
77 maxLookAhead = x.maxLookAhead;
78
79 auxGPS = x.auxGPS;
80 if (auxGPS) auxGPS = auxGPS->clone();
81
82 initDone = x.initDone;
83}
84
85TOIIter::~TOIIter() {
86 delete file;
87 delete rawIter;
88 delete[] interp;
89 delete auxGPS;
90}
91
92#ifdef __MWERKS__
93#define filesep ':'
94#else
95#define filesep '/'
96#endif
97
98void TOIIter::Init() {
99 if (initDone) return;
100 initDone = true;
101 // On a soit un repertoire, soit une liste de fichiers....
102 if (directory == "") {
103 if (files.empty()) { // Ni repertoire, ni fichiers
104 cerr << "toiiter : pas de repertoire, pas de fichiers" << endl;
105 exit(-1);
106 } else {
107 // On a deja une liste de fichiers
108 }
109 } else { // On a un repertoire a explorer
110 // On cherche soit les fichiers dans le repertoire donne, soit des fichiers
111 // dans un sous-repertoire "arch-YY_MM_DD". Les fichiers ont un nom en
112 // "hYY_MM_DD-hh_mm_ss".
113 // Pour l'enregistreur de vol, les fichiers ont un nom en ARKxxxxxx.DAT
114 if (directory[directory.length()-1] != filesep)
115 directory += filesep;
116 DIR* dir = opendir(directory.c_str());
117 struct dirent* ent;
118 while ((ent = readdir(dir)) != NULL) {
119 // si c'est un repertoire, avec un nom de jour, il faut l'explorer...
120 if (!strncmp(ent->d_name, "arch-", 5)) {
121 double mjd = ArcheopsFile::decodeMJD(ent->d_name+5) - 2./24.; // ENTIER + .5 en temps local!
122 if (mjd >= tStart - 1. && mjd <= tEnd) {
123 string direc2 = directory + ent->d_name + filesep;
124 DIR* dir2 = opendir(direc2.c_str());
125 struct dirent* ent2;
126 while ((ent2 = readdir(dir2)) != NULL) {
127 if (*ent2->d_name == 'h') {
128 double mjd2 = ArcheopsFile::decodeMJD(ent->d_name+1) - 2./24.;
129 if (mjd2 >= tStart - 1./24. && mjd2 <= tEnd) {
130 files.insert(direc2 + ent2->d_name);
131 }
132 }
133 }
134 }
135 } else {
136 if (!isOnBoardRecorder && *ent->d_name == 'h') {
137 double mjd = ArcheopsFile::decodeMJD(ent->d_name+1) - 2./24.; // $CHECK$ UTCOffset
138 if (mjd >= tStart - 1./24. && mjd <= tEnd) {
139 files.insert(directory + ent->d_name);
140 }
141 } else if (isOnBoardRecorder) {
142 if (strncmp(ent->d_name, "ARK", 3) && strncmp(ent->d_name, "ark", 3)) continue;
143 char * sfx = ent->d_name + strlen(ent->d_name) - 4;
144 if (strcmp(sfx, ".DAT") && strcmp(sfx, ".dat")) continue;
145 files.insert(directory + ent->d_name);
146 }
147 }
148 }
149 closedir(dir);
150 }
151
152 ScanFiles();
153
154 curFile = files.begin();
155 file = new ArcheopsFile((*curFile).c_str());
156 cout << "opening file " << (*curFile).c_str() << endl;
157
158 // On avance jusqu'a avoir au moins un bloc param et un bloc reglage,
159 // car on ne peut rien faire sans...
160 // Si on a des donnees de l'enregistreur de vol, pas de bloc param, et
161 // on en simule un
162 double oldTStart = tStart;
163 tStart = -9.e99; // pour init, on accepte des blocs avant tstart....
164
165 if (!file->lastParam()) {
166 if (isOnBoardRecorder) {
167 extern param_bolo parametr;
168 block_type_param block;
169 block.param = parametr;
170 valide_block((block_type_modele*)&block, block_param, 0);
171 file->forceBlock((block_type_modele*)&block);
172 } else {
173 file->nextBlock(block_param_mask);
174 }
175 }
176 if (!file->lastReglage()) file->nextBlock(block_reglage_mask);
177
178 // On cherche un bloc GPS pour avoir la correspondance timestamp/UTC.
179 // Pour le moment, on se fonde sur le premier bloc GPS. On pourra faire
180 // mieux en prenant le min de tous les delta_T, a condition d'avoir un
181 // peu plus de details sur la facon dont le GPS est lu.
182
183 if (archParam.acq.tBlock0 < 0) {
184 archParam.acq.tBlock0 = file->getStartMJD();
185
186 file->pushMark();
187 if (file->lastGPS() || file->nextBlock(block_gps_mask)) {
188 // le temps du bloc courant, en secondes
189 double dt = file->blockNum() * file->perBlock();
190 archParam.acq.tBlock0 = file->getGPSMJD() - dt/86400.;
191 } else { // pas de bloc GPS...
192 archParam.acq.tBlock0 = file->getStartMJD();
193 }
194 file->popMark();
195 }
196 tStart = oldTStart; // on restaure
197
198 if (archParam.acq.perEch < 0)
199 archParam.acq.perEch = file->perEchant();
200
201 bool hasInterp = false;
202
203 trigMask = 0;
204 for (vector<info>::iterator i = infos.begin(); i != infos.end(); i++) {
205 if ((*i).interpolated) hasInterp = true;
206 if ((*i).triggering) {
207 switch ((*i).kind) {
208 case boloTens:
209 case boloRaw:
210 trigMask |= block_bolo_mask;
211 break;
212 case gpsTime:
213 case longitude:
214 case latitude:
215 case altitude:
216 trigMask |= block_gps_mask;
217 break;
218 case azimut:
219 file->needSSTProcessMask(SSTHandler::findPeriod);
220 trigMask |= block_sst_mask;
221 break;
222 case sstStarZ:
223 case sstStarF:
224 case sstStarT:
225 file->needSSTProcessMask(SSTHandler::findStars);
226 trigMask |= block_sst_mask;
227 break;
228 case sstDiode:
229 file->needSSTProcessMask(SSTHandler::permDiode);
230 trigMask |= block_sst_mask;
231 break;
232 case sstChannel:
233 trigMask |= block_sst_mask;
234 break;
235 case gyroRaw:
236 trigMask |= block_gyro_mask;
237 break;
238 case alphaAxis:
239 case deltaAxis:
240 case alphaBolo:
241 case deltaBolo:
242 file->needSSTProcessMask(SSTHandler::findAxis);
243 trigMask |= block_sst_mask;
244 break;
245 }
246 }
247 }
248
249 if (trigMask & (block_bolo_mask | block_sst_mask)) {
250 imes = 9999;
251 } else {
252 imes = 0;
253 }
254
255 if (hasInterp) {
256 rawIter = new TOIIter(*this);
257 interp = new TOIInterpolator[infos.size()];
258 for (int i=0; i<infos.size(); i++) {
259 rawIter->infos[i].interpolated = false;
260 }
261 delete file; file = NULL; // on ne travaille plus sur le fichier directement...
262 }
263}
264
265void TOIIter::ScanFiles() {
266 file1stSamp.clear();
267 cout << "Scanning all files" << endl;
268 // Petite astuce pour les STL non conformes comme celles de digital
269 // qui ne supportent pas files.erase(i) suivi de i++....
270 set<string> copy = files;
271 for (set<string>::iterator i = copy.begin(); i != copy.end(); i++) {
272 ArcheopsFile fich((*i).c_str());
273 if (fich.nextBlock()) {
274 file1stSamp[*i] = fich.blockNum()*72; // premier numsample
275 cout << "File " << *i << " 1st sample = " << fich.blockNum()*72 << endl;
276 } else {
277 cout << "File " << *i << " unrecoverable, skipping" << endl;
278 files.erase(*i);
279 }
280 }
281 cout << "Scan done" << endl;
282
283 // Et maintenant, on ne garde que ceux qui tombent dans l'intervalle...
284 copy = files;
285 string prev="";
286 for (set<string>::iterator i = copy.begin(); i != copy.end(); i++) {
287 double smp = file1stSamp[*i];
288 double t = archParam.acq.tBlock0 + smp * archParam.acq.perEch/86400.;
289 if (t>tEnd) { // premier echantillon apres tEnd
290 files.erase(*i);
291 prev = "";
292 continue;
293 }
294 if (t<tStart) { // premier echantillon avant tStart -> on vire le precedent si existe
295 if (prev != "") {
296 files.erase(prev);
297 }
298 }
299 prev = *i;
300 }
301}
302
303bool TOIIter::NextFile() {
304 if (rawIter)
305 return rawIter->NextFile();
306
307 if (files.empty()) {
308 return false;
309 } else {
310 if (curFile == files.end()) return false;
311 curFile++;
312 if (curFile == files.end()) return false;
313 cout << "opening file " << (*curFile).c_str() << endl;
314 ArcheopsFile* newfile = new ArcheopsFile((*curFile).c_str());
315 newfile->grabLastBlocs(*file);
316 delete file;
317 file = newfile;
318 return true;
319 }
320}
321
322
323bool TOIIter::Next() {
324 if (!initDone) Init();
325 while (1) {
326 if (!NextSample()) return false; // end of files
327 double t = getMJD();
328 if (t < tStart) continue;
329 if (t > tEnd) return false;
330 return true;
331 }
332}
333
334bool TOIIter::NextSample() {
335 if (rawIter) { // Delegation pour interpolation
336 // Trouve prochain sample disponible
337 for (int k=0; k<2; k++) {
338 long smp = 2147483647L;
339 for (int i=0; i<infos.size(); i++) {
340 long ss = interp[i].nextSample(lastSample+1);
341 if (ss > 0 && ss < smp) smp=ss;
342 }
343 if (smp != 2147483647L) {
344 lastSample = smp;
345 break;
346 }
347 if (!fetchAhead()) // tout le monde etait en bout de course,
348 return false; // on lit un echantillon, ca suffit, d'ou le k<2
349 }
350 // Verifie que tous les interpolateurs ont assez de donnees pour
351 // trouver la valeur correspondante
352 for (int i=0; i<infos.size(); i++) {
353 //rif (infos[i].interpolated)
354 while (interp[i].needMoreDataFor(lastSample) &&
355 rawIter->getSampleIndex() - lastSample < maxLookAhead)
356 if (!fetchAhead()) return false;
357 }
358
359 // On est pret...
360 return true;
361 }
362
363 // trigger sur info indexee dans bloc bolo ou bloc sst ?
364 if (trigMask & (block_bolo_mask | block_sst_mask)) {
365 imes++;
366 if (imes < file->nEchBlock()) return true;
367 imes = 0;
368 }
369
370 // soit pas d'info indexee, soit fin bloc courant...
371 while (1) {
372 if (file->nextBlock(trigMask)) {
373 while (file->sameBlockNumAhead()) { // tant que meme numero de bloc, on lit
374 if (!file->nextBlock()) { // fin de fichier ?
375 if (NextFile()) file->nextBlock(); // fichier suivant
376 else return false; // tout fini
377 }
378 }
379 return true;
380 }
381 if (!NextFile()) return false;
382 }
383}
384
385/* double TOIIter::getTime() { // MJD
386 // le temps du bloc courant, en secondes
387 double dt = file->blockNum() * file->perBlock();
388 return tBlock0 + dt/86400. + imes*file->perEchant()/86400.;
389 }
390 */
391
392bool TOIIter::canGetValue(int column) {
393 if (!initDone) Init();
394 if (column < 0 || column >= infos.size()) return false;
395 TOIKind kind = infos[column].kind;
396 if (auxGPS &&
397 (kind == longitude || kind == latitude || kind == altitude || kind == tsid)) {
398 double dummy;
399 return auxGPS->getLocation(getMJD(), dummy, dummy, dummy) == 0;
400 }
401 if (rawIter) {
402 return interp[column].canGet(lastSample);
403 }
404 int index = infos[column].index;
405 switch (kind) {
406 case sampleNum:
407 case internalTime:
408 case mjd:
409 return true;
410 case boloTens:
411 if (imes==0 && file->llastBolo()==NULL) return false;
412 return file->lastBolo() != NULL;
413 case boloRaw:
414 return file->lastBolo() != NULL;
415 case sstDiode:
416 case sstChannel:
417 return file->lastSST() != NULL;
418 case sstStarZ:
419 case sstStarF:
420 case sstStarT:{
421 if (file->lastSST() == NULL) return false;
422 int n = file->getNumbStar(imes);
423 return (n > 0 && imes < n);
424 }
425 case gyroRaw:
426 return (file->lastGyro() != NULL);
427 case gpsTime:
428 return file->hasGPSTime();
429 case longitude:
430 case latitude:
431 return file->hasGPSPos();
432 case altitude:
433 return file->hasGPSAlt();
434 case tsid:
435 return file->hasGPSPos();
436 case azimut:
437 case alphaAxis:
438 case deltaAxis:
439 return false;
440 //return (file->lastGPS() != NULL && file->lastSST() != NULL);
441 case alphaSst:
442 case deltaSst:
443 case alphaBolo:
444 case deltaBolo:
445 return false;
446
447 case boloTemp:
448 return false;
449 }
450 return false;
451}
452
453double TOIIter::getValue(int column) {
454 if (!initDone) Init();
455 if (column < 0 || column >= infos.size()) return -1;
456 TOIKind kind = infos[column].kind;
457 if (auxGPS &&
458 (kind == longitude || kind == latitude || kind == altitude || kind == tsid)) {
459 double lat,lon,alt;
460 if (auxGPS->getLocation(getMJD(), lat, lon, alt)) return -99999;
461 if (kind == longitude) return lon;
462 if (kind == latitude) return lat;
463 if (kind == altitude) return alt;
464 if (kind == tsid) {
465 tSid.setLongitude(lon);
466 return tSid.getLST(getMJD());
467 }
468 }
469 if (rawIter) {
470 if (infos[column].interpolated)
471 return interp[column].getIValue(lastSample);
472 else
473 return interp[column].getEValue(lastSample);
474 }
475 int index = infos[column].index;
476 switch (kind) {
477 case sampleNum:
478 return getSampleIndex();
479 case internalTime:
480 return getSampleIndex() * archParam.acq.perEch;
481 case mjd:
482 return getMJD();
483 case boloTens:
484 return file->getMuVBolo(index, imes);
485 case boloRaw:
486 return file->getRawBolo(index, imes);
487 case sstDiode:
488 return file->getSSTSignal(index, imes);
489 case sstChannel:
490 return file->getSSTRawSignal(index, imes);
491 case sstStarZ:
492 return file->getSSTStarZ(index, imes);
493 case sstStarF:
494 return file->getSSTStarF(index, imes);
495 case sstStarT:
496 return file->getSSTStarT(index, imes);
497 case gyroRaw:
498 return file->getGyro(index, imes);
499 case gpsTime:
500 return file->getGPSUTC();
501 case longitude:
502 return file->getGPSLong();
503 case latitude:
504 return file->getGPSLat();
505 case altitude:
506 return file->getGPSAlt();
507 case tsid:
508 tSid.setLongitude(file->getGPSLong());
509 return tSid.getLST(getMJD());
510 case azimut:
511 return file->getAzimut(imes);
512 case alphaAxis:
513 return file->getAlpha(imes);
514 case deltaAxis:
515 return file->getDelta(imes);
516 }
517 return -1;
518}
519
520bool TOIIter::newValue(int column) {
521 if (!initDone) Init();
522 if (column < 0 || column >= infos.size()) return false;
523 TOIKind kind = infos[column].kind;
524 if (auxGPS &&
525 (kind == longitude || kind == latitude || kind == altitude)) {
526 return true;
527 }
528 if (rawIter) {
529 return interp[column].isNewValue(lastSample);
530 }
531 switch (kind) {
532 case sampleNum:
533 case internalTime:
534 case mjd:
535 case tsid:
536 return true;
537 case boloTens:
538 return file->blockNum() == file->getBoloBlockNum();
539 case boloRaw:
540 return file->blockNum() == file->getBoloBlockNum();
541 case sstChannel:
542 case sstDiode:
543 case sstStarZ:
544 case sstStarF:
545 case sstStarT:
546 return file->blockNum() == file->getSSTBlockNum();
547 case gyroRaw:
548 return file->blockNum() == file->getGyroBlockNum();
549 case gpsTime:
550 return file->blockNum() == file->getGPSBlockNum() && imes==0;
551 case longitude:
552 return file->blockNum() == file->getGPSBlockNum() && imes==0;
553 case latitude:
554 return file->blockNum() == file->getGPSBlockNum() && imes==0;
555 case altitude:
556 return file->blockNum() == file->getGPSBlockNum() && imes==0;
557 case azimut:
558 return true; // $CHECK$ with SSTHandler
559 case alphaAxis:
560 return true; // $CHECK$ with SSTHandler
561 case deltaAxis:
562 return true; // $CHECK$ with SSTHandler
563 }
564 return false;
565}
566
567bool TOIIter::extendValue(int column) {
568 return (!infos[column].interpolated && !newValue(column));
569}
570
571bool TOIIter::interpValue(int column) {
572 return (infos[column].interpolated && !newValue(column));
573}
574
575bool TOIIter::isTrig(int column) {
576 if (column < 0 || column >= infos.size()) return false;
577 return infos[column].triggering;
578}
579
580
581TOIKind TOIIter::getKind(int column) {
582 if (column < 0 || column >= infos.size()) return (TOIKind)-1;
583 return infos[column].kind;
584}
585
586int TOIIter::getIndex(int column) {
587 if (column < 0 || column >= infos.size()) return (TOIKind)-1;
588 return infos[column].index;
589}
590
591int TOIIter::getColTOI(TOIKind kind, int index) {
592 for (int i=0; i<infos.size(); i++)
593 if (infos[i].kind == kind && infos[i].index == index)
594 return i;
595 return -1;
596}
597
598bool TOIIter::canGetTOI(TOIKind kind, int index) {
599 int col = getColTOI(kind, index);
600 if (col<0) return false;
601 return canGetValue(col);
602}
603
604double TOIIter::getTOI(TOIKind kind, int index) {
605 int col = getColTOI(kind, index);
606 if (col<0) return -9.e99;
607 return getValue(col);
608}
609
610
611int TOIIter::getBlockSampleIndex() {
612 return imes;
613}
614
615int TOIIter::getSampleIndex() {
616 if (!initDone) Init();
617 if (file) {
618 return file->blockNum() * file->nEchBlock() + imes;
619 } else {
620 return lastSample;
621 }
622}
623
624double TOIIter::getMJD() {
625 if (!initDone) Init();
626 int sample = getSampleIndex();
627 return archParam.acq.tBlock0 + sample*archParam.acq.perEch/86400.;
628}
629
630bool TOIIter::fetchAhead() { // Seulement si delegation
631 if (!rawIter) return false;
632 if (!rawIter->Next()) return false;
633 long sample = rawIter->getSampleIndex();
634 for (int i=0; i<infos.size(); i++) {
635 if (rawIter->canGetValue(i) && rawIter->newValue(i)) {
636 interp[i].enterValue(rawIter->getValue(i), sample);
637 }
638 }
639 return true;
640}
641
642block_type_param* TOIIter::lastParam() {
643 if (!initDone) Init();
644 if (file) return file->lastParam();
645 else return rawIter->lastParam();
646}
647
648
Note: See TracBrowser for help on using the repository browser.