source: Sophya/trunk/SophyaExt/FitsIOServer/datacirclefits.cc@ 4050

#include <iostream>
#include "sopnamsp.h"
#include "datacirclefits.h"

/*
   Class used to recover the circles and measures involved in a mission. The parameters needful to create a circle and the number of samples on this circle are read from a FITS data array. In addition, a function returns the measurement at a given angular position psi on the circle.
*/

DataCircleFits::DataCircleFits() : 
  _fptr(NULL)
, _NMeasurements(0)
, _ICircle(0)
, _mesures(NULL)
, _stored_data(false)
{;}

DataCircleFits::DataCircleFits(fitsfile *is,int hdunum,bool sdata) :
  _fptr(NULL)
, _NMeasurements(0)
, _ICircle(hdunum)
, _mesures(NULL)
, _stored_data(sdata)
{

  // pointer to the FITS file
  _fptr= is;
  int status = 0;

  // move to the HDU containing a circle
  fits_movabs_hdu(_fptr,hdunum,NULL,&status);

  // angles of the circle
  double theta;
  fits_read_key(_fptr,TDOUBLE,"CIRTHETA",&theta,NULL,&status);

  double phi;
  fits_read_key(_fptr,TDOUBLE,"CIRPHI",&phi,NULL,&status);

  double aperture;
  fits_read_key(_fptr,TDOUBLE,"CIRAPER",&aperture,NULL,&status);

  UnitVector temp(theta,phi);
  SetCircle(temp,aperture);

  fits_read_key(_fptr,TINT,"NSAMPLES",&_NMeasurements,NULL,&status);

  if(_stored_data) {
    // vector repeat value for the column (there is only 1 column)
    long repeat;
    fits_get_coltype(_fptr,1,NULL,&repeat,NULL,&status);

    double dnull= 0.;
    int    anull;
    _mesures= new double[_NMeasurements];
   int i;
   for(i = 0; i < _NMeasurements; i++) {
     // row index: range 1 to NAXIS2
     int index= i+1;
     
     // starting line index
     int stline= index/repeat; 
     int stelem= index-stline*repeat;
     if(stelem != 0) {
       stline++;
     } else {
       stelem= repeat;
     }
     fits_read_col(_fptr,TDOUBLE,1,stline,stelem,1,&dnull,&_mesures[i],&anull,&status);
   }
  }
}
  
DataCircleFits::~DataCircleFits() {

  if(_mesures != NULL) delete [] _mesures;
}

int DataCircleFits::NMeasurements() const {

  return _NMeasurements;
}

double DataCircleFits::getData(double psi) const {

  double dtab;
  int ibin= l_ft_nint(psi*_NMeasurements/(2.*M_PI));
  //cout << " ibin= " << ibin << ", " << psi << endl;

  if(_stored_data) {
    dtab= _mesures[ibin];
  } else {
    int status= 0;
    fits_movabs_hdu(_fptr,_ICircle,NULL,&status);

    // vector repeat value for the column (there is only 1 column)
    long repeat;
    fits_get_coltype(_fptr,1,NULL,&repeat,NULL,&status);
   
    double dnull= 0.;
    int    anull;

    int index= ibin+1;
    // starting line index
    int stline= index/repeat; 
    int stelem= index-stline*repeat;
    if(stelem != 0) {
      stline++;
    } else {
      stelem= repeat;
    }
    fits_read_col(_fptr,TDOUBLE,1,stline,stelem,1,&dnull,&dtab,&anull,&status);
  }
  //cout << "DataCircleFits:: bin= " << ibin << ", " << dtab << endl;
  return dtab;
}

void DataCircleFits::print(ostream& out) const {

  out << " Circle:: ApertureAngle= " << ApertureAngle() << ", Theta= " 
      << Theta() << ", Phi= " << Phi() << ", NSamples= " 
      << NMeasurements() << endl;
}

double DataCircleFits::getTMeasure(double psi) const {

  if(_mesures == NULL) {
    cout << "  DataCircleFits::getTMeasure data must be stored"
         << " in an array... verify the option" << endl;
    exit(0);
  }
  int ibin= l_ft_nint(psi*_NMeasurements/(2.*M_PI));
  return _mesures[ibin];
}

double DataCircleFits::getTOffset() const {

  int status= 0;
  fits_movabs_hdu(_fptr,_ICircle,NULL,&status);

  double omes= 0.0;
  fits_read_key_dbl(_fptr,"OFFSET",&omes,NULL,&status);
  return omes;
}