source: Sophya/trunk/SophyaProg/Tests/tfitsdt.cc@ 4036

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <iostream>

#include "sopnamsp.h"
#include "histinit.h"
#include "datatable.h"
#include "fitshdtable.h"
#include "swfitsdtable.h"
#include "swppfdtable.h"
#include "ctimer.h"

/* ----------------------------------------------------------
   Programme de test lecture-ecriture FITS de DataTable
    + test classe SwFitsDataTable
      Oct 2005 - Jan 2006

   --- This code is part of the SOPHYA library ---
  (C) Univ. Paris-Sud   (C) LAL-IN2P3/CNRS   (C) IRFU-CEA
  (C) R. Ansari, C.Magneville    2005 - 2010   

    ----------------------------------------------------------    */

//-----------------------------------------------------------
int test_fdtable() ; 
void tdfill_dtable(BaseDataTable& dtr, uint_8 sz);
void tdcheck_dtable(BaseDataTable& dtr, uint_8& errcntI, uint_8& errcntF);
int test_memdtable(uint_8 sz, bool fgppf=true) ; 
int test_swppfdtable(uint_8 sz) ; 
int test_swfitsdtable(uint_8 sz) ; 
#define _TDTSEGSIZ_  8192
//------------------------------------------------------------

//--------------   MAIN PROGRAM  ------------ 
int main(int narg, char *arg[])
{
  SophyaInit();
// 0 -> test_fdtable , 1 -> test_memdtable/PPF , 2 ->  test_memdtable/FITS 
// 3 -> test_swppfdtable , 4 -> test_swfitsdtable
  int opt = 0;   
  if (narg<2) {
    cout << " Usage: tfitsdt OPT=0/1/2/3/4 [ Size=NRows in units of 1000] \n" 
         << "  0 -> test_fdtable \n  1 -> test_memdtable/PPF , 2 ->  test_memdtable/FITS \n" 
         << "  3 -> test_swppfdtable , 4 -> test_swfitsdtable \n" 
         << "  Rc=0 OK for 1/2/3/4 , Def.Size=2048 (*1000) " << endl;
    return 1; 
  }
  opt = atoi(arg[1]);
  uint_8 sz = 2048*1000;
  if (narg>2)  sz = (uint_8)(atof(arg[2])*1000.);  
  int rc=0;
  try {
    if (opt==1) rc=test_memdtable(sz,true);
    else if (opt==2) rc=test_memdtable(sz,false);
    else if (opt==3) rc=test_swppfdtable(sz);
    else if (opt==4) rc=test_swfitsdtable(sz);
    else test_fdtable();
  }
  catch(PThrowable exc ) {
    cerr << "tfitsdt-main() , Catched exception: \n" << exc.Msg() << endl;
  }
  catch(std::exception ex) {
    cerr << "tfitsdt-main() , Catched exception ! " << (string)(ex.what()) << endl;
  }
  catch(...) {
    cerr << "tfitsdt-main() , Catched ... ! " << endl;
  }
  cout << " ========== End of tfitsdt.cc , Rc= " << rc << " ============= " << endl;
  return rc;
}


/* --Nouvelle-Fonction-- */
int test_fdtable() 
{
  int NL = 600;
  DataTable refdt(100); // Reference DataTable - for tests 

  cout << "======= test_fdtable: Simple DataTable+FITS test ======= " << endl;
  {
    cout << "1/ Creating DataTable / Writing to FITS " << endl;
    DataTable dt(64);
    dt.AddIntegerColumn("line");
    dt.AddDoubleColumn("x");
    dt.AddFloatColumn("f_sin");
    dt.AddDoubleColumn("f_x2");
    dt.AddStringColumn("str_line");
    dt.AddComplexColumn("cmplx_cos_sin");
    MuTyV rec[10];
    char sbuff[32];
    cout << " Filling ... (NLines=" << NL << ")" << endl;
    for(int k = 0; k<NL; k++) {
      rec[0] = k;
      double x = M_PI*k/100.;
      double fx = sin(x)*cos(x);
      rec[1] = x+1;
      rec[2] = sin(x);
      rec[3] = x*x;
      sprintf(sbuff, "SL-%d", k);
      rec[4] = sbuff;
      double sx = sin(x);
      double cx = cos(x);
      rec[5] = complex<r_4>(cx, sx);
      dt.AddLine(rec);
    }
    cout << dt;
    cout << " Copying dt to refdt (refdt = dt) " << endl;
    refdt = dt;
    cout << " Writing dt to fits file dtable.fits ... " << endl;
    FitsInOutFile fios("!dtable.fits", FitsInOutFile::Fits_Create);
    fios << dt;
    cout << " Writing dt to fits as ASCII table with extname= ASC_DTable ... " << endl;
    fios.SetDef_AscTable();
    fios.SetNextExtensionName("ASC_DTable");
    fios << dt;    
    cout << "1.b/ Creating SwFitsDataTable (file swdtable.fits) " << endl;
    FitsInOutFile swf("!swdtable.fits", FitsInOutFile::Fits_Create);
    SwFitsDataTable swdt(swf, 64);
    cout << " Copying from DataTable dt ..." << endl;
    swdt = dt;
    cout << swdt;
  }
  {
    cout << "2/ Reading DataTable from FITS " << endl;
    FitsInOutFile fios("dtable.fits", FitsInOutFile::Fits_RO);
    cout << fios ;
    fios.MoveToNextHDU();
    DataTable dtr;
    fios >> dtr;
    cout << dtr;
    cout << "2.b/  dtr.LineHeaderToString() dtr.LineToString(k)   : " << endl;
    cout << dtr.LineHeaderToString() << endl;
    for(int k = 0; k<NL; k+=NL/12) 
      cout << "Line[" << k << "] " << dtr.LineToString(k) << endl ;  
    cout << "2.b/  Reading from ASCII table   " << endl;
    DataTable dtra;
    fios >> dtra;
    cout << dtra;
    cout << dtra.LineHeaderToString() << endl;
    for(int k = 0; k<NL; k+=NL/12) 
      cout << "Line[" << k << "] " << dtra.LineToString(k) << endl ;  

    cout << "2.c/  Reading SwFitsDataTable from  swdtable.fits  " << endl;   
    SwFitsDataTable swdtr;
    FitsInOutFile swf("swdtable.fits", FitsInOutFile::Fits_RO);   
    swf.MoveAbsToHDU(2);  
    swf >> swdtr;
    cout << swdtr;
    cout << swdtr.LineHeaderToString() << endl;
    for(int k = 0; k<NL; k+=NL/12) 
      cout << "Line[" << k << "] " << swdtr.LineToString(k) << endl ;  
  }
  cout << "============ FIN  test_fdtable =============== " << endl;  
  return 0;
}



/* --Nouvelle-Fonction-- */
void tdfill_dtable(BaseDataTable& dt, uint_8 sz)
{
  dt.AddIntegerColumn("ic1");
  dt.AddFloatColumn("fc2");
  int ic1,ric1;
  float fc2,rfc2;
  MuTyV rec[5];
  uint_8 modprt=sz/10;
  for(uint_8 k=0; k<sz; k++) {
    ic1 = (int)(k%1000000000LL);
    fc2 = (float)ic1+0.5;  fc2 *= (fc2+0.2); 
    rec[0]=ic1;  rec[1]=fc2;
    dt.AddLine(rec);
    if (k%modprt==0)  cout << " tdfill_dtable/ done k=" << k << " : max=" << sz << endl; 
  }
  return;
}

/* --Nouvelle-Fonction-- */
void tdcheck_dtable(BaseDataTable& dtr, uint_8& errcntI, uint_8& errcntF)
{
  errcntI=0;
  errcntF=0;
  int ic1,ric1;
  float fc2,rfc2;
  MuTyV rec[5];
  uint_8 modprt=dtr.NRows()/10;

  for(sa_size_t k=0; k<dtr.NRows(); k++) {
    ic1 = (int)(k%1000000000LL);
    fc2 = (float)ic1+0.5;  fc2 *= (fc2+0.2); 
    dtr.GetRow(k,rec);
    ric1=rec[0];  rfc2=rec[1];
    if (ric1!=ic1)  errcntI++;
    if (fabs(fc2-rfc2)>1.e-9) errcntF++;
    if (k%modprt==0)  
      cout << " tdcheck_dtable/ done k=" << k << " : max=" << dtr.NRows() 
           << " ErrCntI=" << errcntI << " ErrCntF=" << errcntF << endl; 
  }
  return;
}

/* --Nouvelle-Fonction-- */
int test_memdtable(uint_8 sz, bool fgppf) 
{
  Timer tm("Tst_MEM_DataTable");  
  cout << " ------------------- Tst_MEM_DataTable -------------------------- " << endl;
  cout << " --- test_memdtable[1]: creating DataTable with " << sz << " rows (=" << sz/1000 << " kilo) ... " << endl;
  string dtfilename="tstmemdtable.ppf";
  {
    DataTable dt(_TDTSEGSIZ_);
    tdfill_dtable(dt, sz);
    cout << dt;
    tm.Split();
    if (fgppf)   {
      POutPersist po(dtfilename);
      po << dt;
      cout << " --- test_memdtable[2]: datatable written to PPF file " << dtfilename << endl;
    }
    else {
      dtfilename="!tstmemdtable.fits";
      FitsInOutFile fios(dtfilename, FitsInOutFile::Fits_Create);
      fios << dt;
      cout << " --- test_memdtable[2]: datatable written to FITS file " << dtfilename << endl;
    }
  }

  uint_8 errcntI=0;
  uint_8 errcntF=0;
  {
    DataTable dtr;
    if (fgppf)   {
      cout << " --- test_memdtable[3]: reading DataTable from PPF file " << dtfilename << endl;
      PInPersist pin(dtfilename);
      pin >> dtr;
    }
    else {
      dtfilename="tstmemdtable.fits";
      FitsInOutFile fios(dtfilename, FitsInOutFile::Fits_RO);
      fios.MoveToNextHDU();
      cout << " --- test_memdtable[3]: reading DataTable from FITS file " << dtfilename << endl;
      fios >> dtr;
    }
    tm.Split();
    tdcheck_dtable(dtr, errcntI, errcntF);
  } 
  cout << " --- test_memdtable[4]: End Read/Check ErrCntI=" << errcntI << " ErrCntF=" << errcntF << endl; 
  if ((errcntI==0)&&(errcntF==0))  return 0;
  else if (errcntI==0)  return 5;
  else return 9;
}

/* --Nouvelle-Fonction-- */
int test_swppfdtable(uint_8 sz) 
{
  Timer tm("Tst_SwPPF_DataTable");  
  cout << " ------------------- Tst_SwPPF_DataTable ------------------------- " << endl;
  cout << " --- test_swppfdtable[1]: creating DataTable with " << sz << " rows (=" << sz/1000 << " kilo) ... " << endl;
  string dtswfilename="tstswppfdtable.ppf";
  {
    POutPersist po(dtswfilename);
    SwPPFDataTable dt(po, _TDTSEGSIZ_);
    tdfill_dtable(dt, sz);    
    po << dt;
    cout << " --- test_swppfdtable[2]: SwPPFDataTable written to PPF file " << dtswfilename << endl;
    tm.Split();
  }

  uint_8 errcntI=0;
  uint_8 errcntF=0;
  {
    PInPersist pin(dtswfilename);
    SwPPFDataTable dtr;
    cout << " --- test_swppfdtable[3]: reading SwPPFDataTable from PPF file " << dtswfilename << endl;
    pin >> dtr;
    tdcheck_dtable(dtr, errcntI, errcntF);
  } 
  cout << " --- test_swppfdtable[4]: End Read/Check ErrCntI=" << errcntI << " ErrCntF=" << errcntF << endl; 
  if ((errcntI==0)&&(errcntF==0))  return 0;
  else if (errcntI==0)  return 5;
  else return 9;
}

/* --Nouvelle-Fonction-- */ 
int test_swfitsdtable(uint_8 sz) 
{
  Timer tm("Tst_SwFITS_DataTable");  
  cout << " ------------------- Tst_SwFITS_DataTable ----------------------- " << endl;
  cout << " --- test_swfitstable[1]: creating SwFitsDataTable with " << sz << " rows (=" << sz/1000 << " kilo) ... " << endl;
  string dtswfilename="!tstswfitsdtable.fits";
  {
    FitsInOutFile swf(dtswfilename, FitsInOutFile::Fits_Create);
    SwFitsDataTable swdt(swf, _TDTSEGSIZ_);
    tdfill_dtable(swdt, sz);    
    swf << swdt;
    cout << " --- test_swfitstable[2]: SwFitsDataTable written to FITS file " << dtswfilename << endl;
    tm.Split();
  }

  uint_8 errcntI=0;
  uint_8 errcntF=0;
  {
    dtswfilename="tstswfitsdtable.fits";
    SwFitsDataTable swdtr;
    FitsInOutFile swif(dtswfilename, FitsInOutFile::Fits_RO);   
    swif.MoveAbsToHDU(2);  
    cout << " --- test_swfitstable[3]: reading SwFitsDataTable from FITS file " << dtswfilename << endl;
    swif >> swdtr;
    tdcheck_dtable(swdtr, errcntI, errcntF);
  } 
  cout << " --- test_swfitstable[4]: End Read/Check ErrCntI=" << errcntI << " ErrCntF=" << errcntF << endl; 
  if ((errcntI==0)&&(errcntF==0))  return 0;
  else if (errcntI==0)  return 5;
  else return 9;
}