source: Sophya/trunk/SophyaProg/Tests/tnfits.cc@ 3094

#include "sopnamsp.h"
#include "machdefs.h"
#include <stdio.h>
#include <stdlib.h>
#include "math.h"
#include <iostream>
#include <string>

#include <typeinfo>
#include "timing.h"
#include "histinit.h"
#include "array.h"
#include "fitsinoutfile.h"
#include "fitsblkrw.h"
#include "fitsarrhand.h"
#include "fiosinit.h"

/*  Programme test des classes NFits                  */
/*  SOPHYA - R. Ansari (LAL)  -   Sep 2005           */

int tsfitsReadFile(string& flnm);

int main(int narg, char* arg[])
{
  cout << " ---- Programme tnfits - Test classes NFits --- " << endl;
  if ( (narg > 1) && (strcmp(arg[1],"-h")==0) ) {
    cout << " Usage : tnfits [InputFitsFileName] " << endl;
    return 1;
  }
  int rc = 0;
  try {
    SophyaInit();  
    FitsIOServerInit();
    InitTim(); 
    if (narg > 1) {   // Lecture fichier entree 
      string flnm = arg[1];
      cout << "1/=== Lecture fichier FITS : " << flnm << endl;
      FitsInOutFile fios(flnm, FitsInOutFile::Fits_RO);
      fios.Print();
      bool encore = true;
      while (encore) {
        cout << " HDU No" << fios.CurrentHDU() << " Type= " 
             << fios.CurrentHDUTypeStr() << endl;
        DVList dvl;
        fios.GetHeaderRecords(dvl);
        cout << dvl;
        if (fios.CurrentHDUType() == IMAGE_HDU) {
          long naxes[5] = {0,0,0,0,0};
          int naxis=5;
          fios.GetImageHDUInfo(naxis, naxes);
          cout << "--GetImageHDUInfo() naxis= " << naxis << " : "
               << naxes[0] << 'x' << naxes[1] << 'x' << naxes[2] << endl;
          if (naxis == 2) {
            cout << "---Reading Matrix ... mx= " << endl;
            Matrix mx(naxes[1], naxes[0]);
            FitsBlockRW<double>::ReadImageData(fios, mx.Data(), mx.Size());
            cout << mx;
          }
        }
        else {
          vector<string> colnames;
          vector<int> coltypes;
          vector<long> repcnt;
          vector<long> width;
          int ncols = fios.GetColInfo(colnames, coltypes, repcnt, width);
          cout << " NCols from GetColInfo = " << ncols << endl;
          for(int kk=0; kk<colnames.size(); kk++) {
            cout << " kk=" << kk << " ColName= " << colnames[kk] 
                 << " Type= " << coltypes[kk] 
                 << " Repeat= " << repcnt[kk]
                 << " W= " << width[kk] << endl;
          }
        }
        if (fios.CurrentHDU() < fios.NbHDUs()) 
          fios.MoveRelToHDU(1);
        else encore = false;
      }
    }  // ---- Fin lecture fichier entree 
    {
      cout << "2/=== Ecriture BinTable dans toto.fits  " << endl;
      FitsInOutFile fiosc("!toto.fits", FitsInOutFile::Fits_Create);
      //     cout << " CurrHDU-1 : " << fiosc.CurrentHDU() 
      //           << " Type=" << fiosc.CurrentHDUType() <<endl;
      long naxes[5] = {5,5,0,0,0};
      fiosc.CreateImageHDU(FLOAT_IMG, 2, naxes);
      float data[100];
      double data2[100];
      string sdata[100];
      MuTyV mtv;
      for(int kk=0; kk<100; kk++) { 
        data[kk] = kk/5.;
        data2[kk] = cos(data[kk]);
        mtv = kk;
        sdata[kk] = "Str_KK=";
        sdata[kk] += (string)(mtv);
      }
      cout << " CurrHDU-2 : " << fiosc.CurrentHDU() 
           << " Type=" << fiosc.CurrentHDUType() <<endl;
      FitsBlockRW<float>::WriteImageData(fiosc, data, 25);
      fiosc.Print(); 
      /*
      char * ttype[2] = {"Col1", "Col2"};
      char * tform[2] = {"E", "D"};
      char * tunit[2] = {"", ""};
      fiosc.CreateTable(BINARY_TBL, "rzbintbl", 2, ttype, tform, tunit);
      */
      vector<string> ttype, tform, tunit;
      ttype.push_back("Col1");
      ttype.push_back("Col2");
      ttype.push_back("Col3");
      tform.push_back("E");
      tform.push_back("D");
      tform.push_back("18A");
      tunit.push_back("");
      tunit.push_back("");
      tunit.push_back("");
      string tblnm = "rzbintbl";
      fiosc.CreateTable(BINARY_TBL, "rzbintbl", ttype, tform, tunit);

      //      cout << " MovRelToHDU , Type=" << fiosc.MoveRelToHDU(1) << endl;
      cout << " CurrHDU-3 : " << fiosc.CurrentHDU() 
           << " Type=" << fiosc.CurrentHDUType() <<endl;
      fiosc.Print(); 
      /*
      cout << " DBG- fios.InsertColumn(1,Col1,E) " << endl;
      fios.InsertColumn(1,"Col1","E");
      cout << " DBG- fios.InsertColumn(2,Col2,D) " << endl;
      fios.InsertColumn(2,"Col2","D");
      */
      FitsBlockRW<float>::WriteColumnData(fiosc, 1, 1, 1, data, 100);
      FitsBlockRW<double>::WriteColumnData(fiosc, 2, 1, 1, data2, 100);
      FitsBlockRW<string>::WriteColumnData(fiosc, 3, 1, 1, sdata, 100);
            
      // Ecriture entete
      DVList  dvl;
      dvl("toto") = 14;    dvl("titi") = 25.5;    
      dvl("tata") = "Bonjour !"; 
      dvl("hello") = 88;  dvl("Hello") = 77.77;    
      dvl.Comment() = "Blabla ... zone commentaire ds DVList";
      fiosc.WriteHeaderRecords(dvl);
    }

    {
      cout << "3/=== Ecriture tableaux dans ta.fits  " << endl;
      FitsInOutFile fiosc("!ta.fits", FitsInOutFile::Fits_Create);
      cout << " TArray<r_4> mxf(100,50,3) , TMatrix<int_4> mxi(60,30) " << endl;
      TArray<r_4> mxf(100,50,3);
      mxf = RegularSequence(0., 0.2);
      mxf.Info()["info1"] = "info1: TArray<r_4>";
      mxf.Info()["info2"] = "info2: (100,50,3)";
      mxf.Info()["info3 "] = "info3: mot-cle-avec-blanc";
      cout << mxf.Info();
      TMatrix<int_4> mxi(60,30);
      mxi = RegularSequence();
      mxi.Info()["MtxNRow"] = mxi.NRows();
      mxi.Info()["MtxNCol"] = mxi.NCols();
      mxi.Info()["MtxStr"] = "TMatrix<int_4> mxi(60,30)";
      fiosc << mxf << mxi ;
    }
    {
      cout << "4/=== Lecture tableaux depuis ta.fits  " << endl;
      FitsInOutFile fiosc("ta.fits", FitsInOutFile::Fits_RO);
      TArray<r_4> mxf;
      TMatrix<int_4> mxi;
      fiosc >> mxf >> mxi ;
      cout << mxf;
      cout << " mxf.Info()= " << mxf.Info() << mxi;
      cout << " mxi.Info()= " << mxi.Info() << mxi;
    }

  }
  catch (PThrowable & exc) {
    cerr << " tnfits.cc: Catched Exception (PThrowable)" << (string)typeid(exc).name() 
         << "\n...exc.Msg= " << exc.Msg() << endl;
    rc = 99;
  }
  catch (std::exception & e) {
    cerr << " tnfits.cc: Catched std::exception "  
         << " - what()= " << e.what() << endl;
    rc = 98;
  }
  catch (...) {
    cerr << " tnfits.cc: some other exception (...) was caught ! " << endl;
    rc = 97;
  }
  PrtTim("End tnfits " );
  cout << " ---- Programme tnfits.cc-  FIN  (Rc=" << rc << ") --- " << endl;
  return rc;
}