Main Page | Modules | Alphabetical List | Class List | File List | Class Members | File Members

fitshandle.cc

00001 /*
00002  *  This file is part of Healpix_cxx.
00003  *
00004  *  Healpix_cxx is free software; you can redistribute it and/or modify
00005  *  it under the terms of the GNU General Public License as published by
00006  *  the Free Software Foundation; either version 2 of the License, or
00007  *  (at your option) any later version.
00008  *
00009  *  Healpix_cxx is distributed in the hope that it will be useful,
00010  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
00011  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00012  *  GNU General Public License for more details.
00013  *
00014  *  You should have received a copy of the GNU General Public License
00015  *  along with Healpix_cxx; if not, write to the Free Software
00016  *  Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
00017  *
00018  *  For more information about HEALPix, see http://healpix.jpl.nasa.gov
00019  */
00020 
00021 /*
00022  *  Healpix_cxx is being developed at the Max-Planck-Institut fuer Astrophysik
00023  *  and financially supported by the Deutsches Zentrum fuer Luft- und Raumfahrt
00024  *  (DLR).
00025  */
00026 
00027 /*
00028  *  This file contains the implementation of the FITS I/O helper class
00029  *  used by the Planck LevelS package.
00030  *
00031  *  Copyright (C) 2002, 2003, 2004 Max-Planck-Society
00032  *  Author: Martin Reinecke
00033  */
00034 
00035 #include <iostream>
00036 #include <string>
00037 #include <sstream>
00038 #include <vector>
00039 #include <map>
00040 #include "fitsio.h"
00041 #include "fitshandle.h"
00042 #include "cxxutils.h"
00043 
00044 using namespace std;
00045 
00046 string type2char (int type)
00047   {
00048   if (type==TFLOAT) return "E";
00049   if (type==TDOUBLE) return "D";
00050   if (type==TINT) return "J";
00051   if (type==TINT32BIT) return "J";
00052   if (type==TSTRING) return "A";
00053   throw Message_error("wrong datatype in type2char()");
00054   }
00055 
00056 string type2asciiform (int type)
00057   {
00058   if (type==TFLOAT) return "E14.7";
00059   if (type==TDOUBLE) return "E23.15";
00060   if (type==TINT) return "I11";
00061   if (type==TINT32BIT) return "I11";
00062   throw Message_error("wrong datatype in type2asciiform()");
00063   }
00064 
00065 void fitshandle::check_errors()
00066   {
00067   if (status==0) return;
00068   char msg[81];
00069   fits_get_errstatus (status, msg);
00070   cerr << msg << endl;
00071   while (fits_read_errmsg(msg)) cerr << msg << endl;
00072   throw Message_error("FITS error");
00073   }
00074 
00075 void fitshandle::clean_data()
00076   {
00077   if (!fptr) return;
00078   axes_.clear();
00079   columns_.clear();
00080   hdutype_=INVALID;
00081   bitpix_=INVALID;
00082   nrows_=0;
00083   }
00084 
00085 void fitshandle::clean_all()
00086   {
00087   if (!fptr) return;
00088   clean_data();
00089   fits_close_file (fptr, &status);
00090   check_errors();
00091   fptr=0;
00092   }
00093 
00094 void fitshandle::init_image()
00095   {
00096   int naxis;
00097   fits_get_img_type (fptr, &bitpix_, &status);
00098   fits_get_img_dim (fptr, &naxis, &status);
00099   check_errors();
00100   arr<long> naxes(naxis);
00101   fits_get_img_size (fptr, naxis, naxes, &status);
00102   for (long m=0; m<naxis; ++m) axes_.push_back(naxes[naxis-m-1]);
00103   check_errors();
00104   }
00105 
00106 void fitshandle::init_asciitab()
00107   {
00108   char ttype[81], tunit[81], tform[81];
00109   int ncol, typecode;
00110   fits_get_num_cols (fptr, &ncol, &status);
00111   fits_get_num_rows (fptr, &nrows_, &status);
00112   check_errors();
00113   for (int m=1; m<=ncol; ++m)
00114     {
00115     fits_get_acolparms (fptr, m, ttype, 0, tunit, tform,
00116                         0, 0, 0, 0, &status);
00117     fits_ascii_tform (tform, &typecode, 0,0, &status);
00118     check_errors();
00119     columns_.push_back (fitscolumn (ttype,tunit,1,typecode));
00120     }
00121   }
00122 
00123 void fitshandle::init_bintab()
00124   {
00125   char ttype[81], tunit[81], tform[81];
00126   long repc;
00127   int ncol, typecode;
00128   fits_get_num_cols (fptr, &ncol, &status);
00129   fits_get_num_rows (fptr, &nrows_, &status);
00130   check_errors();
00131   for (int m=1; m<=ncol; ++m)
00132     {
00133     fits_get_bcolparms (fptr, m, ttype, tunit, tform, &repc,
00134                         0, 0, 0, 0, &status);
00135     fits_binary_tform (tform, &typecode, 0,0, &status);
00136     check_errors();
00137     columns_.push_back (fitscolumn (ttype,tunit,repc,typecode));
00138     }
00139   }
00140 
00141 void fitshandle::init_data()
00142   {
00143   switch(hdutype_)
00144     {
00145     case IMAGE_HDU:
00146       init_image();
00147       break;
00148     case ASCII_TBL:
00149       init_asciitab();
00150       break;
00151     case BINARY_TBL:
00152       init_bintab();
00153       break;
00154     default:
00155       throw Message_error("init_data(): wrong HDU type");
00156       break;
00157     }
00158   }
00159 
00160 void fitshandle::read_col (int colnum, void *data, long ndata, int dtype,
00161   long offset)
00162   {
00163   assert_table_hdu("fitshandle::read_column()",colnum);
00164   long repc = columns_[colnum-1].repcount();
00165   planck_assert (ndata<=(repc*nrows_-offset),
00166                  "read_column(): array too large");
00167   long frow = offset/repc+1;
00168   long felem = offset%repc+1;
00169   fits_read_col (fptr, dtype, colnum, frow, felem, ndata, 0, data, 0, &status);
00170   check_errors();
00171   }
00172 void fitshandle::write_col (int colnum, const void *data, long ndata,
00173   int dtype, long offset)
00174   {
00175   assert_table_hdu("fitshandle::write_column()",colnum);
00176   long repc = columns_[colnum-1].repcount();
00177   long frow = offset/repc+1;
00178   long felem = offset%repc+1;
00179   fits_write_col (fptr, dtype, colnum, frow, felem, ndata,
00180     const_cast<void *>(data), &status);
00181   fits_get_num_rows (fptr, &nrows_, &status);
00182   check_errors();
00183   }
00184 
00185 void fitshandle::open (const string &fname, int rwmode)
00186   {
00187   clean_all();
00188   fits_open_file(&fptr, fname.c_str(), rwmode, &status);
00189   check_errors();
00190   }
00191 
00192 void fitshandle::create (const string &fname)
00193   {
00194   clean_all();
00195   fits_create_file(&fptr, fname.c_str(), &status);
00196   fits_write_imghdr(fptr, 8, 0, 0, &status); // insert empty PHDU
00197   fits_write_date(fptr, &status);
00198   check_errors();
00199   }
00200 
00201 void fitshandle::goto_hdu (int hdu)
00202   {
00203   clean_data();
00204   fits_movabs_hdu(fptr, hdu, &hdutype_, &status);
00205   check_errors();
00206   init_data();
00207   }
00208 
00209 void fitshandle::insert_bintab (const vector<fitscolumn> &cols)
00210   {
00211   clean_data();
00212   int ncol=cols.size();
00213   arr2b<char> ttype(ncol,81), tform(ncol,81), tunit(ncol,81);
00214 
00215   for (long m=0; m<ncol; ++m)
00216     {
00217     strcpy (ttype[m], cols[m].name().c_str());
00218     strcpy (tunit[m], cols[m].unit().c_str());
00219     ostringstream x;
00220     x << cols[m].repcount() << type2char (cols[m].type());
00221     strcpy (tform[m], x.str().c_str());
00222     }
00223   fits_insert_btbl (fptr, nrows_, ncol, ttype, tform, tunit,
00224     const_cast<char *>("xtension"), 0, &status);
00225   check_errors();
00226   hdutype_=BINARY_TBL;
00227   columns_=cols;
00228   }
00229 
00230 void fitshandle::insert_asctab (const vector<fitscolumn> &cols)
00231   {
00232   clean_data();
00233   int ncol=cols.size();
00234   arr2b<char> ttype(ncol,81), tform(ncol,81), tunit(ncol,81);
00235 
00236   for (long m=0; m<ncol; ++m)
00237     {
00238     strcpy (ttype[m], cols[m].name().c_str());
00239     strcpy (tunit[m], cols[m].unit().c_str());
00240     ostringstream x;
00241     if (cols[m].type()!=TSTRING)
00242       {
00243       planck_assert (cols[m].repcount()==1,"bad repcount for ASCII table");
00244       x << type2asciiform (cols[m].type());
00245       }
00246     else
00247       {
00248       x << "A" << dataToString(cols[m].repcount());
00249       }
00250     strcpy (tform[m], x.str().c_str());
00251     }
00252   fits_insert_atbl (fptr, 0, nrows_, ncol, ttype, 0, tform, tunit,
00253     const_cast<char *>("xtension"), &status);
00254   check_errors();
00255   hdutype_=ASCII_TBL;
00256   columns_=cols;
00257   }
00258 
00259 void fitshandle::insert_image (int btpx, const vector<long> &Axes)
00260   {
00261   clean_data();
00262   hdutype_=IMAGE_HDU;
00263   bitpix_=btpx;
00264   axes_=Axes;
00265   arr<long> tmpax(Axes.size());
00266   for (long m=0; m<long(Axes.size()); m++) tmpax[m]=Axes[Axes.size()-1-m];
00267   fits_insert_img (fptr, btpx, Axes.size(), tmpax, &status);
00268   check_errors();
00269   }
00270 
00271 template<typename T>
00272   void fitshandle::insert_image (int btpx, const arr2<T> &data)
00273   {
00274   clean_data();
00275   hdutype_=IMAGE_HDU;
00276   bitpix_=btpx;
00277   axes_.push_back(data.size1());
00278   axes_.push_back(data.size2());
00279   arr<long> tmpax(axes_.size());
00280   for (long m=0; m<long(axes_.size()); m++) tmpax[m]=axes_[axes_.size()-1-m];
00281   fits_insert_img (fptr, btpx, axes_.size(), tmpax, &status);
00282   arr2<T> &tmparr = const_cast<arr2<T> &> (data);
00283   fits_write_img (fptr, FITSUTIL<T>::DTYPE, 1, axes_[0]*axes_[1],
00284     &tmparr[0][0], &status);
00285   check_errors();
00286   }
00287 
00288 template void fitshandle::insert_image (int btpx, const arr2<double> &data);
00289 template void fitshandle::insert_image (int btpx, const arr2<float> &data);
00290 template void fitshandle::insert_image (int btpx, const arr2<int> &data);
00291 
00292 void fitshandle::copy_historified_header (const fitshandle &orig)
00293   {
00294   const char *inclist[] = { "*" };
00295   const char *exclist[] = {
00296       "SIMPLE","BITPIX","NAXIS","NAXIS#","PCOUNT","GCOUNT",
00297       "EXTEND","ORIGIN","DATE*","TFIELDS","TTYPE#","TFORM#",
00298       "TUNIT#","EXTNAME","CTYPE#","CRVAL#","CRPIX#","CDELT#",
00299       "XTENSION","INSTRUME","TELESCOP","PDMTYPE","TBCOL#" };
00300   char card[81];
00301   string card2;
00302   orig.assert_connected("fitshandle::copy_historified_header()");
00303   assert_connected("fitshandle::copy_historified_header()");
00304   fits_read_record (orig.fptr, 0, card, &status);
00305   check_errors();
00306   while (true)
00307     {
00308     fits_find_nextkey (orig.fptr, const_cast<char **>(inclist), 1,
00309       const_cast<char **>(exclist), 23, card, &status);
00310     if (status!=0) break;
00311     card2=trim(card);
00312     if (card2!="END" && card2!="COMMENT" && card2!="HISTORY")
00313       {
00314       if (card2.find("COMMENT")==0)
00315         card2.replace(0,7,"HISTORY");
00316       if (card2.find("HISTORY")!=0)
00317         card2.insert(0,"HISTORY ");
00318       if (card2.length()<=80)
00319         fits_write_record (fptr, card2.c_str(), &status);
00320       else
00321         {
00322         fits_write_record (fptr, card2.substr(0,80).c_str(), &status);
00323         card2=card2.substr(80,string::npos);
00324         card2.insert(0,"HISTORY ");
00325         fits_write_record (fptr, card2.c_str(),&status);
00326         }
00327       }
00328     check_errors();
00329     }
00330   if (status==KEY_NO_EXIST) { fits_clear_errmsg(); status=0; }
00331   check_errors();
00332   }
00333 
00334 void fitshandle::copy_header (const fitshandle &orig)
00335   {
00336   const char *inclist[] = { "*" };
00337   const char *exclist[] = {
00338       "SIMPLE","BITPIX","NAXIS","NAXIS#","PCOUNT","GCOUNT",
00339       "EXTEND","ORIGIN","DATE*","TFIELDS","TTYPE#","TFORM#",
00340       "TUNIT#","EXTNAME","CTYPE#","CRVAL#","CRPIX#","CDELT#",
00341       "XTENSION","INSTRUME","TELESCOP","PDMTYPE","TBCOL#" };
00342   char card[81];
00343   string card2;
00344   orig.assert_connected("fitshandle::copy_header()");
00345   assert_connected("fitshandle::copy_header()");
00346   fits_read_record (orig.fptr, 0, card, &status);
00347   check_errors();
00348   while (true)
00349     {
00350     fits_find_nextkey (orig.fptr, const_cast<char **>(inclist), 1,
00351       const_cast<char **>(exclist), 23, card, &status);
00352     if (status!=0) break;
00353     card2=trim(card);
00354     if (card2!="END" && card2!="COMMENT" && card2!="HISTORY")
00355       fits_write_record (fptr, card, &status);
00356     check_errors();
00357     }
00358   if (status==KEY_NO_EXIST) { fits_clear_errmsg(); status=0; }
00359   check_errors();
00360   }
00361 
00362 void fitshandle::check_key_present(const string &name)
00363   {
00364   char card[81];
00365   fits_read_card(fptr, const_cast<char *>(name.c_str()), card, &status);
00366   if (status==KEY_NO_EXIST)
00367     { fits_clear_errmsg(); status=0; return; }
00368   check_errors();
00369 // FIXME: disabled for now; but this issue needs to be resolved!
00370 //  cerr << "Warning: key " << name << " set more than once!" << endl;
00371   }
00372 
00373 template<typename T> void fitshandle::add_key (const string &name,
00374   const T &value, const string &comment)
00375   {
00376   assert_connected("fitshandle::add_key()");
00377   check_key_present (name);
00378   fits_write_key (fptr, FITSUTIL<T>::DTYPE, const_cast<char *>(name.c_str()),
00379     const_cast<T *>(&value), const_cast<char *>(comment.c_str()), &status);
00380   check_errors();
00381   }
00382 
00383 template void fitshandle::add_key(const string &name,
00384   const int &value, const string &comment);
00385 template void fitshandle::add_key(const string &name,
00386   const long &value, const string &comment);
00387 template void fitshandle::add_key(const string &name,
00388   const float &value, const string &comment);
00389 template void fitshandle::add_key(const string &name,
00390   const double &value, const string &comment);
00391 template<> void fitshandle::add_key(const string &name,
00392   const bool &value, const string &comment)
00393   {
00394   assert_connected("fitshandle::add_key()");
00395   check_key_present (name);
00396   int val=value;
00397   fits_write_key (fptr, TLOGICAL, const_cast<char *>(name.c_str()),
00398     &val, const_cast<char *>(comment.c_str()),
00399     &status);
00400   check_errors();
00401   }
00402 template<> void fitshandle::add_key (const string &name,
00403   const string &value, const string &comment)
00404   {
00405   assert_connected("fitshandle::add_key()");
00406   check_key_present (name);
00407   fits_write_key (fptr, TSTRING, const_cast<char *>(name.c_str()),
00408     const_cast<char *>(value.c_str()), const_cast<char *>(comment.c_str()),
00409     &status);
00410   check_errors();
00411   }
00412 
00413 template<typename T> void fitshandle::update_key (const string &name,
00414   const T &value, const string &comment)
00415   {
00416   assert_connected("fitshandle::update_key()");
00417   fits_update_key (fptr, FITSUTIL<T>::DTYPE, const_cast<char *>(name.c_str()),
00418     const_cast<T *>(&value), const_cast<char *>(comment.c_str()), &status);
00419   check_errors();
00420   }
00421 
00422 template void fitshandle::update_key(const string &name,
00423   const int &value, const string &comment);
00424 template void fitshandle::update_key(const string &name,
00425   const long &value, const string &comment);
00426 template void fitshandle::update_key(const string &name,
00427   const float &value, const string &comment);
00428 template void fitshandle::update_key(const string &name,
00429   const double &value, const string &comment);
00430 template<> void fitshandle::update_key(const string &name,
00431   const bool &value, const string &comment)
00432   {
00433   assert_connected("fitshandle::update_key()");
00434   int val=value;
00435   fits_update_key (fptr, TLOGICAL, const_cast<char *>(name.c_str()),
00436     &val, const_cast<char *>(comment.c_str()),
00437     &status);
00438   check_errors();
00439   }
00440 template<> void fitshandle::update_key (const string &name,
00441   const string &value, const string &comment)
00442   {
00443   assert_connected("fitshandle::update_key()");
00444   fits_update_key (fptr, TSTRING, const_cast<char *>(name.c_str()),
00445     const_cast<char *>(value.c_str()), const_cast<char *>(comment.c_str()),
00446     &status);
00447   check_errors();
00448   }
00449 
00450 void fitshandle::delete_key (const string &name)
00451   {
00452   assert_connected("fitshandle::delete_key()");
00453   fits_delete_key (fptr, const_cast<char *>(name.c_str()), &status);
00454   check_errors();
00455   }
00456 
00457 void fitshandle::add_comment(const string &comment)
00458   {
00459   assert_connected("fitshandle::add_comment()");
00460   fits_write_comment(fptr,const_cast<char *>(comment.c_str()),&status);
00461   check_errors();
00462   }
00463 
00464 template<typename T> void fitshandle::get_key (const string &name, T &value)
00465   {
00466   assert_connected("fitshandle::get_key()");
00467   fits_read_key (fptr, FITSUTIL<T>::DTYPE, const_cast<char *>(name.c_str()),
00468     &value, 0, &status);
00469   if (status==KEY_NO_EXIST) throw Message_error
00470     ("Fitshandle::get_key(): key "+name+" not found");
00471   check_errors();
00472   }
00473 template void fitshandle::get_key(const string &name,int &value);
00474 template void fitshandle::get_key(const string &name,long &value);
00475 template void fitshandle::get_key(const string &name,float &value);
00476 template void fitshandle::get_key(const string &name,double &value);
00477 template<> void fitshandle::get_key(const string &name,bool &value)
00478   {
00479   assert_connected("fitshandle::get_key()");
00480   int val;
00481   fits_read_key (fptr, TLOGICAL, const_cast<char *>(name.c_str()), &val, 0,
00482     &status);
00483   if (status==KEY_NO_EXIST) throw Message_error
00484     ("Fitshandle::get_key(): key "+name+" not found");
00485   check_errors();
00486   value=val;
00487   }
00488 template<> void fitshandle::get_key (const string &name,string &value)
00489   {
00490   char tmp[2049];
00491   assert_connected("fitshandle::get_key()");
00492   fits_read_key (fptr, TSTRING, const_cast<char *>(name.c_str()), tmp, 0,
00493     &status);
00494   if (status==KEY_NO_EXIST) throw Message_error
00495     ("Fitshandle::get_key(): key "+name+" not found");
00496   check_errors();
00497   value=tmp;
00498   }
00499 
00500 bool fitshandle::key_present(const string &name)
00501   {
00502   char card[81];
00503   fits_read_card(fptr, const_cast<char *>(name.c_str()), card, &status);
00504   if (status==KEY_NO_EXIST)
00505     { fits_clear_errmsg(); status=0; return false; }
00506   check_errors();
00507   return true;
00508   }
00509 
00510 void fitshandle::assert_pdmtype (const string &pdmtype)
00511   {
00512   string type;
00513   get_key("PDMTYPE",type);
00514   if (pdmtype==type) return;
00515   cerr << "PDMTYPE " << pdmtype << " expected, but found " << type << endl;
00516   }
00517 
00518 void fitshandle::read_column_raw_void
00519   (int colnum, void *data, int type, long num, long offset)
00520   {
00521   switch (type)
00522     {
00523     case PLANCK_INT32:
00524       read_col (colnum, data, num, TINT, offset); break;
00525     case PLANCK_FLOAT32:
00526       read_col (colnum, data, num, TFLOAT, offset); break;
00527     case PLANCK_FLOAT64:
00528       read_col (colnum, data, num, TDOUBLE, offset); break;
00529     case PLANCK_STRING:
00530       {
00531       string *data2 = static_cast<string *> (data);
00532       assert_table_hdu("fitshandle::read_column()",colnum);
00533       planck_assert (num<=(nrows_-offset),
00534         "read_column(): array too large");
00535       arr2b<char> tdata(num, columns_[colnum-1].repcount()+1);
00536       fits_read_col (fptr, TSTRING, colnum, offset+1, 1, num,
00537         0, tdata, 0, &status);
00538       check_errors();
00539       for (long m=0;m<num;++m) data2[m]=tdata[m];
00540       break;
00541       }
00542     default:
00543       throw Message_error ("unsupported data type in read_column_raw_void()");
00544     }
00545   }
00546 
00547 void fitshandle::write_column_raw_void
00548   (int colnum, const void *data, int type, long num, long offset)
00549   {
00550   switch (type)
00551     {
00552     case PLANCK_INT32:
00553       write_col (colnum, data, num, TINT, offset); break;
00554     case PLANCK_FLOAT32:
00555       write_col (colnum, data, num, TFLOAT, offset); break;
00556     case PLANCK_FLOAT64:
00557       write_col (colnum, data, num, TDOUBLE, offset); break;
00558     case PLANCK_STRING:
00559       {
00560       const string *data2 = static_cast<const string *> (data);
00561       assert_table_hdu("fitshandle::write_column()",colnum);
00562       int stringlen = columns_[colnum-1].repcount()+1;
00563       arr2b<char> tdata(num, stringlen);
00564       for (long m=0;m<num;++m)
00565         {
00566         strncpy(tdata[m],data2[m].c_str(),stringlen-1);
00567         tdata[m][stringlen-1] = '\0';
00568         }
00569       fits_write_col (fptr, TSTRING, colnum, offset+1, 1, num,
00570         tdata, &status);
00571       fits_get_num_rows (fptr, &nrows_, &status);
00572       check_errors();
00573       break;
00574       }
00575     default:
00576       throw Message_error ("unsupported data type in write_column_raw_void()");
00577     }
00578   }
00579 
00580 template<typename T> void fitshandle::write_image (const arr2<T> &data)
00581   {
00582   assert_image_hdu("fitshandle::write_image()");
00583   planck_assert (axes_.size()==2, "wrong number of dimensions");
00584   planck_assert (axes_[0]==data.size1(), "wrong size of dimension 1");
00585   planck_assert (axes_[1]==data.size2(), "wrong size of dimension 2");
00586 
00587   fits_write_img (fptr, FITSUTIL<T>::DTYPE, 1, axes_[0]*axes_[1],
00588     const_cast<T *>(&data[0][0]), &status);
00589   check_errors();
00590   }
00591 
00592 template void fitshandle::write_image (const arr2<float> &data);
00593 template void fitshandle::write_image (const arr2<double> &data);
00594 template void fitshandle::write_image (const arr2<int> &data);
00595 
00596 template<typename T> void fitshandle::write_subimage
00597   (const arr<T> &data, long offset)
00598   {
00599   assert_image_hdu("fitshandle::write_subimage()");
00600   fits_write_img (fptr, FITSUTIL<T>::DTYPE, 1+offset,
00601       data.size(), const_cast<T *>(&data[0]), &status);
00602   check_errors();
00603   }
00604 
00605 template void fitshandle::write_subimage (const arr<float> &data, long offset);
00606 template void fitshandle::write_subimage
00607   (const arr<double> &data, long offset);
00608 template void fitshandle::write_subimage (const arr<int> &data, long offset);
00609 
00610 template<typename T> void fitshandle::read_image (arr2<T> &data)
00611   {
00612   assert_image_hdu("fitshandle::read_image()");
00613   planck_assert (axes_.size()==2, "wrong number of dimensions");
00614   data.alloc(axes_[0], axes_[1]);
00615   fits_read_img (fptr, FITSUTIL<T>::DTYPE, 1, axes_[0]*axes_[1], 0,
00616     &data[0][0], 0, &status);
00617   check_errors();
00618   }
00619 
00620 template void fitshandle::read_image (arr2<float> &data);
00621 template void fitshandle::read_image (arr2<double> &data);
00622 template void fitshandle::read_image (arr2<int> &data);
00623 
00624 template<typename T> void fitshandle::read_subimage
00625   (arr2<T> &data, int xl, int yl)
00626   {
00627   assert_image_hdu("fitshandle::read_subimage()");
00628   planck_assert (axes_.size()==2, "wrong number of dimensions");
00629   for (int m=0; m<data.size1(); ++m)
00630     fits_read_img (fptr, FITSUTIL<T>::DTYPE, (xl+m)*axes_[1]+yl+1,
00631       data.size2(), 0, &data[m][0], 0, &status);
00632   check_errors();
00633   }
00634 
00635 template void fitshandle::read_subimage (arr2<float> &data, int xl, int yl);
00636 template void fitshandle::read_subimage (arr2<double> &data, int xl, int yl);
00637 template void fitshandle::read_subimage (arr2<int> &data, int xl, int yl);
00638 
00639 template<typename T> void fitshandle::read_subimage (arr<T> &data, long offset)
00640   {
00641   assert_image_hdu("fitshandle::read_subimage()");
00642   fits_read_img (fptr, FITSUTIL<T>::DTYPE, 1+offset,
00643       data.size(), 0, &data[0], 0, &status);
00644   check_errors();
00645   }
00646 
00647 template void fitshandle::read_subimage (arr<float> &data, long offset);
00648 template void fitshandle::read_subimage (arr<double> &data, long offset);
00649 template void fitshandle::read_subimage (arr<int> &data, long offset);
00650 
00651 void fitshandle::add_healpix_keys (int datasize)
00652   {
00653   int nside = isqrt(datasize/12);
00654   planck_assert (12*nside*nside==datasize, "Wrong Healpix map size");
00655 
00656   add_key ("PIXTYPE",string("HEALPIX"),"HEALPIX pixelisation");
00657   add_key ("ORDERING",string("RING"),
00658            "Pixel ordering scheme, either RING or NESTED");
00659   add_key ("NSIDE",nside,"Resolution parameter for HEALPIX");
00660   add_key ("FIRSTPIX",0,"First pixel # (0 based)");
00661   add_key ("LASTPIX",datasize-1,"Last pixel # (0 based)");
00662   add_key ("INDXSCHM",string("IMPLICIT"), "Indexing : IMPLICIT or EXPLICIT");
00663   add_key ("GRAIN",0,"Grain of pixel indexing");
00664   }

Generated on Fri Jul 8 09:37:14 2005 for LevelS C++ support library