source: Sophya/trunk/Poubelle/archTOI.old/archeopsfile.cc@ 310

#define utilitaires_de_block_archeops
#include "archeopsfile.h"
#include "compress.h"
#include <iostream.h>

class BlockSet {
public:
  BlockSet();
  BlockSet(BlockSet const&);
  ~BlockSet();
  void setBloc(block_type_modele const& blk);
  void setRawBloc(block_type_modele const& blk);

  block_type_param*         lastParam;
  block_type_journal*       lastJournal;
  block_type_reglage*       lastReglage;
  block_type_dilution*      lastDilution;
  block_type_gps*           lastGPS;
  block_type_une_periode*   lastUnePeriode;
  block_type_synchro_sol*   lastSynchroSol;
  block_type_pointage_sol*  lastPointageSol;
  block_type_bolo*          lastBolo;   // can be uncompressed bolo_comp
  block_type_bolo*          llastBolo;   
  block_type_gyro*          lastGyro;
  block_type_sst*           lastSST;
  block_type_bolo_comprime* lastBoloComp;   // can be uncompressed bolo_comp
  block_type_gyro_comprime* lastGyroComp;
  block_type_sst_comprime*  lastSSTComp;
};

BlockSet::BlockSet() {
  lastParam      = NULL;
  lastJournal    = NULL;
  lastReglage    = NULL;
  lastDilution   = NULL;
  lastGPS        = NULL;
  lastUnePeriode = NULL;
  lastSynchroSol = NULL;
  lastPointageSol= NULL;
  lastBolo       = NULL;   
  llastBolo      = NULL;   
  lastGyro       = NULL;
  lastSST        = NULL;
  lastBoloComp   = NULL;   
  lastGyroComp   = NULL;
  lastSSTComp    = NULL;
}

BlockSet::~BlockSet() {
  delete lastParam;
  delete lastJournal;
  delete lastReglage;
  delete lastDilution;
  delete lastGPS;
  delete lastUnePeriode;
  delete lastSynchroSol;
  delete lastPointageSol;
  delete lastBolo;   
  delete llastBolo;   
  delete lastGyro;
  delete lastSST;
  delete lastBoloComp;   
  delete lastGyroComp;
  delete lastSSTComp;
}

BlockSet::BlockSet(BlockSet const& x) {
  lastParam      = NULL;
  lastJournal    = NULL;
  lastReglage    = NULL;
  lastDilution   = NULL;
  lastGPS        = NULL;
  lastUnePeriode = NULL;
  lastSynchroSol = NULL;
  lastPointageSol= NULL;
  lastBolo       = NULL;   
  llastBolo      = NULL;   
  lastGyro       = NULL;
  lastSST        = NULL;
  lastBoloComp   = NULL;   
  lastGyroComp   = NULL;
  lastSSTComp    = NULL;
  if (x.lastParam) {
    lastParam = new block_type_param;
    *lastParam = *x.lastParam;
  }
  if (x.lastJournal) {
    lastJournal = new block_type_journal;
    *lastJournal = *x.lastJournal;
  }
  if (x.lastReglage) {
    lastReglage = new block_type_reglage;
    *lastReglage = *x.lastReglage;
  }
  if (x.lastDilution) {
    lastDilution = new block_type_dilution;
    *lastDilution = *x.lastDilution;
  }
  if (x.lastGPS) {
    lastGPS = new block_type_gps;
    *lastGPS = *x.lastGPS;
  }
  if (x.lastUnePeriode) {
    lastUnePeriode = new block_type_une_periode;
    *lastUnePeriode = *x.lastUnePeriode;
  }
  if (x.lastSynchroSol) {
    lastSynchroSol = new block_type_synchro_sol;
    *lastSynchroSol = *x.lastSynchroSol;
  }
  if (x.lastPointageSol) {
    lastPointageSol = new block_type_pointage_sol;
    *lastPointageSol = *x.lastPointageSol;
  }
  if (x.lastBolo) {
    lastBolo = new block_type_bolo;
    *lastBolo = *x.lastBolo;
  }
  if (x.llastBolo) {
    llastBolo = new block_type_bolo;
    *llastBolo = *x.llastBolo;
  }
  if (x.lastGyro) {
    lastGyro = new block_type_gyro;
    *lastGyro = *x.lastGyro;
  }
  if (x.lastSST) {
    lastSST = new block_type_sst;
    *lastSST = *x.lastSST;
  }
  if (x.lastBoloComp) {
    lastBoloComp = new block_type_bolo_comprime;
    *lastBoloComp = *x.lastBoloComp;
  }
  if (x.lastGyroComp) {
    lastGyroComp = new block_type_gyro_comprime;
    *lastGyroComp = *x.lastGyroComp;
  }
  if (x.lastSSTComp) {
    lastSSTComp = new block_type_sst_comprime;
    *lastSSTComp = *x.lastSSTComp;
  }
}

void BlockSet::setBloc(block_type_modele const& blk)
{
  int kind = type_block(&blk);
  switch (kind) {
    case block_param:
      if (!lastParam) lastParam = new block_type_param;
      memcpy(lastParam, &blk, sizeof(block_type_param));
      // Les fichiers fournis contiennent des valeurs debiles...
      if (lastParam->param.n_max_bolo < 0)
        lastParam->param.n_max_bolo = nb_max_bolo;
      if (lastParam->param.n_per_block < 0)
        lastParam->param.n_per_block = nb_per_block;
      if (lastParam->param.n_max_mes_per < 0)
        lastParam->param.n_max_mes_per = nb_max_mes_per;
      break;
    case block_journal:
      if (!lastJournal) lastJournal = new block_type_journal;
      memcpy(lastJournal, &blk, sizeof(block_type_journal));
      break;
    case block_reglage:
      if (!lastReglage) lastReglage = new block_type_reglage;
      memcpy(lastReglage, &blk, sizeof(block_type_reglage));
      break;
    case block_dilution:
      if (!lastDilution) lastDilution = new block_type_dilution;
      memcpy(lastDilution, &blk, sizeof(block_type_dilution));
      break;
    case block_gps:
      if (!lastGPS) lastGPS = new block_type_gps;
      memcpy(lastGPS, &blk, sizeof(block_type_gps));
      break;
    case block_une_periode:
      if (!lastUnePeriode) lastUnePeriode = new block_type_une_periode;
      memcpy(lastUnePeriode, &blk, sizeof(block_type_une_periode));
      break;
    case block_synchro_sol:
      if (!lastSynchroSol) lastSynchroSol = new block_type_synchro_sol;
      memcpy(lastSynchroSol, &blk, sizeof(block_type_synchro_sol));
      break;
    case block_pointage_sol:
      if (!lastPointageSol) lastPointageSol = new block_type_pointage_sol;
      memcpy(lastPointageSol, &blk, sizeof(block_type_pointage_sol));
      break;
    case block_bolo:
      if (!lastBolo) lastBolo = new block_type_bolo;
      else if (!llastBolo) {
        llastBolo = new block_type_bolo;
        memcpy(llastBolo, lastBolo, sizeof(block_type_bolo));
      }
      memcpy(lastBolo, &blk, sizeof(block_type_bolo));
      break;
    case block_gyro:
      if (!lastGyro) lastGyro = new block_type_gyro;
      memcpy(lastGyro, &blk, sizeof(block_type_gyro));
      break;
    case block_sst:
      if (!lastSST) lastSST = new block_type_sst;
      memcpy(lastSST, &blk, sizeof(block_type_sst));
      break;
  }
}

void BlockSet::setRawBloc(block_type_modele const& blk)
{
  int kind = type_block(&blk);
  switch (kind) {
    case block_bolo_comprime:
      if (!lastBoloComp) lastBoloComp = new block_type_bolo_comprime;
      memcpy(lastBoloComp, &blk, sizeof(block_type_bolo_comprime));
      break;
    case block_gyro_comprime:
      if (!lastGyroComp) lastGyroComp = new block_type_gyro_comprime;
      memcpy(lastGyroComp, &blk, sizeof(block_type_gyro_comprime));
      break;
    case block_sst_comprime:
      if (!lastSSTComp) lastSSTComp = new block_type_sst_comprime;
      memcpy(lastSSTComp, &blk, sizeof(block_type_sst_comprime));
      break;
  }
}



ArcheopsFile::ArcheopsFile(string const& fname) {
  f = fopen(fname.c_str(), "rb");
  fn = fname;
  if (!f) throw ArchExc("file not found");
  fseek(f,0,SEEK_END);
  fLen = ftell(f);
  curPos = 0;
  curKind = -1;
  curRawKind = -1;
  blockSet = new BlockSet;
  utcOffset=2;
  computeMJD(fname);
}

ArcheopsFile::ArcheopsFile(ArcheopsFile const& x) {
  blockSet = x.blockSet ? new BlockSet(*x.blockSet) : NULL;
  fposStack = x.fposStack;
  
  stack<BlockSet*> tmp;
  ArcheopsFile& y = (ArcheopsFile&) x;
  while (!y.blockStack.empty()) {
    tmp.push(y.blockStack.top());
    y.blockStack.pop();
  }
  while (!tmp.empty()) {
    y.blockStack.push(tmp.top());
    blockStack.push(new BlockSet(*tmp.top()));
    tmp.pop();
  }

  curBlock   = x.curBlock;
  curKind    = x.curKind;
  curRawKind = x.curRawKind;
  curPos     = x.curPos;
  fLen       = x.fLen;
  fn         = x.fn;
  f = fopen(fn.c_str(), "rb");
}

ArcheopsFile::~ArcheopsFile() {
  if (f) fclose(f);
  delete blockSet;
  while (!blockStack.empty())
    { delete blockStack.top(); blockStack.pop();}
}

void ArcheopsFile::grabLastBlocs(ArcheopsFile const& oldFile) {
  delete blockSet;
  blockSet = new BlockSet(*oldFile.blockSet);
  setUTCOffset(oldFile.utcOffset);
}


def_nom_block

bool ArcheopsFile::nextBlock() {
  //if (curPos>0) saveCurBlock();
  //if (curPos<0) curPos = 0;
  if (curPos+12 > fLen) return false;
  fseek(f,curPos,SEEK_SET);
  size_t read = fread(&curBlock,1,sizeof(curBlock),f);
  if (read < longueur_block(&curBlock)) return false;
  if (verifie_block(&curBlock) != block_correct) {
  printf("block invalide\n"); throw ArchExc("invalid block");
  }
  curRawKind = curKind = type_block(&curBlock);
  curPos += longueur_block(&curBlock);
  //printf("block %d  :  %s \n",numero_block(&curBlock),nom_block[curKind]);
  postProcessBlock();
  saveCurBlock();
  return true;
}

bool ArcheopsFile::nextBlock(long mask) {
  if (!mask) return false;
  while (1) {
    if (!nextBlock()) return false;
    if (( 1 << curRawKind) & mask) return true;
  }
}

int  ArcheopsFile::blockKind() {
  return curKind;
}

int  ArcheopsFile::blockRawKind() {
  return curRawKind;
}

int  ArcheopsFile::blockNum() {
  return numero_block(&curBlock);
}

#define  bitmot                 24              // nb de bit horloge dans un mot ADC

double ArcheopsFile::perEchant() { // periode d'echantillonage pour le dernier bloc reglage
  double p,f1,f2,f3;
  int pp;
  pp=lastReglage()->reglage.horloge.periode;
  p=pp/5.;
  f1=1000/p;f2=f1/bitmot;f3=f2*1000./(double)(lastReglage()->reglage.horloge.nb_mesures);
  return 0.5/f3;                                        //  2 fois la frequence de modulation
}

double ArcheopsFile::perBlock() {  // duree (en secondes) correspondant a un bloc bolo
  return perEchant() * (double)lastParam()->param.n_per_block*2.;
}

int ArcheopsFile::nEchBlock() {  // Nb d'echantillons dans un bloc
  return lastParam()->param.n_per_block*2.;
}

string  ArcheopsFile::blockKdName() {
  return nom_block[curKind];
}

string  ArcheopsFile::blockRawKdName() {
  return nom_block[curRawKind];
}


block_type_modele*  ArcheopsFile::currentBlock() {
  if (curPos<0) return NULL;
  return &curBlock;
}

block_type_param* ArcheopsFile::lastParam() {
  return blockSet->lastParam;
}
block_type_journal* ArcheopsFile::lastJournal() {
  return blockSet->lastJournal;
}
block_type_reglage* ArcheopsFile::lastReglage() {
  return blockSet->lastReglage;
}
block_type_dilution* ArcheopsFile::lastDilution() {
  return blockSet->lastDilution;
}
block_type_gps* ArcheopsFile::lastGPS() {
  return blockSet->lastGPS;
}
block_type_une_periode* ArcheopsFile::lastUnePeriode() {
  return blockSet->lastUnePeriode;
}
block_type_synchro_sol* ArcheopsFile::lastSynchroSol() {
  return blockSet->lastSynchroSol;
}
block_type_pointage_sol* ArcheopsFile::lastPointageSol() {
  return blockSet->lastPointageSol;
}
block_type_bolo* ArcheopsFile::lastBolo() {
  return blockSet->lastBolo;
} 
block_type_bolo* ArcheopsFile::llastBolo() {
  return blockSet->llastBolo;
} 
block_type_gyro* ArcheopsFile::lastGyro() {
  return blockSet->lastGyro;
}
block_type_sst* ArcheopsFile::lastSST() {
  return blockSet->lastSST;
}
block_type_bolo_comprime* ArcheopsFile::lastBoloComp() {
  return blockSet->lastBoloComp;
} 
block_type_gyro_comprime* ArcheopsFile::lastGyroComp() {
  return blockSet->lastGyroComp;
}
block_type_sst_comprime* ArcheopsFile::lastSSTComp() {
  return blockSet->lastSSTComp;
}

void ArcheopsFile::postProcessBlock() {
  switch (curKind) {
     case block_bolo_comprime: {
       blockSet->setRawBloc(curBlock);
       block_type_bolo blk2;
       block_type_bolo_comprime* blk = (block_type_bolo_comprime*) &curBlock;
       for(int j=0;j<nb_bolo_util;j++)  
                {
                   decompress_7_2(blk->data_bolo[j],blk2.data_bolo[j],nb_per_block*2);
                }
       valide_block((block_type_modele*)&blk2,block_bolo,numero_block(blk));    
       memcpy(&curBlock,&blk2,sizeof(blk2));
       curKind = block_bolo;
       break;
     }
  }
}

void ArcheopsFile::saveCurBlock() {
  blockSet->setBloc(curBlock);
}

void ArcheopsFile::pushMark() {    // push current file position, and "last" blocks`
  fposStack.push(curPos);
  blockStack.push(new BlockSet(*blockSet));
}

void ArcheopsFile::popMark() {    // pops last file position and "last" blocks
  if (! blockStack.empty()) {
    delete blockSet;
    blockSet = blockStack.top();
    curPos = fposStack.top();
    blockStack.pop();
    fposStack.pop();
  }
}

static int mlen[] = {31,28,31,30,31,30,31,31,30,31,30,31};

void ArcheopsFile::computeMJD(string const& fname) {
  //h99_04_29-15h36m22
  short y,m,d,hh,mm,ss;
  char const* p = fname.c_str();
  char const* p2 = p  + fname.length()-1;
  while (*p2 != ':' && *p2 != '/' && p2 != p) p2--;
  if (*p2 == ':' || *p2 == '/') p2++;
  if (*p2 == 'h') p2++;
  y  = atoi(p2); p2+=3;
  m  = atoi(p2); p2+=3;
  d  = atoi(p2); p2+=3;
  hh = atoi(p2); p2+=3;
  mm = atoi(p2); p2+=3;
  ss = atoi(p2);
  
  if (y<50) y += 100;
  // 1. depuis 0/1/97 minuit
  rawMJD = (int) (365.25 * (y-97));
  for (int i=0; i<m-1; i++)
    rawMJD += mlen[i];
  if (y%4 == 0 && m > 2) rawMJD++;
  rawMJD += d;
  rawMJD += hh/24. + mm/24./60. + ss/24./60./60;
  
  rawMJD += 448.5; // 0/1/97 minuit
  startMJD = rawMJD - utcOffset/24.;
}

void ArcheopsFile::setUTCOffset(int UTCOffset) {
  utcOffset = UTCOffset;
  startMJD = rawMJD - utcOffset/24.;
}

double ArcheopsFile::getStartMJD() {
  return startMJD;
}

// GPS
// $GPGGA,hhmmss.ss,ddmm.mmmm,n,dddmm.mmmm,e,q,ss,y.y,a.a,z,


int ArcheopsFile::getGPSBlockNum() {
  if (!lastGPS()) return -1;
  return numero_block(lastGPS());
}

double ArcheopsFile::getGPSUTC() {  // en secondes depuis minuit UTC jour courant
  if (!lastGPS()) return -1;
  char* p = lastGPS()->gps;
  if (strncmp(p, "$GPGGA,", 7)) return -1;
  p += 7;
  double t;
  sscanf(p, "%lg", &t);
  int h = int(t/10000); t -= h*10000;
  int m = int(t/100);   t -= m*100;
  return h*3600. + m*60. + t;
}

double ArcheopsFile::getGPSMJD() {  // modified julian day du dernier bloc GPS, JD - 2450000
  double t = getGPSUTC()/86400.;
  if (t<0) return t;
  if (t > (startMJD - int(startMJD))) {
    return int(startMJD) + t;
  } else {
    return int(startMJD) + 1. + t;
  }
}

double ArcheopsFile::getGPSLat() {  // degres, +  = NORD
  if (!lastGPS()) return 99999;
  char* p = lastGPS()->gps;
  if (strncmp(p, "$GPGGA,", 7)) return 99999;
  char* fence = p+80;
  p += 7;
  while (*p != ',' && p<fence) p++;
  if (*p != ',') return 99999; 
  p++;
  double t;
  sscanf(p, "%lg", &t);
  int d = int(t/100); t -= d*100;
  t = d + t/60;
  while (*p != ',' && p<fence) p++;
  if (*p != ',') return 99999; 
  p++;
  if (*p == 'S') t = -t;
  return t;
}

double ArcheopsFile::getGPSLong() { // degres, +  = EST
  if (!lastGPS()) return -1;
  char* p = lastGPS()->gps;
  if (strncmp(p, "$GPGGA,", 7)) return 99999;
  char* fence = p+80;
  p += 7;
  for (int i=0; i<3; i++) {
    while (*p != ',' && p<fence) p++;
    if (*p != ',') return 99999; 
    p++;
  }
  double t;
  sscanf(p, "%lg", &t);
  int d = int(t/100); t -= d*100;
  t = d + t/60;
  while (*p != ',' && p<fence) p++;
  if (*p != ',') return 99999; 
  p++;
  if (*p == 'W') t = -t;
  return t;
}



// Bolo

long ArcheopsFile::getRawBolo(int ibolo, int imesure) {    // donnee brute, avec seulement soustraction offset
  int nb_coups,aa;  
  block_type_bolo* blk = imesure >= 0 ? lastBolo() : llastBolo();
  if (imesure < 0) imesure += nEchBlock();
  if (!blk) return 0;
  block_type_reglage* reglage = lastReglage();
  block_type_param* param = lastParam();
  
  nb_coups= reglage->reglage.horloge.nb_mesures/2 - reglage->reglage.horloge.temp_mort;
  aa = (nb_coups<<14) + (nb_coups*190) ;
  
  int s = imesure % 2 ? 1 : -1;

  return s*(((blk->data_bolo[ibolo][imesure]-aa)<<1)/nb_coups);
}

def_gains 

double ArcheopsFile::getMuVBolo(int ibolo, int imesure) { // microvolts, filtre avec filtre carre
  double y = (getRawBolo(ibolo, imesure-1) + getRawBolo(ibolo, imesure))/2.;
  block_type_reglage* reglage = lastReglage();
  return ((1e5*y)/(65536.*gain(reglage->reglage.bolo[ibolo])));
}


// SST, gyros...
// $CHECK$ TBD


double ArcheopsFile::getAzimut(int imesure) {return imesure*360/nEchBlock();}
double ArcheopsFile::getPendDirect(int /*imesure*/) {return 0;}
double ArcheopsFile::getPendOrth(int /*imesure*/) {return 0;}
double ArcheopsFile::getAlpha(int /*imesure*/) {return 0;}
double ArcheopsFile::getDelta(int /*imesure*/) {return 0;}