#include "machdefs.h" #include #include "fitsfile.h" #include "pexceptions.h" #include "strutil.h" #include "anydataobj.h" #include "fitsspherehealpix.h" FitsFile::FitsFile() { fptr_= NULL; hdutype_= 0; hdunum_ = 0; nrows_ = 0; } FitsFile::~FitsFile() { int status = 0; if( fptr_ != NULL) { fits_close_file(fptr_,&status); delete fptr_; } if( status ) printerror( status ); } int FitsFile::NbBlocks(char flnm[]) { int status = 0; int nbhdu = 0; fitsfile* fileptr; fits_open_file(&fileptr,flnm,READONLY,&status); if( status ) printerror( status, "NbBlocks: erreur ouverture fichier" ); fits_get_num_hdus(fileptr, &nbhdu, &status); fits_close_file(fileptr,&status); return nbhdu; } void FitsFile::getBlockType(char flnm[], int hdunum, string& typeOfExtension, int& naxis, vector& naxisn, string& dataType, DVList& dvl ) { int status = 0; fitsfile* fileptr; fits_open_file(&fileptr,flnm,READONLY,&status); if( status ) printerror( status, "getBlockType: erreur ouverture fichier" ); // move to the specified HDU number int hdutype = 0; fits_movabs_hdu(fileptr,hdunum,&hdutype,&status); if( status ) printerror( status,"getBlockType: erreur movabs"); if(hdutype == IMAGE_HDU) { typeOfExtension = "IMAGE"; int bitpix; GetImageParameters (fileptr, bitpix, naxis, naxisn); if(bitpix == DOUBLE_IMG) dataType = "double"; else if(bitpix == FLOAT_IMG) dataType = "float"; else if(bitpix == LONG_IMG || bitpix == SHORT_IMG ) dataType = "int"; else { cout << " bitpix= " << bitpix << endl; throw PException(" FitsFile::getBlockType : unsupprted FITS data type"); } } else if(hdutype == ASCII_TBL || hdutype == BINARY_TBL) { int nrows = 0; vector noms; vector types; vector taille_des_chaines; GetBinTabParameters(fileptr, naxis, nrows, naxisn, noms, types, taille_des_chaines); int k; for (k=0; k< naxisn.size(); k++) naxisn[k] *= nrows; if(hdutype == ASCII_TBL) { typeOfExtension = "ASCII_TBL"; dataType = "ASCII"; } else { typeOfExtension = "BINARY_TBL"; if(types[0] == 'D') dataType = "double"; else if(types[0] == 'E') dataType = "float"; else if(types[0] == 'I' ) dataType = "int"; else if(types[0] == 'S' ) dataType = "char*"; else { cout << " types[0]= " << types[0] << endl; throw PException(" FitsFile::getBlockType : unsupprted FITS data type"); } } } else { cout << " hdutype= " << hdutype << endl; throw IOExc("FitsFile::getBlockType: this HDU type is unknown"); } KeywordsIntoDVList(fileptr, dvl, hdunum); fits_close_file(fileptr,&status); } void FitsFile::ReadF(char flnm[],int hdunum) { int status = 0; hdutype_= 0; fits_open_file(&fptr_,flnm,READONLY,&status); if( status ) printerror( status ); // if (hdunum <= 1) { hdunum_ = 1; // presence of image ? int naxis= 0; fits_read_key(fptr_,TINT,"NAXIS",&naxis,NULL,&status); if( status ) printerror( status ); if (naxis > 0 ) // there is an image { hdutype_ = IMAGE_HDU; GetImageParameters (fptr_, bitpix_, naxis_, naxisn_); nbData_ = 1; int k; for (k=0; k 0) nbData_ *= naxisn_[k]; KeywordsIntoDVList(fptr_, dvl_,hdunum_); } else { throw PException(" first header : no image, probably error in hdunum"); } // } else { hdunum_ = hdunum-1; int hdutype; fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status); if( status ) printerror( status,":FitsFile::ReadF : erreur movabs"); moveToFollowingHeader(); } ReadFromFits(*this); } void FitsFile::moveToFollowingHeader() { int status = 0; int hdutype; fits_movrel_hdu(fptr_, 1,&hdutype,&status); if( status ) printerror( status," lecture du header suivant" ); hdunum_++; hdutype_= hdutype; if(hdutype_ == IMAGE_HDU) { GetImageParameters (fptr_, bitpix_, naxis_, naxisn_); nbData_ = 1; int k; for (k=0; k 0) nbData_ *= naxisn_[k]; KeywordsIntoDVList(fptr_, dvl_,hdunum_); } if(hdutype_ == ASCII_TBL || hdutype_ == BINARY_TBL) { GetBinTabParameters(fptr_,nbcols_, nrows_,repeat_, noms_, types_, taille_des_chaines_); KeywordsIntoDVList(fptr_, dvl_, hdunum_); } } void FitsFile::WriteF(char flnm[], bool OldFile) { int status = 0; hdutype_= 0; hdunum_ = 0; // create new FITS file if (!OldFile) { fits_create_file(&fptr_,flnm,&status); if( status ) printerror(status,"file already exists"); } else { fits_open_file(&fptr_,flnm,READWRITE,&status); if( status ) printerror(status,"file does not exist"); fits_get_num_hdus(fptr_, &hdunum_, &status); int hdutype; fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status); if( status ) printerror( status,":FitsFile::WriteF : erreur movabs"); } WriteToFits(*this); } void FitsFile::GetSingleColumn(double* map, int nentries) const { int status = 0; if(hdutype_ == IMAGE_HDU) { if(bitpix_ != DOUBLE_IMG) { cout << " The data type on fits file is not double..."; cout << " Conversion to double achieved by cfitsio lib" << endl; } // no checking for undefined pixels int anull; double dnull= 0.; long nels= nentries; fits_read_img(fptr_,TDOUBLE,1,nels,&dnull,map,&anull,&status); if( status ) printerror( status ); } else if(hdutype_ == ASCII_TBL || hdutype_ == BINARY_TBL) { GetBinTabFCol(map,nentries, 0); } else { cout << " hdutype= " << hdutype_ << endl; throw IOExc("FitsFile::GetSingleColumn this HDU is unknown"); } } void FitsFile::GetSingleColumn(float* map, int nentries) const { int status = 0; if(hdutype_ == IMAGE_HDU) { if(bitpix_ != FLOAT_IMG) { cout << " The data type on fits file is not float "; cout << " Conversion to float achieved by cfitsio lib" << endl; } // no checking for undefined pixels int anull; float fnull= 0.; long nels= nentries; fits_read_img(fptr_,TFLOAT,1,nels,&fnull, map,&anull,&status); if( status ) printerror( status ); } else if(hdutype_ == ASCII_TBL || hdutype_ == BINARY_TBL) { GetBinTabFCol(map,nentries, 0); } else { cout << " hdutype= " << hdutype_ << endl; throw IOExc("FitsFile::GetSingleColumn this HDU is unknown"); } } void FitsFile::GetSingleColumn( int* map, int nentries) const { int status = 0; if(hdutype_ == IMAGE_HDU) { if(bitpix_ != LONG_IMG) { cout << " The data type on fits file is not int "; cout << " Conversion to float achieved by cfitsio lib" << endl; } // no checking for undefined pixels int anull; float fnull= 0.; long nels= nentries; fits_read_img(fptr_,TINT,1,nels,&fnull,map,&anull,&status); if( status ) printerror( status ); } else if(hdutype_ == ASCII_TBL || hdutype_ == BINARY_TBL) { GetBinTabFCol(map,nentries, 0); } else { cout << " hdutype= " << hdutype_ << endl; throw IOExc("FitsFile::GetSingleColumn this HDU is unknown"); } } void FitsFile::makeHeaderImageOnFits(char type, int nbdim, int* naxisn) { int status = 0; long naxis = nbdim; long* naxes = new long[nbdim]; int k; for (k=0; k< nbdim; k++) naxes[k] = (long)naxisn[k]; if (type == 'D') fits_create_img(fptr_,DOUBLE_IMG,naxis,naxes,&status); else if (type == 'E') fits_create_img(fptr_,FLOAT_IMG,naxis,naxes,&status); else if (type == 'I') fits_create_img(fptr_,LONG_IMG,naxis,naxes,&status); else { cout << " type of data: " << type << endl; throw PException("FitsFile:::makeHeaderImageOnFits:unprogrammed type of data "); } hdunum_++; delete [] naxes; if( status ) printerror( status ); } void FitsFile::putImageToFits(int nbData, double* map) const { int status = 0; long npix= nbData; fits_write_img(fptr_,TDOUBLE,1,npix,map,&status); if( status ) printerror( status ); writeSignatureOnFits(); } void FitsFile::putImageToFits(int nbData, float* map) const { int status = 0; long npix= nbData; fits_write_img(fptr_,TFLOAT,1,npix, map,&status); if( status ) printerror( status ); writeSignatureOnFits(); } void FitsFile::putImageToFits( int nbData, int* map) const { int status = 0; long npix= nbData; fits_write_img(fptr_,TINT,1,npix,map,&status); if( status ) printerror( status ); writeSignatureOnFits(); } void FitsFile::GetImageParameters (fitsfile* fileptr,int& bitpix,int& naxis,vector& naxisn) { int hdunum=0; cout << " Reading a FITS image in HDU : " << fits_get_hdu_num(fileptr,&hdunum) << endl; int status= 0; // bits per pixels fits_read_key(fileptr,TINT,"BITPIX",&bitpix,NULL,&status); if( status ) printerror( status ); // number of dimensions in the FITS array naxis= 0; fits_read_key(fileptr,TINT,"NAXIS",&naxis,NULL,&status); if( status ) printerror( status ); // read the NAXISn keywords to get image size long* naxes = new long[naxis] ; int nfound; fits_read_keys_lng(fileptr,"NAXIS",1,naxis,naxes,&nfound,&status); if( status ) printerror( status ); if (nfound != naxis ) cout << " WARNING : " << nfound << " axes found, expected naxis= " << naxis << endl; int k; for (k=0; k& repeat, vector& noms, vector& types, vector& taille_des_chaines) { int status= 0; int hdunum=0; int hdutype=0; fits_get_hdu_num(fileptr,&hdunum); fits_get_hdu_type(fileptr, &hdutype, &status); if(hdutype != ASCII_TBL && hdutype != BINARY_TBL) { throw IOExc("FitsFile::GetBinTabParameters this HDU is not an ASCII table nor a binary table"); } if(hdutype == ASCII_TBL) cout << " Reading a FITS ascii table in HDU : " << hdunum << endl; if(hdutype == BINARY_TBL) cout << " Reading a FITS binary table in HDU : " << hdunum << endl; // get the number of columns fits_get_num_cols(fileptr, &nbcols,&status); if( status ) printerror( status ); // get the number of rows long naxis2= 0; fits_get_num_rows(fileptr,&naxis2,&status); if( status ) printerror( status ); nrows = (int)naxis2; // get the datatype, names and the repeat count noms.clear(); noms.reserve(nbcols); types.clear(); types.reserve(nbcols); repeat.clear(); repeat.reserve(nbcols); taille_des_chaines.clear(); char **ttype = new char*[nbcols]; int ii; for (ii=0; ii < nbcols; ii++) ttype[ii]=new char[FLEN_VALUE]; int nfound; fits_read_keys_str(fileptr, "TTYPE",1,nbcols,ttype,&nfound, &status); if( status ) printerror( status,"erreur lecture des noms de colonne"); int rept=0; for(ii = 0; ii < nbcols; ii++) { int DTYPE; long width; long repete = 0; fits_get_coltype(fileptr,ii+1,&DTYPE,&repete,&width,&status); if( status ) printerror( status,"erreur lecture type de colonne"); rept = repete; noms.push_back(string(ttype[ii])); switch (DTYPE) { case TDOUBLE : types.push_back('D'); break; case TFLOAT : types.push_back('E'); break; case TLONG : types.push_back('I'); break; case TINT : types.push_back('I'); break; case TSHORT : types.push_back('I'); break; case TSTRING : types.push_back('S'); taille_des_chaines.push_back(width); rept/=width; break; default : cout << " field " << ii+1 << " DTYPE= " << DTYPE << endl; throw IOExc("FitsFile:: unknown type of field"); } repeat.push_back(rept); } } void FitsFile::makeHeaderBntblOnFits( char* fieldType, char** Noms, int nentries, int tfields, DVList &dvl, char* extname, vector taille_des_chaines) { int status = 0; long nrows; if (strlen(fieldType) != tfields) { cout << " nombre de champs :" << tfields << "nombre de types: " << strlen(fieldType) << endl; throw ParmError("FitsFile:: fields and types don't match"); } char ** ttype= new char*[tfields]; char ** tform= new char*[tfields]; char largeur[FLEN_VALUE]; int noColString=0; int k; for (k=0; k0) { throw IOExc("FitsFile::putColToFits, this HDU is an ASCII table, nocol>0 forbidden"); } int code; long repeat, width; fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status); if( code != TFLOAT) { cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " (FITS code), to be written= FLOAT " << endl; } fits_write_col(fptr_,TFLOAT,nocol+1,1,1,nentries, donnees ,&status); if( status ) printerror( status,"erreur ecriture du fichier fits" ); } void FitsFile::putColToFits(int nocol, int nentries, int* donnees) const { int status = 0; int hdutype; fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status); if( status ) printerror(status,"putColToFits: le movabs a foire"); fits_get_hdu_type(fptr_, &hdutype, &status); if(hdutype != ASCII_TBL && hdutype != BINARY_TBL) { cout << " hdunum= " << hdunum_ << " hdutype= " << hdutype << endl; throw IOExc("FitsFile::putColToFits, this HDU is not an ASCII table nor a binary table"); } if(hdutype == ASCII_TBL && nocol>0) { throw IOExc("FitsFile::putColToFits, this HDU is an ASCII table, nocol>0 forbidden"); } int code; long repeat, width; fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status); if( code != TLONG && code != TINT && code != TSHORT ) { cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " (FITS code), to be written= INT " << endl; } fits_write_col(fptr_,TINT,nocol+1,1,1,nentries, donnees ,&status); if( status ) printerror( status," ecriture du fichier fits" ); } void FitsFile::putColToFits(int nocol, int nentries, char** donnees) const { int status = 0; int hdutype; fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status); if( status ) printerror(status,"putColToFits: le movabs a foire"); fits_get_hdu_type(fptr_, &hdutype, &status); if(hdutype != ASCII_TBL && hdutype != BINARY_TBL) { cout << " hdunum= " << hdunum_ << " hdutype= " << hdutype << endl; throw IOExc("FitsFile::putColToFits, this HDU is not an ASCII table nor a binary table"); } if(hdutype == ASCII_TBL && nocol>0) { throw IOExc("FitsFile::putColToFits, this HDU is an ASCII table, nocol>0 forbidden"); } int code; long repeat, width; fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status); if( code != TSTRING) { cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " (FITS code), to be written= char** " << endl; } fits_write_col(fptr_,TSTRING,nocol+1,1,1,nentries, donnees ,&status); if( status ) printerror( status,"erreur ecriture du fichier fits" ); } void FitsFile::readheader() //*************************************/ //* Print out all the header keywords */ //*************************************/ { // standard string lengths defined in fitsioc.h char card[FLEN_CARD]; int status = 0; // get the number of keywords int nkeys, keypos; if( fits_get_hdrpos(fptr_,&nkeys,&keypos,&status) ) printerror( status ); cout << " Header listing for HDU : " << hdunum_ << endl; int jj; for(jj = 1; jj <= nkeys; jj++) { if( fits_read_record(fptr_,jj,card,&status) ) printerror( status ); // print the keyword card cout << card << endl; } cout << "END" << endl; } void FitsFile::writeSignatureOnFits() const { int status = 0; char keyname[LEN_KEYWORD]; char strval[FLEN_VALUE]; char comment[FLEN_COMMENT]; strncpy(keyname, "CREATOR", LEN_KEYWORD); keyname[LEN_KEYWORD-1] = '\0'; strcpy(strval, "SOPHYA"); strcpy(comment," SOPHYA Package - FITSIOServer "); fits_write_key(fptr_, TSTRING, keyname, &strval, comment, &status); if( status ) printerror( status ); fits_write_comment(fptr_,"..............................................", &status); fits_write_comment(fptr_, " SOPHYA package - FITSIOSever ", &status); fits_write_comment(fptr_, " (C) LAL/IN2P3-CNRS Orsay, FRANCE 2000", &status); fits_write_comment(fptr_, " (C) DAPNIA/CEA Saclay, FRANCE 2000", &status); fits_write_comment(fptr_,"..............................................", &status); if( status ) printerror( status ); } void FitsFile::printerror(int &status) //*****************************************************/ //* Print out cfitsio error messages and exit program */ //*****************************************************/ { if( status ) { fits_report_error(stderr,status); throw IOExc("FitsFile:: error FITSIO status"); } return; } void FitsFile::printerror(int& status, char* texte) //*****************************************************/ //* Print out cfitsio error messages and exit program */ //*****************************************************/ { // print out cfitsio error messages and exit program // print error report fits_report_error(stderr, status); cout << " erreur:: " << texte << endl; throw IOExc("FitsFile:: error FITSIO status"); }