00001
00012 #include "FitsFile.h"
00013
00014 #include "pattern/string_convert.h"
00015
00016 #include <algorithm>
00017 #include <stdexcept>
00018
00019 #include <cassert>
00020
00021 using std::string;
00022 using std::vector;
00023
00024 using namespace hippodraw;
00025
00026 FitsFile::FitsFile ( const std::string & filename, bool write )
00027 : FitsFileBase ( filename, write )
00028 {
00029 }
00030
00031 void
00032 FitsFile::
00033 fillHDUNames ( std::vector < std::string > & names )
00034 {
00035 names.clear ();
00036 int number = getNumberOfHDU ( );
00037
00038 for ( int i = 0; i < number; i++ ) {
00039 moveToHDU ( i );
00040 string value = stringValueForKey ( "EXTNAME" );
00041
00042 if ( value == "" ) {
00043 if ( i == 0 ) {
00044 int number = intValueForKey ( "NAXIS" );
00045 if ( number == 0 ) {
00046 value = "Empty image";
00047 }
00048 }
00049 }
00050
00051 if ( value == "" ) value = "<no name>";
00052 names.push_back ( value );
00053 }
00054 }
00055
00056 int
00057 FitsFile::
00058 fillColumnNames ( std::vector < std::string > & labels )
00059 {
00060 HduType type = getHduType ();
00061 if ( type == Image ) {
00062 m_status = fillColumnNamesFromImage ( labels );
00063 }
00064 else {
00065 m_status = fillColumnNamesFromTable ( labels );
00066 }
00067
00068 return m_status;
00069 }
00070
00071 int
00072 FitsFile::
00073 fillColumnNamesFromTable ( std::vector < std::string > & labels )
00074 {
00075 labels.clear ();
00076 int columns = getNumberOfColumns ();
00077
00078 if ( columns != 0 ) {
00079 m_status = 0;
00080
00081 for ( int i = 0; i < columns; i++ ) {
00082 char colname [ FLEN_VALUE ];
00083 int col;
00084 fits_get_colname ( m_fptr, CASEINSEN, const_cast< char * > ("*"),
00085 colname, &col, &m_status );
00086 labels.push_back ( string ( colname ) );
00087 }
00088 }
00089 else {
00090 labels.push_back ( string ( "none" ) );
00091 }
00092
00093 return m_status;
00094 }
00095
00096 int
00097 FitsFile::
00098 fillColumnNamesFromImage ( std::vector < std::string > & labels )
00099 {
00100 labels.clear ();
00101 int dimension = getImageDimensions ();
00102 if ( dimension == 2 ) {
00103 labels.push_back ( "Values" );
00104 }
00105 else if ( dimension == 3 ) {
00106 const char type[] = "CTYPE3";
00107 bool yes = hasKey ( type );
00108 if ( yes ) {
00109 const string name = stringValueForKey ( type );
00110 int number = intValueForKey ( "NAXIS3" );
00111 for ( int i = 0; i < number; i++ ) {
00112 const string iv = String::convert ( i );
00113 labels.push_back ( name + " " + iv );
00114 }
00115 }
00116 }
00117
00118 return m_status;
00119 }
00120
00121 int
00122 FitsFile::
00123 fillDoubleVectorFromColumn ( std::vector < double > & vec,
00124 int column )
00125 {
00126 int retval = 0;
00127
00128 HduType type = getHduType ();
00129 if ( type == Btable || type == Atable ) {
00130 retval = fillFromTableColumn ( vec, column );
00131 }
00132 else {
00133 retval = fillFromImage ( vec, column );
00134 }
00135
00136 return retval;
00137 }
00138
00139 int
00140 FitsFile::
00141 fillFromTableColumn ( std::vector < double > & vec, int column )
00142 {
00143 int anynul;
00144 double nulval = 0;
00145 long nelements = vec.size();
00146 vector<double>::iterator it = vec.begin();
00147 double * ptr = &*it;
00148
00149 m_status = 0;
00150 fits_read_col_dbl( m_fptr, column+1, 1, 1, nelements, nulval,
00151 ptr, &anynul, &m_status);
00152
00153 return m_status;
00154 }
00155
00156 void
00157 FitsFile::
00158 fillShape ( std::vector < intptr_t > & shape, int column )
00159 {
00160 m_status = 0;
00161
00162 string tdimn ( "TDIM" );
00163 tdimn += String::convert ( column+1 );
00164 char value[80];
00165 fits_read_keyword ( m_fptr,
00166 const_cast < char * > ( tdimn.c_str() ),
00167 value,
00168 NULL, & m_status );
00169
00170 long rows = getNumberOfRows ();
00171 if ( value[0] != 0 ) {
00172 const string dims = value;
00173 unsigned int count = std::count ( dims.begin(), dims.end(), ',' );
00174 count += 1;
00175 int icount = 0;
00176 vector < long > naxes ( count );
00177 fits_read_tdim ( m_fptr, column+1, count,
00178 &icount, &naxes[0], & m_status );
00179 shape.resize ( count + 1 );
00180 shape [0] = rows;
00181 std::copy ( naxes.begin(), naxes.end(), shape.begin() + 1 );
00182 }
00183 else {
00184 shape.resize ( 1, rows );
00185 }
00186 }
00187
00188 int
00189 FitsFile::
00190 fillAxisSizes ( std::vector < long > & vec ) const
00191 {
00192 vec.clear ();
00193 int naxis = getImageDimensions ();
00194 vec.resize ( naxis );
00195
00196 vector < long > ::iterator first = vec.begin();
00197 long * ptr = & *first;
00198 m_status = 0;
00199 fits_get_img_size ( m_fptr, naxis, ptr, & m_status );
00200 assert ( m_status == 0 );
00201
00202 return m_status;
00203 }
00204
00205 void
00206 FitsFile::
00207 fillImageDeltas ( std::vector < double > & deltas ) const
00208 {
00209 deltas.clear ();
00210 int naxis = getImageDimensions ();
00211 deltas.reserve ( naxis );
00212
00213 char key [ FLEN_KEYWORD];
00214 char * keyroot = const_cast < char * > ("CDELT");
00215 for ( int i = 0; i < naxis; i++ ) {
00216 m_status = 0;
00217 fits_make_keyn ( keyroot, i+1, key, & m_status );
00218 bool yes = hasKey ( key );
00219
00220 if ( yes ) {
00221 double value = doubleValueForKey ( key );
00222 deltas.push_back ( value );
00223 }
00224 else {
00225 deltas.push_back ( 1.0 );
00226 }
00227 }
00228 }
00229
00230 void
00231 FitsFile::
00232 fillRefPixelIndices ( std::vector < int > & indices ) const
00233 {
00234 indices.clear ();
00235 int naxis = getImageDimensions ();
00236 indices.reserve ( naxis );
00237
00238 char key [ FLEN_KEYWORD];
00239 char * keyroot = const_cast < char * > ( "CRPIX" );
00240 for ( int i = 0; i < naxis; i++ ) {
00241 m_status = 0;
00242 fits_make_keyn ( keyroot, i+1, key, & m_status );
00243 bool yes = hasKey ( key );
00244
00245 if ( yes ) {
00246 int value = intValueForKey ( key );
00247 indices.push_back ( value );
00248 }
00249 else {
00250 indices.push_back ( 1 );
00251 }
00252 }
00253 }
00254
00255 void
00256 FitsFile::
00257 fillRefPixelValues ( std::vector < double > & values ) const
00258 {
00259 values.clear ();
00260 int naxis = getImageDimensions ();
00261 values.reserve ( naxis );
00262
00263 char key [ FLEN_KEYWORD];
00264 char * keyroot = const_cast < char * > ( "CRVAL" );
00265 for ( int i = 0; i < naxis; i++ ) {
00266 m_status = 0;
00267 fits_make_keyn ( keyroot, i+1, key, & m_status );
00268 bool yes = hasKey ( key );
00269
00270 if ( yes ) {
00271 double value = doubleValueForKey ( key );
00272 values.push_back ( value );
00273 }
00274 else {
00275 values.push_back ( 0. );
00276 }
00277 }
00278 }
00279
00280 bool
00281 FitsFile::
00282 isHammerAitoff () const
00283 {
00284 bool hammer = true;
00285
00286 char key [FLEN_KEYWORD];
00287 char * keyroot = const_cast < char * > ( "CTYPE" );
00288 for ( int i = 0; i < 2; i++ ) {
00289 fits_make_keyn ( keyroot, i+1, key, & m_status );
00290 bool yes = hasKey ( key );
00291
00292 if ( yes ) {
00293 string value = stringValueForKey ( key );
00294 if ( value.find ( "-AIT" ) == string::npos ) {
00295 hammer = false;
00296 }
00297 } else {
00298 hammer = false;
00299 }
00300 }
00301
00302 return hammer;
00303 }
00304
00305 int
00306 FitsFile::
00307 fillFromImage ( std::vector < double > & vec, unsigned int zplane )
00308 {
00309 vector < long > naxes;
00310 fillAxisSizes ( naxes );
00311
00312 int datatype = Double;
00313 long nelements = naxes[0] * naxes[1];
00314 double nulval = 0;
00315 vector < double >::iterator first = vec.begin();
00316 double * ptr = & *first;
00317 int anynul;
00318 m_status = 0;
00319 if ( naxes.size () == 2 ) {
00320 long fpixel[2] = { 1, 1 };
00321 fits_read_pix ( m_fptr, datatype, fpixel, nelements,
00322 & nulval, ptr, & anynul, & m_status );
00323 } else {
00324 long fpixel[] = { 1, 1, 1 };
00325 fpixel[2] = zplane + 1;
00326 fits_read_pix ( m_fptr, datatype, fpixel, nelements,
00327 & nulval, ptr, & anynul, & m_status );
00328 }
00329
00330 return m_status;
00331 }
00332
00333 int
00334 FitsFile::
00335 fillIntVectorFromColumn( std::vector < int > & vec, int column )
00336 {
00337 int anynul, status = 0;
00338 int nulval = 0;
00339 long nelements = vec.size();
00340 vector<int>::iterator it = vec.begin();
00341 int * ptr = new int [ nelements ];
00342
00343 m_status = 0;
00344 fits_read_col_int( m_fptr, column+1, 1, 1, nelements, nulval,
00345 ptr, &anynul, &status);
00346 copy ( ptr, ptr+nelements, it );
00347
00348 return status;
00349 }
00350
00351 void
00352 FitsFile::
00353 writeHDU ( long rows, int columns,
00354 const std::vector < std::string > & names,
00355 const std::vector < std::vector < int > > & shapes,
00356 const std::string & extname )
00357 {
00358 char ** types = new char * [ columns ];
00359 char ** tform = new char * [ columns ];
00360 char ** tunits = 0;
00361 char * m_extname;
00362
00363 vector < string > forms;
00364 for ( int i = 0; i < columns; i ++ ) {
00365 types[i]=const_cast<char*> (names[i].c_str());
00366 int size = 1;
00367 const vector < int > & shape = shapes [ i ];
00368 unsigned int rank = shape.size();
00369 if ( rank > 1 ) {
00370 for ( unsigned int j = 1; j < rank; j++ ) {
00371 int dim = shape [ j ];
00372 size *= dim;
00373 }
00374 }
00375 string form = String::convert ( size );
00376 form += "D";
00377 forms.push_back ( form );
00378
00379 tform[i] = strdup ( form.c_str() );
00380
00381
00382 }
00383
00384 m_extname = const_cast<char*> (extname.c_str());
00385
00386 fits_create_tbl( m_fptr, BINARY_TBL, 0, columns,
00387 types, tform, tunits, m_extname, &m_status);
00388
00389 for ( int i = 0; i < columns; i ++ ) {
00390 const vector < int > & shape = shapes [ i ];
00391 unsigned int rank = shape.size();
00392 if ( rank > 1 ) {
00393 vector < long > naxes;
00394 for ( unsigned int j = 1; j < rank; j++ ) {
00395 naxes.push_back ( shape [ j ] );
00396 }
00397 fits_write_tdim ( m_fptr, i+1,
00398 rank -1, &naxes [0], &m_status );
00399 }
00400 }
00401
00402 delete types;
00403 delete tform;
00404 }
00405
00406 void
00407 FitsFile::
00408 writeImageHDU ( long x, long y )
00409 {
00410 long naxes[2]={x,y};
00411 fits_create_img ( m_fptr, DOUBLE_IMG, 2, naxes, &m_status );
00412 }
00413
00414
00415 void
00416 FitsFile::
00417 writeColumn ( int col, const std::vector < double > & data )
00418 {
00419 LONGLONG nelements = data.size ();
00420 LONGLONG firstrow = 1;
00421 LONGLONG firstelem = 1;
00422 fits_write_col( m_fptr, TDOUBLE, col+1, firstrow, firstelem, nelements,
00423 const_cast < double * > ( &data[0] ), &m_status);
00424 }
00425
00426 void
00427 FitsFile::
00428 writePix ( long x, long y, const std::vector < double > & data )
00429 {
00430 long fpixel[2]={1,1};
00431 long nelements = x * y;
00432 fits_write_pix ( m_fptr, TDOUBLE, fpixel, nelements,
00433 const_cast < double * > ( &data[0] ), &m_status );
00434 }
00435
00436 void
00437 FitsFile::
00438 writeCloseFile ()
00439 {
00440 fits_close_file( m_fptr, &m_status );
00441 m_fptr = 0;
00442 }
00443
00444 bool
00445 FitsFile::
00446 pixCenter () const
00447 {
00448
00449 char * key = const_cast < char * > ( "PIXCENT" );
00450 bool yes = hasKey ( key );
00451 if ( yes )
00452 {
00453 if ( stringValueForKey ( key ) == "T" ) return true;
00454 }
00455
00456 return false;
00457
00458 }
00459
00460 void
00461 FitsFile::
00462 writeRefPixelValues ( double value1, double value2 )
00463 {
00464 char * key1 = const_cast < char * > ( "CRVAL1" );
00465 char * key2 = const_cast < char * > ( "CRVAL2" );
00466
00467 fits_write_key( m_fptr, TDOUBLE, key1, &value1, NULL, &m_status );
00468 fits_write_key( m_fptr, TDOUBLE, key2, &value2, NULL, &m_status );
00469 }