// @(#)root/graf2d:$Id$
// Author:  Claudi Martinez, July 19th 2010

/*************************************************************************
 * Copyright (C) 1995-2010, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/

// TFITS is an interface that lets you reading Flexible Image Transport System
// (FITS) files, which are generally used in astronomy. This file format
// was standardized 1981 and today is still widely used among professional
// and amateur astronomers. FITS is not only an image file, but also
// it can contain spectrums, data tables, histograms, and multidimensional
// data. Furthermore, FITS data can be described itself by containing
// human-readable information that let us to interpret the data within
// the FITS file. For example, a FITS could contain a 3D data cube,
// but an additional description would tell us that we must read it, for
// example, as a 3-layer image.
//
// TFITS requires CFITSIO library to be installed on your system. It
// is currently maintained by NASA/GSFC and can be downloaded from
// http://fits.gsfc.nasa.gov, as well as documentation.
//
// Using this interface is easy and straightforward. There is only 1 class
// called "TFITSHDU" which has several methods to extract data from a
// FITS file, more specifically, from an HDU within the file. An HDU, or
// Header Data Unit, is a chunk of data with a header containing several
// "keyword = value" tokens. The header describes the structure of data
// within the HDU. An HDU can be of two types: an "image HDU" or a "table
// HDU". The former can be any kind of multidimensional array of real numbers,
// by which the name "image" may be confusing: you can store an image, but
// you can also store a N-dimensional data cube. On the other hand, table
// HDUs are sets of several rows and columns (a.k.a fields) which contain
// generic data, as strings, real or complex numbers and even arrays.
//
// Please have a look to the tutorials ($ROOTSYS/tutorials/fitsio/) to see
// some examples. IMPORTANT: to run tutorials it is required that
// you change the current working directory of ROOT (CINT) shell to the
// tutorials directory. Example:
//
// root [1] gSystem->ChangeDirectory("tutorials/fitsio")
// root [1] .x FITS_tutorial1.C

// LIST OF TODO
// - Support for complex values within data tables
// - Support for reading arrays from table cells
// - Support for grouping


// IMPLEMENTATION NOTES:
// CFITSIO library uses standard C types ('int', 'long', ...). To avoid
// confusion, the same types are used internally by the access methods.
// However, class's fields are ROOT-defined types.

//______________________________________________________________________________
/* Begin_Html
<center><h2>FITS file interface class</h2></center>
TFITS is a class that allows extracting images and data from FITS files and contains
several methods to manage them.
End_Html */

#include "TFITS.h"
#include "TROOT.h"
#include "TImage.h"
#include "TArrayI.h"
#include "TArrayD.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TH3D.h"
#include "TVectorD.h"
#include "TMatrixD.h"
#include "TObjArray.h"
#include "TObjString.h"
#include "TCanvas.h"

#include "fitsio.h"
#include <stdlib.h>


ClassImp(TFITSHDU)

/**************************************************************************/
// TFITSHDU
/**************************************************************************/

//______________________________________________________________________________
void TFITSHDU::CleanFilePath(const char *filepath_with_filter, TString &dst)
{
   // Clean path from possible filter and put the result in 'dst'.

   dst = filepath_with_filter;

   Ssiz_t ndx = dst.Index("[", 1, 0, TString::kExact);
   if (ndx != kNPOS) {
      dst.Resize(ndx);
   }
}


//______________________________________________________________________________
TFITSHDU::TFITSHDU(const char *filepath_with_filter)
{
   // TFITSHDU constructor from file path with HDU selection filter.
   // Please refer to CFITSIO manual for more information about
   // HDU selection filters.
   // Examples:
   // - TFITSHDU("/path/to/myfile.fits"): just open the PRIMARY HDU
   // - TFITSHDU("/path/to/myfile.fits[1]"): open HDU #1
   // - TFITSHDU("/path/to/myfile.fits[PICS]"): open HDU called 'PICS'
   // - TFITSHDU("/path/to/myfile.fits[ACQ][EXPOSURE > 5]"): open the (table) HDU called 'ACQ' and
   //                                                        selects the rows that have column 'EXPOSURE'
   //                                                        greater than 5.

   _initialize_me();

   fFilePath = filepath_with_filter;
   CleanFilePath(filepath_with_filter, fBaseFilePath);

   if (kFALSE == LoadHDU(fFilePath)) {
      _release_resources();
      throw -1;
   }
}

//______________________________________________________________________________
TFITSHDU::TFITSHDU(const char *filepath, Int_t extension_number)
{
   // TFITSHDU constructor from filepath and extension number.

   _initialize_me();
   CleanFilePath(filepath, fBaseFilePath);

   //Add "by extension number" filter
   fFilePath.Form("%s[%d]", fBaseFilePath.Data(), extension_number);

   if (kFALSE == LoadHDU(fFilePath)) {
      _release_resources();
      throw -1;
   }
}

//______________________________________________________________________________
TFITSHDU::TFITSHDU(const char *filepath, const char *extension_name)
{
   // TFITSHDU constructor from filepath and extension name.

   _initialize_me();
   CleanFilePath(filepath, fBaseFilePath);

   //Add "by extension number" filter
   fFilePath.Form("%s[%s]", fBaseFilePath.Data(), extension_name);


   if (kFALSE == LoadHDU(fFilePath)) {
      _release_resources();
      throw -1;
   }
}

//______________________________________________________________________________
TFITSHDU::~TFITSHDU()
{
   // TFITSHDU destructor.

   _release_resources();
}

//______________________________________________________________________________
void TFITSHDU::_release_resources()
{
   // Release internal resources.

   if (fRecords) delete [] fRecords;

   if (fType == kImageHDU) {
      if (fSizes)  delete fSizes;
      if (fPixels) delete fPixels;
   } else {
      if (fColumnsInfo) {
         if (fCells) {
            for (Int_t i = 0; i < fNColumns; i++) {
               if (fColumnsInfo[i].fType == kString) {
                  //Deallocate character arrays allocated for kString columns
                  Int_t offset = i * fNRows;
                  for (Int_t row = 0; row < fNRows; row++) {
                     delete [] fCells[offset+row].fString;
                  }
               } else if (fColumnsInfo[i].fType == kRealVector) {
                  //Deallocate character arrays allocated for kString columns
                  Int_t offset = i * fNRows;
                  for (Int_t row = 0; row < fNRows; row++) {
                     delete [] fCells[offset+row].fRealVector;
                  }
               }
            }

            delete [] fCells;
         }

         delete [] fColumnsInfo;
      }


   }
}

//______________________________________________________________________________
void TFITSHDU::_initialize_me()
{
   // Do some initializations.

   fRecords = 0;
   fPixels = 0;
   fSizes = 0;
   fColumnsInfo = 0;
   fNColumns = fNRows = 0;
   fCells = 0;
}

//______________________________________________________________________________
Bool_t TFITSHDU::LoadHDU(TString& filepath_filter)
{
   // Load HDU from fits file satisfying the specified filter.
   // Returns kTRUE if success. Otherwise kFALSE.
   // If filter == "" then the primary array is selected

   fitsfile *fp=0;
   int status = 0;
   char errdescr[FLEN_STATUS+1];

   // Open file with filter
   fits_open_file(&fp, filepath_filter.Data(), READONLY, &status);
   if (status) goto ERR;

   // Read HDU number
   int hdunum;
   fits_get_hdu_num(fp, &hdunum);
   fNumber = Int_t(hdunum);

   // Read HDU type
   int hdutype;
   fits_get_hdu_type(fp, &hdutype, &status);
   if (status) goto ERR;
   fType = (hdutype == IMAGE_HDU) ? kImageHDU : kTableHDU;

   //Read HDU header records
   int nkeys, morekeys;
   char keyname[FLEN_KEYWORD+1];
   char keyvalue[FLEN_VALUE+1];
   char comment[FLEN_COMMENT+1];

   fits_get_hdrspace(fp, &nkeys, &morekeys, &status);
   if (status) goto ERR;

   fRecords = new struct HDURecord[nkeys];

   for (int i = 1; i <= nkeys; i++) {
      fits_read_keyn(fp, i, keyname, keyvalue, comment, &status);
      if (status) goto ERR;
      fRecords[i-1].fKeyword = keyname;
      fRecords[i-1].fValue = keyvalue;
      fRecords[i-1].fComment = comment;
   }

   fNRecords = Int_t(nkeys);

   //Set extension name
   fExtensionName = "PRIMARY"; //Default
   for (int i = 0; i < nkeys; i++) {
      if (fRecords[i].fKeyword == "EXTNAME") {
         fExtensionName = fRecords[i].fValue;
         break;
      }
   }

   //Read HDU's data
   if (fType == kImageHDU) {
      //Image
      int param_ndims=0;
      long *param_dimsizes;

      //Read image number of dimensions
      fits_get_img_dim(fp, &param_ndims, &status);
      if (status) goto ERR;
      if (param_ndims > 0) {
         //Read image sizes in each dimension
         param_dimsizes = new long[param_ndims];
         fits_get_img_size(fp, param_ndims, param_dimsizes, &status);
         if (status) {
            delete [] param_dimsizes;
            goto ERR;
         }

         fSizes = new TArrayI(param_ndims);
         fSizes = new TArrayI(param_ndims);
         for (int i = 0; i < param_ndims; i++) { //Use for loop to copy values instead of passing array to constructor, since 'Int_t' size may differ from 'long' size
            fSizes->SetAt(param_dimsizes[i], i);
         }

         delete [] param_dimsizes;

         //Read pixels
         int anynul;
         long *firstpixel = new long[param_ndims];
         double nulval = 0;
         long npixels = 1;

         for (int i = 0; i < param_ndims; i++) {
            npixels *= (long) fSizes->GetAt(i); //Compute total number of pixels
            firstpixel[i] = 1; //Set first pixel to read from.
         }

         double *pixels = new double[npixels];

         fits_read_pix(fp, TDOUBLE, firstpixel, npixels,
                     (void *) &nulval, (void *) pixels, &anynul, &status);

         if (status) {
            delete [] firstpixel;
            delete [] pixels;
            goto ERR;
         }

         fPixels = new TArrayD(npixels, pixels);

         delete [] firstpixel;
         delete [] pixels;

      } else {
         //Null array
         fSizes = new TArrayI();
         fPixels = new TArrayD();
      }
   } else {
      // Table

      // Get table's number of rows and columns
      long table_rows;
      int  table_cols;

      fits_get_num_rows(fp, &table_rows, &status);
      if (status) goto ERR;

      fNRows = Int_t(table_rows);

      fits_get_num_cols(fp, &table_cols, &status);
      if (status) goto ERR;

      fNColumns = Int_t(table_cols);

      // Allocate column info array
      fColumnsInfo = new struct Column[table_cols];

      // Read column names
      char colname[80];
      int colnum;

      fits_get_colname(fp, CASEINSEN, (char*) "*", colname, &colnum, &status);
      while (status == COL_NOT_UNIQUE)
      {
         fColumnsInfo[colnum-1].fName = colname;
         fits_get_colname(fp, CASEINSEN, (char*) "*", colname, &colnum, &status);
      }
      if (status != COL_NOT_FOUND) goto ERR;
      status = 0;

      //Allocate cells
      fCells = new union Cell [table_rows * table_cols];

      // Read columns
      int typecode;
      long repeat, width;
      Int_t cellindex;


      for (colnum = 0, cellindex = 0; colnum < fNColumns; colnum++) {
         fits_get_coltype(fp, colnum+1, &typecode, &repeat, &width, &status);
         if (status) goto ERR;

         if ((typecode == TDOUBLE) || (typecode == TSHORT) || (typecode == TLONG)
                                   || (typecode == TFLOAT) || (typecode == TLOGICAL) || (typecode == TBIT)
                                   || (typecode == TBYTE)  || (typecode == TSTRING)) {

            fColumnsInfo[colnum].fType = (typecode == TSTRING) ? kString : kRealNumber;

            if (typecode == TSTRING) {
               // String column
               int dispwidth=0;
               fits_get_col_display_width(fp, colnum+1, &dispwidth, &status);
               if (status) goto ERR;


               char *nulval = (char*) "";
               int anynul=0;
               char **array;

               if (dispwidth <= 0) {
                  dispwidth = 1;
               }

               array = new char* [table_rows];
               for (long row = 0; row < table_rows; row++) {
                  array[row] = new char[dispwidth+1]; //also room for end null!
               }

               if (repeat > 0) {
                  fits_read_col(fp, TSTRING, colnum+1, 1, 1, table_rows, nulval, array, &anynul, &status);
                  if (status) {
                     for (long row = 0; row < table_rows; row++) {
                        delete [] array[row];
                     }
                     delete [] array;
                     goto ERR;
                  }

               } else {
                  //No elements: set dummy
                  for (long row = 0; row < table_rows; row++) {
                     strlcpy(array[row], "-",dispwidth+1);
                  }
               }


               //Save values
               for (long row = 0; row < table_rows; row++) {
                  fCells[cellindex++].fString = array[row];
               }

               delete [] array; //Delete temporal array holding pointer to strings, but not delete strings themselves!


            } else {
               //Numeric or vector column
               double nulval = 0;
               int anynul=0;

               fColumnsInfo[colnum].fDim = (Int_t) repeat;

               double *array;

               if (repeat > 0) {
                  array = new double [table_rows * repeat]; //Hope you got a big machine! Ask China otherwise :-)
                  fits_read_col(fp, TDOUBLE, colnum+1, 1, 1, table_rows * repeat, &nulval, array, &anynul, &status);

                  if (status) {
                     delete [] array;
                     goto ERR;
                  }
               } else {
                  //No elements: set dummy
                  array = new double [table_rows];
                  for (long row = 0; row < table_rows; row++) {
                     array[row] = 0.0;
                  }
               }

               //Save values
               if (repeat == 1) {
                  //Scalar
                  for (long row = 0; row < table_rows; row++) {
                     fCells[cellindex++].fRealNumber = array[row];
                  }
               } else if (repeat > 1) {
                  //Vector
                  for (long row = 0; row < table_rows; row++) {
                     double *vec = new double [repeat];
                     long offset = row * repeat;
                     for (long component = 0; component < repeat; component++) {
                        vec[component] = array[offset++];
                     }
                     fCells[cellindex++].fRealVector = vec;
                  }
               }

               delete [] array;

            }

         } else {
            Warning("LoadHDU", "error opening FITS file. Column type %d is currently not supported", typecode);
         }
      }



      if (hdutype == ASCII_TBL) {
         //ASCII table

      } else {
         //Binary table
      }

   }

   // Close file
   fits_close_file(fp, &status);
   return kTRUE;

ERR:
   fits_get_errstatus(status, errdescr);
   Warning("LoadHDU", "error opening FITS file. Details: %s", errdescr);
   status = 0;
   if (fp) fits_close_file(fp, &status);
   return kFALSE;
}

//______________________________________________________________________________
struct TFITSHDU::HDURecord* TFITSHDU::GetRecord(const char *keyword)
{
   // Get record by keyword.

   for (int i = 0; i < fNRecords; i++) {
      if (fRecords[i].fKeyword == keyword) {
         return &fRecords[i];
      }
   }
   return 0;
}

//______________________________________________________________________________
TString& TFITSHDU::GetKeywordValue(const char *keyword)
{
   // Get the value of a given keyword. Return "" if not found.

   HDURecord *rec = GetRecord(keyword);
   if (rec) {
      return rec->fValue;
   } else {
      return *(new TString(""));
   }
}

//______________________________________________________________________________
void TFITSHDU::PrintHDUMetadata(const Option_t *) const
{
   // Print records.

   for (int i = 0; i < fNRecords; i++) {
      if (fRecords[i].fComment.Length() > 0) {
         printf("%-10s = %s / %s\n", fRecords[i].fKeyword.Data(), fRecords[i].fValue.Data(), fRecords[i].fComment.Data());
      } else {
         printf("%-10s = %s\n", fRecords[i].fKeyword.Data(), fRecords[i].fValue.Data());
      }
   }
}

//______________________________________________________________________________
void TFITSHDU::PrintFileMetadata(const Option_t *opt) const
{
   // Print HDU's parent file's metadata.

   fitsfile *fp=0;
   int status = 0;
   char errdescr[FLEN_STATUS+1];
   int hducount, extnum;
   int hdutype = IMAGE_HDU;
   const char *exttype;
   char extname[FLEN_CARD]="PRIMARY"; //room enough
   int verbose = (opt[0] ? 1 : 0);

   // Open file with no filters: current HDU will be the primary one.
   fits_open_file(&fp, fBaseFilePath.Data(), READONLY, &status);
   if (status) goto ERR;

   // Read HDU count
   fits_get_num_hdus(fp, &hducount, &status);
   if (status) goto ERR;
   printf("Total: %d HDUs\n", hducount);

   extnum = 0;
   while(hducount) {
      // Read HDU type
      fits_get_hdu_type(fp, &hdutype, &status);
      if (status) goto ERR;

      if (hdutype == IMAGE_HDU) {
         exttype="IMAGE";
      } else if (hdutype == ASCII_TBL) {
         exttype="ASCII TABLE";
      } else {
         exttype="BINARY TABLE";
      }

      //Read HDU header records
      int nkeys, morekeys;
      char keyname[FLEN_KEYWORD+1];
      char keyvalue[FLEN_VALUE+1];
      char comment[FLEN_COMMENT+1];

      fits_get_hdrspace(fp, &nkeys, &morekeys, &status);
      if (status) goto ERR;

      struct HDURecord *records = new struct HDURecord[nkeys];

      for (int i = 1; i <= nkeys; i++) {
         fits_read_keyn(fp, i, keyname, keyvalue, comment, &status);
         if (status) {
            delete [] records;
            goto ERR;
         }

         records[i-1].fKeyword = keyname;
         records[i-1].fValue = keyvalue;
         records[i-1].fComment = comment;

         if (strcmp(keyname, "EXTNAME") == 0) {
            //Extension name
            strlcpy(extname, keyvalue,FLEN_CARD);
         }
      }

      //HDU info
      printf("   [%d] %s (%s)\n", extnum, exttype, extname);

      //HDU records
      if (verbose) {
         for (int i = 0; i < nkeys; i++) {
            if (comment[0]) {
               printf("      %-10s = %s / %s\n", records[i].fKeyword.Data(), records[i].fValue.Data(), records[i].fComment.Data());
            } else {
               printf("      %-10s = %s\n", records[i].fKeyword.Data(), records[i].fValue.Data());
            }
         }
      }
      printf("\n");

      delete [] records;

      hducount--;
      extnum++;
      if (hducount){
         fits_movrel_hdu(fp, 1, &hdutype, &status);
         if (status) goto ERR;
      }
   }

   // Close file
   fits_close_file(fp, &status);
   return;

ERR:
   fits_get_errstatus(status, errdescr);
   Warning("PrintFileMetadata", "error opening FITS file. Details: %s", errdescr);
   status = 0;
   if (fp) fits_close_file(fp, &status);
}

//______________________________________________________________________________
void TFITSHDU::PrintColumnInfo(const Option_t *) const
{
   // Print column information

   if (fType != kTableHDU) {
      Warning("PrintColumnInfo", "this is not a table HDU.");
      return;
   }

   for (Int_t i = 0; i < fNColumns; i++) {
      printf("%-20s : %s\n", fColumnsInfo[i].fName.Data(), (fColumnsInfo[i].fType == kRealNumber) ? "REAL NUMBER" : "STRING");
   }
}

//______________________________________________________________________________
void TFITSHDU::PrintFullTable(const Option_t *) const
{
   // Print full table contents

   int printed_chars;

   if (fType != kTableHDU) {
      Warning("PrintColumnInfo", "this is not a table HDU.");
      return;
   }

   // Dump header
   putchar('\n');
   printed_chars = 0;
   for (Int_t col = 0; col < fNColumns; col++) {
      printed_chars += printf("%-10s| ", fColumnsInfo[col].fName.Data());
   }
   putchar('\n');
   while(printed_chars--) {
      putchar('-');
   }
   putchar('\n');

   // Dump rows
   for (Int_t row = 0; row < fNRows; row++) {
      for (Int_t col = 0; col < fNColumns; col++) {
         if (fColumnsInfo[col].fType == kString) {
            printf("%-10s", fCells[col * fNRows + row].fString);
         } else {
            printed_chars = printf("%.2lg", fCells[col * fNRows + row].fRealNumber);
            printed_chars -= 10;
            while (printed_chars < 0) {
               putchar(' ');
               printed_chars++;
            }
         }

         if (col <= fNColumns - 1) printf("| ");
      }
      printf("\n");
   }
}

//______________________________________________________________________________
void TFITSHDU::Print(const Option_t *opt) const
{
   // Print metadata.
   // Currently supported options:
   // ""  :  print HDU record data
   // "F" :  print FITS file's extension names, numbers and types
   // "F+":  print FITS file's extension names and types and their record data
   // "T" :  print column information when HDU is a table
   // "T+" : print full table (columns header and rows)

   if ((opt[0] == 'F') || (opt[0] == 'f')) {
      PrintFileMetadata((opt[1] == '+') ? "+" : "");
   } else if ((opt[0] == 'T') || (opt[0] == 't')) {
      if (opt[1] == '+') {
         PrintFullTable("");
      } else {
         PrintColumnInfo("");
      }

   } else {
      PrintHDUMetadata("");
   }
}


//______________________________________________________________________________
TImage *TFITSHDU::ReadAsImage(Int_t layer, TImagePalette *pal)
{
   // Read image HDU as a displayable image. Return 0 if conversion cannot be done.
   // If the HDU seems to be a multilayer image, 'layer' parameter can be used
   // to retrieve the specified layer (starting from 0)

   if (fType != kImageHDU) {
      Warning("ReadAsImage", "this is not an image HDU.");
      return 0;
   }


   if (((fSizes->GetSize() != 2) && (fSizes->GetSize() != 3) && (fSizes->GetSize() != 4)) || ((fSizes->GetSize() == 4) && (fSizes->GetAt(3) > 1))) {
      Warning("ReadAsImage", "could not convert image HDU to image because it has %d dimensions.", fSizes->GetSize());
      return 0;
   }

   Int_t width, height;
   UInt_t pixels_per_layer;

   width  = Int_t(fSizes->GetAt(0));
   height = Int_t(fSizes->GetAt(1));

   pixels_per_layer = UInt_t(width) * UInt_t(height);

   if (  ((fSizes->GetSize() == 2) && (layer > 0))
      || (((fSizes->GetSize() == 3) || (fSizes->GetSize() == 4)) && (layer >= fSizes->GetAt(2)))) {

      Warning("ReadAsImage", "layer out of bounds.");
      return 0;
   }

   // Get the maximum and minimum pixel values in the layer to auto-stretch pixels
   Double_t maxval = 0, minval = 0;
   UInt_t i;
   Double_t pixvalue;
   Int_t offset = layer * pixels_per_layer;

   for (i = 0; i < pixels_per_layer; i++) {
      pixvalue = fPixels->GetAt(offset + i);

      if (pixvalue > maxval) {
         maxval = pixvalue;
      }

      if ((i == 0) || (pixvalue < minval)) {
         minval = pixvalue;
      }
   }

   //Build the image stretching pixels into a range from 0.0 to 255.0
   //TImage *im = new TImage(width, height);
   TImage *im = TImage::Create();
   if (!im) return 0;
   TArrayD *layer_pixels = new TArrayD(pixels_per_layer);


   if (maxval == minval) {
      //plain image
      for (i = 0; i < pixels_per_layer; i++) {
         layer_pixels->SetAt(255.0, i);
      }
   } else {
      Double_t factor = 255.0 / (maxval-minval);
      for (i = 0; i < pixels_per_layer; i++) {
         pixvalue = fPixels->GetAt(offset + i);
         layer_pixels->SetAt(factor * (pixvalue-minval), i) ;
      }
   }

   if (pal == 0) {
      // Default palette: grayscale palette
      pal = new TImagePalette(256);
      for (i = 0; i < 256; i++) {
         pal->fPoints[i] = ((Double_t)i)/255.0;
         pal->fColorAlpha[i] = 0xffff;
         pal->fColorBlue[i] = pal->fColorGreen[i] = pal->fColorRed[i] = i << 8;
      }
      pal->fPoints[0] = 0;
      pal->fPoints[255] = 1.0;

      im->SetImage(*layer_pixels, UInt_t(width), pal);

      delete pal;

   } else {
      im->SetImage(*layer_pixels, UInt_t(width), pal);
   }

   delete layer_pixels;

   return im;
}

//______________________________________________________________________________
void TFITSHDU::Draw(Option_t *)
{
   // If the HDU is an image, draw the first layer of the primary array
   // To set a title to the canvas, pass it in "opt"

   if (fType != kImageHDU) {
      Warning("Draw", "cannot draw. This is not an image HDU.");
      return;
   }

   TImage *im = ReadAsImage(0, 0);
   if (im) {
      Int_t width = Int_t(fSizes->GetAt(0));
      Int_t height = Int_t(fSizes->GetAt(1));
      TString cname, ctitle;
      cname.Form("%sHDU", this->GetName());
      ctitle.Form("%d x %d", width, height);
      new TCanvas(cname, ctitle, width, height);
      im->Draw();
   }
}


//______________________________________________________________________________
TMatrixD* TFITSHDU::ReadAsMatrix(Int_t layer, Option_t *opt)
{
   // Read image HDU as a matrix. Return 0 if conversion cannot be done
   // If the HDU seems to be a multilayer image, 'layer' parameter can be used
   // to retrieve the specified layer (starting from 0) in matrix form.
   // Options (value of 'opt'):
   // "S": stretch pixel values to a range from 0.0 to 1.0

   if (fType != kImageHDU) {
      Warning("ReadAsMatrix", "this is not an image HDU.");
      return 0;
   }

   if (((fSizes->GetSize() != 2) && (fSizes->GetSize() != 3) && (fSizes->GetSize() != 4)) || ((fSizes->GetSize() == 4) && (fSizes->GetAt(3) > 1))) {
      Warning("ReadAsMatrix", "could not convert image HDU to image because it has %d dimensions.", fSizes->GetSize());
      return 0;
   }


   if (   ((fSizes->GetSize() == 2) && (layer > 0))
       || (((fSizes->GetSize() == 3) || (fSizes->GetSize() == 4)) && (layer >= fSizes->GetAt(2)))) {
      Warning("ReadAsMatrix", "layer out of bounds.");
      return 0;
   }

   Int_t width, height;
   UInt_t pixels_per_layer;
   Int_t offset;
   UInt_t i;
   TMatrixD *mat=0;

   width  = Int_t(fSizes->GetAt(0));
   height = Int_t(fSizes->GetAt(1));

   pixels_per_layer = UInt_t(width) * UInt_t(height);
   offset = layer * pixels_per_layer;

   double *layer_pixels = new double[pixels_per_layer];

   if ((opt[0] == 'S') || (opt[0] == 's')) {
      //Stretch
      // Get the maximum and minimum pixel values in the layer to auto-stretch pixels
      Double_t factor, maxval=0, minval=0;
      Double_t pixvalue;
      for (i = 0; i < pixels_per_layer; i++) {
         pixvalue = fPixels->GetAt(offset + i);

         if (pixvalue > maxval) {
            maxval = pixvalue;
         }

         if ((i == 0) || (pixvalue < minval)) {
            minval = pixvalue;
         }
      }

      if (maxval == minval) {
         //plain image
         for (i = 0; i < pixels_per_layer; i++) {
            layer_pixels[i] = 1.0;
         }
      } else {
         factor = 1.0 / (maxval-minval);
         mat = new TMatrixD(height, width);

         for (i = 0; i < pixels_per_layer; i++) {
            layer_pixels[i] = factor * (fPixels->GetAt(offset + i) - minval);
         }
      }

   } else {
      //No stretching
      mat = new TMatrixD(height, width);

      for (i = 0; i < pixels_per_layer; i++) {
         layer_pixels[i] = fPixels->GetAt(offset + i);
      }
   }

   if (mat) {
      // mat->Use(height, width, layer_pixels);
      memcpy(mat->GetMatrixArray(), layer_pixels, pixels_per_layer*sizeof(double));
   }

   delete [] layer_pixels;
   return mat;
}


//______________________________________________________________________________
TH1 *TFITSHDU::ReadAsHistogram()
{
   // Read image HDU as a histogram. Return 0 if conversion cannot be done.
   // The returned object can be TH1D, TH2D or TH3D depending on data dimensionality.
   // Please, check condition (returnedValue->IsA() == TH*D::Class()) to
   // determine the object class.
   // NOTE: do not confuse with image histogram! This function interprets
   // the array as a histogram. It does not compute the histogram of pixel
   // values of an image! Here "pixels" are interpreted as number of entries.

   if (fType != kImageHDU) {
      Warning("ReadAsHistogram", "this is not an image HDU.");
      return 0;
   }

   TH1 *result=0;

   if ((fSizes->GetSize() != 1) && (fSizes->GetSize() != 2) && (fSizes->GetSize() != 3)) {
      Warning("ReadAsHistogram", "could not convert image HDU to histogram because it has %d dimensions.", fSizes->GetSize());
      return 0;
   }

   if (fSizes->GetSize() == 1) {
      //1D
      UInt_t Nx = UInt_t(fSizes->GetAt(0));
      UInt_t x;

      TH1D *h = new TH1D("", "", Int_t(Nx), 0, Nx-1);

      for (x = 0; x < Nx; x++) {
         Int_t nentries = Int_t(fPixels->GetAt(x));
         if (nentries < 0) nentries = 0; //Crop negative values
         h->Fill(x, nentries);
      }

      result = h;

   } else if (fSizes->GetSize() == 2) {
      //2D
      UInt_t Nx = UInt_t(fSizes->GetAt(0));
      UInt_t Ny = UInt_t(fSizes->GetAt(1));
      UInt_t x,y;

      TH2D *h = new TH2D("", "", Int_t(Nx), 0, Nx-1, Int_t(Ny), 0, Ny-1);

      for (y = 0; y < Ny; y++) {
         UInt_t offset = y * Nx;
         for (x = 0; x < Nx; x++) {
            Int_t nentries = Int_t(fPixels->GetAt(offset + x));
            if (nentries < 0) nentries = 0; //Crop negative values
            h->Fill(x,y, nentries);
         }
      }

      result = h;

   } else if (fSizes->GetSize() == 3) {
      //3D
      UInt_t Nx = UInt_t(fSizes->GetAt(0));
      UInt_t Ny = UInt_t(fSizes->GetAt(1));
      UInt_t Nz = UInt_t(fSizes->GetAt(2));
      UInt_t x,y,z;

      TH3D *h = new TH3D("", "", Int_t(Nx), 0, Nx-1, Int_t(Ny), 0, Ny-1, Int_t(Nz), 0, Nz-1);


      for (z = 0; z < Nz; z++) {
         UInt_t offset1 = z * Nx * Ny;
         for (y = 0; y < Ny; y++) {
            UInt_t offset2 = y * Nx;
            for (x = 0; x < Nx; x++) {
               Int_t nentries = Int_t(fPixels->GetAt(offset1 + offset2 + x));
               if (nentries < 0) nentries = 0; //Crop negative values
               h->Fill(x, y, z, nentries);
            }
         }
      }

      result = h;
   }

   return result;
}

//______________________________________________________________________________
TVectorD* TFITSHDU::GetArrayRow(UInt_t row)
{
   // Get a row from the image HDU when it's a 2D array.

   if (fType != kImageHDU) {
      Warning("GetArrayRow", "this is not an image HDU.");
      return 0;
   }

   if (fSizes->GetSize() != 2) {
      Warning("GetArrayRow", "could not get row from HDU because it has %d dimensions.", fSizes->GetSize());
      return 0;
   }

   UInt_t i, offset, W,H;

   W =  UInt_t(fSizes->GetAt(0));
   H =  UInt_t(fSizes->GetAt(1));


   if (row >= H) {
      Warning("GetArrayRow", "index out of bounds.");
      return 0;
   }

   offset = W * row;
   double *v = new double[W];

   for (i = 0; i < W; i++) {
      v[i] = fPixels->GetAt(offset+i);
   }

   TVectorD *vec = new TVectorD(W, v);

   delete [] v;

   return vec;
}

//______________________________________________________________________________
TVectorD* TFITSHDU::GetArrayColumn(UInt_t col)
{
   // Get a column from the image HDU when it's a 2D array.

   if (fType != kImageHDU) {
      Warning("GetArrayColumn", "this is not an image HDU.");
      return 0;
   }

   if (fSizes->GetSize() != 2) {
      Warning("GetArrayColumn", "could not get row from HDU because it has %d dimensions.", fSizes->GetSize());
      return 0;
   }

   UInt_t i, W,H;

   W =  UInt_t(fSizes->GetAt(0));
   H =  UInt_t(fSizes->GetAt(1));


   if (col >= W) {
      Warning("GetArrayColumn", "index out of bounds.");
      return 0;
   }

   double *v = new double[H];

   for (i = 0; i < H; i++) {
      v[i] = fPixels->GetAt(W*i+col);
   }

   TVectorD *vec = new TVectorD(H, v);

   delete [] v;

   return vec;
}


//______________________________________________________________________________
Int_t TFITSHDU::GetColumnNumber(const char *colname)
{
   //Get column number given its name

   Int_t colnum;
   for (colnum = 0; colnum < fNColumns; colnum++) {
      if (fColumnsInfo[colnum].fName == colname) {
         return colnum;
      }
   }
   return -1;
}

//______________________________________________________________________________
TObjArray* TFITSHDU::GetTabStringColumn(Int_t colnum)
{
   // Get a string-typed column from a table HDU given its column index (>=0).

   if (fType != kTableHDU) {
      Warning("GetTabStringColumn", "this is not a table HDU.");
      return 0;
   }

   if ((colnum < 0) || (colnum >= fNColumns)) {
      Warning("GetTabStringColumn", "column index out of bounds.");
      return 0;
   }

   if (fColumnsInfo[colnum].fType != kString) {
      Warning("GetTabStringColumn", "attempting to read a column that is not of type 'kString'.");
      return 0;
   }

   Int_t offset = colnum * fNRows;

   TObjArray *res = new TObjArray();
   for (Int_t row = 0; row < fNRows; row++) {
      res->Add(new TObjString(fCells[offset + row].fString));
   }

   return res;
}

//______________________________________________________________________________
TObjArray* TFITSHDU::GetTabStringColumn(const char *colname)
{
   // Get a string-typed column from a table HDU given its name

   if (fType != kTableHDU) {
      Warning("GetTabStringColumn", "this is not a table HDU.");
      return 0;
   }


   Int_t colnum = GetColumnNumber(colname);

   if (colnum == -1) {
      Warning("GetTabStringColumn", "column not found.");
      return 0;
   }

   if (fColumnsInfo[colnum].fType != kString) {
      Warning("GetTabStringColumn", "attempting to read a column that is not of type 'kString'.");
      return 0;
   }

   Int_t offset = colnum * fNRows;

   TObjArray *res = new TObjArray();
   for (Int_t row = 0; row < fNRows; row++) {
      res->Add(new TObjString(fCells[offset + row].fString));
   }

   return res;
}

//______________________________________________________________________________
TVectorD* TFITSHDU::GetTabRealVectorColumn(Int_t colnum)
{
   // Get a real number-typed column from a table HDU given its column index (>=0).

   if (fType != kTableHDU) {
      Warning("GetTabRealVectorColumn", "this is not a table HDU.");
      return 0;
   }

   if ((colnum < 0) || (colnum >= fNColumns)) {
      Warning("GetTabRealVectorColumn", "column index out of bounds.");
      return 0;
   }

   if (fColumnsInfo[colnum].fType != kRealNumber) {
      Warning("GetTabRealVectorColumn", "attempting to read a column that is not of type 'kRealNumber'.");
      return 0;
   } else if (fColumnsInfo[colnum].fDim > 1) {
      Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded vectors, not real scalars. Use GetTabRealVectorCells() instead.");
      return 0;
   }

   Int_t offset = colnum * fNRows;

   Double_t *arr = new Double_t[fNRows];

   for (Int_t row = 0; row < fNRows; row++) {
      arr[row] = fCells[offset + row].fRealNumber;
   }

   TVectorD *res = new TVectorD();
   res->Use(fNRows, arr);

   return res;
}

//______________________________________________________________________________
TVectorD* TFITSHDU::GetTabRealVectorColumn(const char *colname)
{
   // Get a real number-typed column from a table HDU given its name

   if (fType != kTableHDU) {
      Warning("GetTabRealVectorColumn", "this is not a table HDU.");
      return 0;
   }

   Int_t colnum = GetColumnNumber(colname);

   if (colnum == -1) {
      Warning("GetTabRealVectorColumn", "column not found.");
      return 0;
   }

   if (fColumnsInfo[colnum].fType != kRealNumber) {
      Warning("GetTabRealVectorColumn", "attempting to read a column that is not of type 'kRealNumber'.");
      return 0;
   } else if (fColumnsInfo[colnum].fDim > 1) {
      Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded vectors, not real scalars. Use GetTabRealVectorCells() instead.");
      return 0;
   }

   Int_t offset = colnum * fNRows;

   Double_t *arr = new Double_t[fNRows];

   for (Int_t row = 0; row < fNRows; row++) {
      arr[row] = fCells[offset + row].fRealNumber;
   }

   TVectorD *res = new TVectorD();
   res->Use(fNRows, arr);

   return res;
}

//______________________________________________________________________________
Bool_t TFITSHDU::Change(const char *filter)
{
   // Change to another HDU given by "filter".
   // The parameter "filter" will be appended to the
   // FITS file's base path. For example:
   // hduObject.Change("[EVENTS][TIME > 5]");
   // Please, see documentation of TFITSHDU(const char *filepath_with_filter) constructor
   // for further information.

   TString tmppath;
   tmppath.Form("%s%s", fBaseFilePath.Data(), filter);

   _release_resources();
   _initialize_me();

   if (kFALSE == LoadHDU(tmppath)) {
      //Failed! Restore previous hdu
      Warning("Change", "error changing HDU. Restoring the previous one...");

      _release_resources();
      _initialize_me();

      if (kFALSE == LoadHDU(fFilePath)) {
         Warning("Change", "could not restore previous HDU. This object is no longer reliable!!");
      }
      return kFALSE;
   }

   //Set new full path
   fFilePath = tmppath;
   return kTRUE;
}

//______________________________________________________________________________
Bool_t TFITSHDU::Change(Int_t extension_number)
{
   // Change to another HDU given by extension_number

   TString tmppath;
   tmppath.Form("[%d]", extension_number);

   return Change(tmppath.Data());
}

//______________________________________________________________________________
TObjArray *TFITSHDU::GetTabRealVectorCells(Int_t colnum)
{
   // Get a collection of real vectors embedded in cells along a given column from a table HDU. colnum >= 0.
   if (fType != kTableHDU) {
      Warning("GetTabRealVectorCells", "this is not a table HDU.");
      return 0;
   }

   if ((colnum < 0) || (colnum >= fNColumns)) {
      Warning("GetTabRealVectorCells", "column index out of bounds.");
      return 0;
   }

   if (fColumnsInfo[colnum].fType != kRealNumber) {
      Warning("GetTabRealVectorCells", "attempting to read a column that is not of type 'kRealNumber'.");
      return 0;
   }

   Int_t offset = colnum * fNRows;

   TObjArray *res = new TObjArray();
   Int_t dim = fColumnsInfo[colnum].fDim;

   for (Int_t row = 0; row < fNRows; row++) {
      TVectorD *v = new TVectorD();
      v->Use(dim, fCells[offset + row].fRealVector);
      res->Add(v);
   }

   //Make the collection to own the allocated TVectorD objects, so when
   //destroying the collection, the vectors will be destroyed too.
   res->SetOwner(kTRUE);

   return res;
}

//______________________________________________________________________________
TObjArray *TFITSHDU::GetTabRealVectorCells(const char *colname)
{
   // Get a collection of real vectors embedded in cells along a given column from a table HDU by name
   if (fType != kTableHDU) {
      Warning("GetTabRealVectorCells", "this is not a table HDU.");
      return 0;
   }

   Int_t colnum = GetColumnNumber(colname);

   if (colnum == -1) {
      Warning("GetTabRealVectorCells", "column not found.");
      return 0;
   }

   return GetTabRealVectorCells(colnum);
}

//______________________________________________________________________________
TVectorD *TFITSHDU::GetTabRealVectorCell(Int_t rownum, Int_t colnum)
{
   // Get a real vector embedded in a cell given by (row>=0, column>=0)
   if (fType != kTableHDU) {
      Warning("GetTabRealVectorCell", "this is not a table HDU.");
      return 0;
   }

   if ((colnum < 0) || (colnum >= fNColumns)) {
      Warning("GetTabRealVectorCell", "column index out of bounds.");
      return 0;
   }

   if ((rownum < 0) || (rownum >= fNRows)) {
      Warning("GetTabRealVectorCell", "row index out of bounds.");
      return 0;
   }

   if (fColumnsInfo[colnum].fType != kRealNumber) {
      Warning("GetTabRealVectorCell", "attempting to read a column that is not of type 'kRealNumber'.");
      return 0;
   }

   TVectorD *v = new TVectorD();
   v->Use(fColumnsInfo[colnum].fDim, fCells[(colnum * fNRows) + rownum].fRealVector);

   return v;
}

//______________________________________________________________________________
TVectorD *TFITSHDU::GetTabRealVectorCell(Int_t rownum, const char *colname)
{
   // Get a real vector embedded in a cell given by (row>=0, column name)
   if (fType != kTableHDU) {
      Warning("GetTabRealVectorCell", "this is not a table HDU.");
      return 0;
   }

   Int_t colnum = GetColumnNumber(colname);

   if (colnum == -1) {
      Warning("GetTabRealVectorCell", "column not found.");
      return 0;
   }

   return GetTabRealVectorCell(rownum, colnum);
}


//______________________________________________________________________________
const TString& TFITSHDU::GetColumnName(Int_t colnum)
{
   // Get the name of a column given its index (column>=0).
   // In case of error the column name is "".

   static TString noName;
   if (fType != kTableHDU) {
      Error("GetColumnName", "this is not a table HDU.");
      return noName;
   }

   if ((colnum < 0) || (colnum >= fNColumns)) {
      Error("GetColumnName", "column index out of bounds.");
      return noName;
   }
   return fColumnsInfo[colnum].fName;
}
 TFITS.cxx:1
 TFITS.cxx:2
 TFITS.cxx:3
 TFITS.cxx:4
 TFITS.cxx:5
 TFITS.cxx:6
 TFITS.cxx:7
 TFITS.cxx:8
 TFITS.cxx:9
 TFITS.cxx:10
 TFITS.cxx:11
 TFITS.cxx:12
 TFITS.cxx:13
 TFITS.cxx:14
 TFITS.cxx:15
 TFITS.cxx:16
 TFITS.cxx:17
 TFITS.cxx:18
 TFITS.cxx:19
 TFITS.cxx:20
 TFITS.cxx:21
 TFITS.cxx:22
 TFITS.cxx:23
 TFITS.cxx:24
 TFITS.cxx:25
 TFITS.cxx:26
 TFITS.cxx:27
 TFITS.cxx:28
 TFITS.cxx:29
 TFITS.cxx:30
 TFITS.cxx:31
 TFITS.cxx:32
 TFITS.cxx:33
 TFITS.cxx:34
 TFITS.cxx:35
 TFITS.cxx:36
 TFITS.cxx:37
 TFITS.cxx:38
 TFITS.cxx:39
 TFITS.cxx:40
 TFITS.cxx:41
 TFITS.cxx:42
 TFITS.cxx:43
 TFITS.cxx:44
 TFITS.cxx:45
 TFITS.cxx:46
 TFITS.cxx:47
 TFITS.cxx:48
 TFITS.cxx:49
 TFITS.cxx:50
 TFITS.cxx:51
 TFITS.cxx:52
 TFITS.cxx:53
 TFITS.cxx:54
 TFITS.cxx:55
 TFITS.cxx:56
 TFITS.cxx:57
 TFITS.cxx:58
 TFITS.cxx:59
 TFITS.cxx:60
 TFITS.cxx:61
 TFITS.cxx:62
 TFITS.cxx:63
 TFITS.cxx:64
 TFITS.cxx:65
 TFITS.cxx:66
 TFITS.cxx:67
 TFITS.cxx:68
 TFITS.cxx:69
 TFITS.cxx:70
 TFITS.cxx:71
 TFITS.cxx:72
 TFITS.cxx:73
 TFITS.cxx:74
 TFITS.cxx:75
 TFITS.cxx:76
 TFITS.cxx:77
 TFITS.cxx:78
 TFITS.cxx:79
 TFITS.cxx:80
 TFITS.cxx:81
 TFITS.cxx:82
 TFITS.cxx:83
 TFITS.cxx:84
 TFITS.cxx:85
 TFITS.cxx:86
 TFITS.cxx:87
 TFITS.cxx:88
 TFITS.cxx:89
 TFITS.cxx:90
 TFITS.cxx:91
 TFITS.cxx:92
 TFITS.cxx:93
 TFITS.cxx:94
 TFITS.cxx:95
 TFITS.cxx:96
 TFITS.cxx:97
 TFITS.cxx:98
 TFITS.cxx:99
 TFITS.cxx:100
 TFITS.cxx:101
 TFITS.cxx:102
 TFITS.cxx:103
 TFITS.cxx:104
 TFITS.cxx:105
 TFITS.cxx:106
 TFITS.cxx:107
 TFITS.cxx:108
 TFITS.cxx:109
 TFITS.cxx:110
 TFITS.cxx:111
 TFITS.cxx:112
 TFITS.cxx:113
 TFITS.cxx:114
 TFITS.cxx:115
 TFITS.cxx:116
 TFITS.cxx:117
 TFITS.cxx:118
 TFITS.cxx:119
 TFITS.cxx:120
 TFITS.cxx:121
 TFITS.cxx:122
 TFITS.cxx:123
 TFITS.cxx:124
 TFITS.cxx:125
 TFITS.cxx:126
 TFITS.cxx:127
 TFITS.cxx:128
 TFITS.cxx:129
 TFITS.cxx:130
 TFITS.cxx:131
 TFITS.cxx:132
 TFITS.cxx:133
 TFITS.cxx:134
 TFITS.cxx:135
 TFITS.cxx:136
 TFITS.cxx:137
 TFITS.cxx:138
 TFITS.cxx:139
 TFITS.cxx:140
 TFITS.cxx:141
 TFITS.cxx:142
 TFITS.cxx:143
 TFITS.cxx:144
 TFITS.cxx:145
 TFITS.cxx:146
 TFITS.cxx:147
 TFITS.cxx:148
 TFITS.cxx:149
 TFITS.cxx:150
 TFITS.cxx:151
 TFITS.cxx:152
 TFITS.cxx:153
 TFITS.cxx:154
 TFITS.cxx:155
 TFITS.cxx:156
 TFITS.cxx:157
 TFITS.cxx:158
 TFITS.cxx:159
 TFITS.cxx:160
 TFITS.cxx:161
 TFITS.cxx:162
 TFITS.cxx:163
 TFITS.cxx:164
 TFITS.cxx:165
 TFITS.cxx:166
 TFITS.cxx:167
 TFITS.cxx:168
 TFITS.cxx:169
 TFITS.cxx:170
 TFITS.cxx:171
 TFITS.cxx:172
 TFITS.cxx:173
 TFITS.cxx:174
 TFITS.cxx:175
 TFITS.cxx:176
 TFITS.cxx:177
 TFITS.cxx:178
 TFITS.cxx:179
 TFITS.cxx:180
 TFITS.cxx:181
 TFITS.cxx:182
 TFITS.cxx:183
 TFITS.cxx:184
 TFITS.cxx:185
 TFITS.cxx:186
 TFITS.cxx:187
 TFITS.cxx:188
 TFITS.cxx:189
 TFITS.cxx:190
 TFITS.cxx:191
 TFITS.cxx:192
 TFITS.cxx:193
 TFITS.cxx:194
 TFITS.cxx:195
 TFITS.cxx:196
 TFITS.cxx:197
 TFITS.cxx:198
 TFITS.cxx:199
 TFITS.cxx:200
 TFITS.cxx:201
 TFITS.cxx:202
 TFITS.cxx:203
 TFITS.cxx:204
 TFITS.cxx:205
 TFITS.cxx:206
 TFITS.cxx:207
 TFITS.cxx:208
 TFITS.cxx:209
 TFITS.cxx:210
 TFITS.cxx:211
 TFITS.cxx:212
 TFITS.cxx:213
 TFITS.cxx:214
 TFITS.cxx:215
 TFITS.cxx:216
 TFITS.cxx:217
 TFITS.cxx:218
 TFITS.cxx:219
 TFITS.cxx:220
 TFITS.cxx:221
 TFITS.cxx:222
 TFITS.cxx:223
 TFITS.cxx:224
 TFITS.cxx:225
 TFITS.cxx:226
 TFITS.cxx:227
 TFITS.cxx:228
 TFITS.cxx:229
 TFITS.cxx:230
 TFITS.cxx:231
 TFITS.cxx:232
 TFITS.cxx:233
 TFITS.cxx:234
 TFITS.cxx:235
 TFITS.cxx:236
 TFITS.cxx:237
 TFITS.cxx:238
 TFITS.cxx:239
 TFITS.cxx:240
 TFITS.cxx:241
 TFITS.cxx:242
 TFITS.cxx:243
 TFITS.cxx:244
 TFITS.cxx:245
 TFITS.cxx:246
 TFITS.cxx:247
 TFITS.cxx:248
 TFITS.cxx:249
 TFITS.cxx:250
 TFITS.cxx:251
 TFITS.cxx:252
 TFITS.cxx:253
 TFITS.cxx:254
 TFITS.cxx:255
 TFITS.cxx:256
 TFITS.cxx:257
 TFITS.cxx:258
 TFITS.cxx:259
 TFITS.cxx:260
 TFITS.cxx:261
 TFITS.cxx:262
 TFITS.cxx:263
 TFITS.cxx:264
 TFITS.cxx:265
 TFITS.cxx:266
 TFITS.cxx:267
 TFITS.cxx:268
 TFITS.cxx:269
 TFITS.cxx:270
 TFITS.cxx:271
 TFITS.cxx:272
 TFITS.cxx:273
 TFITS.cxx:274
 TFITS.cxx:275
 TFITS.cxx:276
 TFITS.cxx:277
 TFITS.cxx:278
 TFITS.cxx:279
 TFITS.cxx:280
 TFITS.cxx:281
 TFITS.cxx:282
 TFITS.cxx:283
 TFITS.cxx:284
 TFITS.cxx:285
 TFITS.cxx:286
 TFITS.cxx:287
 TFITS.cxx:288
 TFITS.cxx:289
 TFITS.cxx:290
 TFITS.cxx:291
 TFITS.cxx:292
 TFITS.cxx:293
 TFITS.cxx:294
 TFITS.cxx:295
 TFITS.cxx:296
 TFITS.cxx:297
 TFITS.cxx:298
 TFITS.cxx:299
 TFITS.cxx:300
 TFITS.cxx:301
 TFITS.cxx:302
 TFITS.cxx:303
 TFITS.cxx:304
 TFITS.cxx:305
 TFITS.cxx:306
 TFITS.cxx:307
 TFITS.cxx:308
 TFITS.cxx:309
 TFITS.cxx:310
 TFITS.cxx:311
 TFITS.cxx:312
 TFITS.cxx:313
 TFITS.cxx:314
 TFITS.cxx:315
 TFITS.cxx:316
 TFITS.cxx:317
 TFITS.cxx:318
 TFITS.cxx:319
 TFITS.cxx:320
 TFITS.cxx:321
 TFITS.cxx:322
 TFITS.cxx:323
 TFITS.cxx:324
 TFITS.cxx:325
 TFITS.cxx:326
 TFITS.cxx:327
 TFITS.cxx:328
 TFITS.cxx:329
 TFITS.cxx:330
 TFITS.cxx:331
 TFITS.cxx:332
 TFITS.cxx:333
 TFITS.cxx:334
 TFITS.cxx:335
 TFITS.cxx:336
 TFITS.cxx:337
 TFITS.cxx:338
 TFITS.cxx:339
 TFITS.cxx:340
 TFITS.cxx:341
 TFITS.cxx:342
 TFITS.cxx:343
 TFITS.cxx:344
 TFITS.cxx:345
 TFITS.cxx:346
 TFITS.cxx:347
 TFITS.cxx:348
 TFITS.cxx:349
 TFITS.cxx:350
 TFITS.cxx:351
 TFITS.cxx:352
 TFITS.cxx:353
 TFITS.cxx:354
 TFITS.cxx:355
 TFITS.cxx:356
 TFITS.cxx:357
 TFITS.cxx:358
 TFITS.cxx:359
 TFITS.cxx:360
 TFITS.cxx:361
 TFITS.cxx:362
 TFITS.cxx:363
 TFITS.cxx:364
 TFITS.cxx:365
 TFITS.cxx:366
 TFITS.cxx:367
 TFITS.cxx:368
 TFITS.cxx:369
 TFITS.cxx:370
 TFITS.cxx:371
 TFITS.cxx:372
 TFITS.cxx:373
 TFITS.cxx:374
 TFITS.cxx:375
 TFITS.cxx:376
 TFITS.cxx:377
 TFITS.cxx:378
 TFITS.cxx:379
 TFITS.cxx:380
 TFITS.cxx:381
 TFITS.cxx:382
 TFITS.cxx:383
 TFITS.cxx:384
 TFITS.cxx:385
 TFITS.cxx:386
 TFITS.cxx:387
 TFITS.cxx:388
 TFITS.cxx:389
 TFITS.cxx:390
 TFITS.cxx:391
 TFITS.cxx:392
 TFITS.cxx:393
 TFITS.cxx:394
 TFITS.cxx:395
 TFITS.cxx:396
 TFITS.cxx:397
 TFITS.cxx:398
 TFITS.cxx:399
 TFITS.cxx:400
 TFITS.cxx:401
 TFITS.cxx:402
 TFITS.cxx:403
 TFITS.cxx:404
 TFITS.cxx:405
 TFITS.cxx:406
 TFITS.cxx:407
 TFITS.cxx:408
 TFITS.cxx:409
 TFITS.cxx:410
 TFITS.cxx:411
 TFITS.cxx:412
 TFITS.cxx:413
 TFITS.cxx:414
 TFITS.cxx:415
 TFITS.cxx:416
 TFITS.cxx:417
 TFITS.cxx:418
 TFITS.cxx:419
 TFITS.cxx:420
 TFITS.cxx:421
 TFITS.cxx:422
 TFITS.cxx:423
 TFITS.cxx:424
 TFITS.cxx:425
 TFITS.cxx:426
 TFITS.cxx:427
 TFITS.cxx:428
 TFITS.cxx:429
 TFITS.cxx:430
 TFITS.cxx:431
 TFITS.cxx:432
 TFITS.cxx:433
 TFITS.cxx:434
 TFITS.cxx:435
 TFITS.cxx:436
 TFITS.cxx:437
 TFITS.cxx:438
 TFITS.cxx:439
 TFITS.cxx:440
 TFITS.cxx:441
 TFITS.cxx:442
 TFITS.cxx:443
 TFITS.cxx:444
 TFITS.cxx:445
 TFITS.cxx:446
 TFITS.cxx:447
 TFITS.cxx:448
 TFITS.cxx:449
 TFITS.cxx:450
 TFITS.cxx:451
 TFITS.cxx:452
 TFITS.cxx:453
 TFITS.cxx:454
 TFITS.cxx:455
 TFITS.cxx:456
 TFITS.cxx:457
 TFITS.cxx:458
 TFITS.cxx:459
 TFITS.cxx:460
 TFITS.cxx:461
 TFITS.cxx:462
 TFITS.cxx:463
 TFITS.cxx:464
 TFITS.cxx:465
 TFITS.cxx:466
 TFITS.cxx:467
 TFITS.cxx:468
 TFITS.cxx:469
 TFITS.cxx:470
 TFITS.cxx:471
 TFITS.cxx:472
 TFITS.cxx:473
 TFITS.cxx:474
 TFITS.cxx:475
 TFITS.cxx:476
 TFITS.cxx:477
 TFITS.cxx:478
 TFITS.cxx:479
 TFITS.cxx:480
 TFITS.cxx:481
 TFITS.cxx:482
 TFITS.cxx:483
 TFITS.cxx:484
 TFITS.cxx:485
 TFITS.cxx:486
 TFITS.cxx:487
 TFITS.cxx:488
 TFITS.cxx:489
 TFITS.cxx:490
 TFITS.cxx:491
 TFITS.cxx:492
 TFITS.cxx:493
 TFITS.cxx:494
 TFITS.cxx:495
 TFITS.cxx:496
 TFITS.cxx:497
 TFITS.cxx:498
 TFITS.cxx:499
 TFITS.cxx:500
 TFITS.cxx:501
 TFITS.cxx:502
 TFITS.cxx:503
 TFITS.cxx:504
 TFITS.cxx:505
 TFITS.cxx:506
 TFITS.cxx:507
 TFITS.cxx:508
 TFITS.cxx:509
 TFITS.cxx:510
 TFITS.cxx:511
 TFITS.cxx:512
 TFITS.cxx:513
 TFITS.cxx:514
 TFITS.cxx:515
 TFITS.cxx:516
 TFITS.cxx:517
 TFITS.cxx:518
 TFITS.cxx:519
 TFITS.cxx:520
 TFITS.cxx:521
 TFITS.cxx:522
 TFITS.cxx:523
 TFITS.cxx:524
 TFITS.cxx:525
 TFITS.cxx:526
 TFITS.cxx:527
 TFITS.cxx:528
 TFITS.cxx:529
 TFITS.cxx:530
 TFITS.cxx:531
 TFITS.cxx:532
 TFITS.cxx:533
 TFITS.cxx:534
 TFITS.cxx:535
 TFITS.cxx:536
 TFITS.cxx:537
 TFITS.cxx:538
 TFITS.cxx:539
 TFITS.cxx:540
 TFITS.cxx:541
 TFITS.cxx:542
 TFITS.cxx:543
 TFITS.cxx:544
 TFITS.cxx:545
 TFITS.cxx:546
 TFITS.cxx:547
 TFITS.cxx:548
 TFITS.cxx:549
 TFITS.cxx:550
 TFITS.cxx:551
 TFITS.cxx:552
 TFITS.cxx:553
 TFITS.cxx:554
 TFITS.cxx:555
 TFITS.cxx:556
 TFITS.cxx:557
 TFITS.cxx:558
 TFITS.cxx:559
 TFITS.cxx:560
 TFITS.cxx:561
 TFITS.cxx:562
 TFITS.cxx:563
 TFITS.cxx:564
 TFITS.cxx:565
 TFITS.cxx:566
 TFITS.cxx:567
 TFITS.cxx:568
 TFITS.cxx:569
 TFITS.cxx:570
 TFITS.cxx:571
 TFITS.cxx:572
 TFITS.cxx:573
 TFITS.cxx:574
 TFITS.cxx:575
 TFITS.cxx:576
 TFITS.cxx:577
 TFITS.cxx:578
 TFITS.cxx:579
 TFITS.cxx:580
 TFITS.cxx:581
 TFITS.cxx:582
 TFITS.cxx:583
 TFITS.cxx:584
 TFITS.cxx:585
 TFITS.cxx:586
 TFITS.cxx:587
 TFITS.cxx:588
 TFITS.cxx:589
 TFITS.cxx:590
 TFITS.cxx:591
 TFITS.cxx:592
 TFITS.cxx:593
 TFITS.cxx:594
 TFITS.cxx:595
 TFITS.cxx:596
 TFITS.cxx:597
 TFITS.cxx:598
 TFITS.cxx:599
 TFITS.cxx:600
 TFITS.cxx:601
 TFITS.cxx:602
 TFITS.cxx:603
 TFITS.cxx:604
 TFITS.cxx:605
 TFITS.cxx:606
 TFITS.cxx:607
 TFITS.cxx:608
 TFITS.cxx:609
 TFITS.cxx:610
 TFITS.cxx:611
 TFITS.cxx:612
 TFITS.cxx:613
 TFITS.cxx:614
 TFITS.cxx:615
 TFITS.cxx:616
 TFITS.cxx:617
 TFITS.cxx:618
 TFITS.cxx:619
 TFITS.cxx:620
 TFITS.cxx:621
 TFITS.cxx:622
 TFITS.cxx:623
 TFITS.cxx:624
 TFITS.cxx:625
 TFITS.cxx:626
 TFITS.cxx:627
 TFITS.cxx:628
 TFITS.cxx:629
 TFITS.cxx:630
 TFITS.cxx:631
 TFITS.cxx:632
 TFITS.cxx:633
 TFITS.cxx:634
 TFITS.cxx:635
 TFITS.cxx:636
 TFITS.cxx:637
 TFITS.cxx:638
 TFITS.cxx:639
 TFITS.cxx:640
 TFITS.cxx:641
 TFITS.cxx:642
 TFITS.cxx:643
 TFITS.cxx:644
 TFITS.cxx:645
 TFITS.cxx:646
 TFITS.cxx:647
 TFITS.cxx:648
 TFITS.cxx:649
 TFITS.cxx:650
 TFITS.cxx:651
 TFITS.cxx:652
 TFITS.cxx:653
 TFITS.cxx:654
 TFITS.cxx:655
 TFITS.cxx:656
 TFITS.cxx:657
 TFITS.cxx:658
 TFITS.cxx:659
 TFITS.cxx:660
 TFITS.cxx:661
 TFITS.cxx:662
 TFITS.cxx:663
 TFITS.cxx:664
 TFITS.cxx:665
 TFITS.cxx:666
 TFITS.cxx:667
 TFITS.cxx:668
 TFITS.cxx:669
 TFITS.cxx:670
 TFITS.cxx:671
 TFITS.cxx:672
 TFITS.cxx:673
 TFITS.cxx:674
 TFITS.cxx:675
 TFITS.cxx:676
 TFITS.cxx:677
 TFITS.cxx:678
 TFITS.cxx:679
 TFITS.cxx:680
 TFITS.cxx:681
 TFITS.cxx:682
 TFITS.cxx:683
 TFITS.cxx:684
 TFITS.cxx:685
 TFITS.cxx:686
 TFITS.cxx:687
 TFITS.cxx:688
 TFITS.cxx:689
 TFITS.cxx:690
 TFITS.cxx:691
 TFITS.cxx:692
 TFITS.cxx:693
 TFITS.cxx:694
 TFITS.cxx:695
 TFITS.cxx:696
 TFITS.cxx:697
 TFITS.cxx:698
 TFITS.cxx:699
 TFITS.cxx:700
 TFITS.cxx:701
 TFITS.cxx:702
 TFITS.cxx:703
 TFITS.cxx:704
 TFITS.cxx:705
 TFITS.cxx:706
 TFITS.cxx:707
 TFITS.cxx:708
 TFITS.cxx:709
 TFITS.cxx:710
 TFITS.cxx:711
 TFITS.cxx:712
 TFITS.cxx:713
 TFITS.cxx:714
 TFITS.cxx:715
 TFITS.cxx:716
 TFITS.cxx:717
 TFITS.cxx:718
 TFITS.cxx:719
 TFITS.cxx:720
 TFITS.cxx:721
 TFITS.cxx:722
 TFITS.cxx:723
 TFITS.cxx:724
 TFITS.cxx:725
 TFITS.cxx:726
 TFITS.cxx:727
 TFITS.cxx:728
 TFITS.cxx:729
 TFITS.cxx:730
 TFITS.cxx:731
 TFITS.cxx:732
 TFITS.cxx:733
 TFITS.cxx:734
 TFITS.cxx:735
 TFITS.cxx:736
 TFITS.cxx:737
 TFITS.cxx:738
 TFITS.cxx:739
 TFITS.cxx:740
 TFITS.cxx:741
 TFITS.cxx:742
 TFITS.cxx:743
 TFITS.cxx:744
 TFITS.cxx:745
 TFITS.cxx:746
 TFITS.cxx:747
 TFITS.cxx:748
 TFITS.cxx:749
 TFITS.cxx:750
 TFITS.cxx:751
 TFITS.cxx:752
 TFITS.cxx:753
 TFITS.cxx:754
 TFITS.cxx:755
 TFITS.cxx:756
 TFITS.cxx:757
 TFITS.cxx:758
 TFITS.cxx:759
 TFITS.cxx:760
 TFITS.cxx:761
 TFITS.cxx:762
 TFITS.cxx:763
 TFITS.cxx:764
 TFITS.cxx:765
 TFITS.cxx:766
 TFITS.cxx:767
 TFITS.cxx:768
 TFITS.cxx:769
 TFITS.cxx:770
 TFITS.cxx:771
 TFITS.cxx:772
 TFITS.cxx:773
 TFITS.cxx:774
 TFITS.cxx:775
 TFITS.cxx:776
 TFITS.cxx:777
 TFITS.cxx:778
 TFITS.cxx:779
 TFITS.cxx:780
 TFITS.cxx:781
 TFITS.cxx:782
 TFITS.cxx:783
 TFITS.cxx:784
 TFITS.cxx:785
 TFITS.cxx:786
 TFITS.cxx:787
 TFITS.cxx:788
 TFITS.cxx:789
 TFITS.cxx:790
 TFITS.cxx:791
 TFITS.cxx:792
 TFITS.cxx:793
 TFITS.cxx:794
 TFITS.cxx:795
 TFITS.cxx:796
 TFITS.cxx:797
 TFITS.cxx:798
 TFITS.cxx:799
 TFITS.cxx:800
 TFITS.cxx:801
 TFITS.cxx:802
 TFITS.cxx:803
 TFITS.cxx:804
 TFITS.cxx:805
 TFITS.cxx:806
 TFITS.cxx:807
 TFITS.cxx:808
 TFITS.cxx:809
 TFITS.cxx:810
 TFITS.cxx:811
 TFITS.cxx:812
 TFITS.cxx:813
 TFITS.cxx:814
 TFITS.cxx:815
 TFITS.cxx:816
 TFITS.cxx:817
 TFITS.cxx:818
 TFITS.cxx:819
 TFITS.cxx:820
 TFITS.cxx:821
 TFITS.cxx:822
 TFITS.cxx:823
 TFITS.cxx:824
 TFITS.cxx:825
 TFITS.cxx:826
 TFITS.cxx:827
 TFITS.cxx:828
 TFITS.cxx:829
 TFITS.cxx:830
 TFITS.cxx:831
 TFITS.cxx:832
 TFITS.cxx:833
 TFITS.cxx:834
 TFITS.cxx:835
 TFITS.cxx:836
 TFITS.cxx:837
 TFITS.cxx:838
 TFITS.cxx:839
 TFITS.cxx:840
 TFITS.cxx:841
 TFITS.cxx:842
 TFITS.cxx:843
 TFITS.cxx:844
 TFITS.cxx:845
 TFITS.cxx:846
 TFITS.cxx:847
 TFITS.cxx:848
 TFITS.cxx:849
 TFITS.cxx:850
 TFITS.cxx:851
 TFITS.cxx:852
 TFITS.cxx:853
 TFITS.cxx:854
 TFITS.cxx:855
 TFITS.cxx:856
 TFITS.cxx:857
 TFITS.cxx:858
 TFITS.cxx:859
 TFITS.cxx:860
 TFITS.cxx:861
 TFITS.cxx:862
 TFITS.cxx:863
 TFITS.cxx:864
 TFITS.cxx:865
 TFITS.cxx:866
 TFITS.cxx:867
 TFITS.cxx:868
 TFITS.cxx:869
 TFITS.cxx:870
 TFITS.cxx:871
 TFITS.cxx:872
 TFITS.cxx:873
 TFITS.cxx:874
 TFITS.cxx:875
 TFITS.cxx:876
 TFITS.cxx:877
 TFITS.cxx:878
 TFITS.cxx:879
 TFITS.cxx:880
 TFITS.cxx:881
 TFITS.cxx:882
 TFITS.cxx:883
 TFITS.cxx:884
 TFITS.cxx:885
 TFITS.cxx:886
 TFITS.cxx:887
 TFITS.cxx:888
 TFITS.cxx:889
 TFITS.cxx:890
 TFITS.cxx:891
 TFITS.cxx:892
 TFITS.cxx:893
 TFITS.cxx:894
 TFITS.cxx:895
 TFITS.cxx:896
 TFITS.cxx:897
 TFITS.cxx:898
 TFITS.cxx:899
 TFITS.cxx:900
 TFITS.cxx:901
 TFITS.cxx:902
 TFITS.cxx:903
 TFITS.cxx:904
 TFITS.cxx:905
 TFITS.cxx:906
 TFITS.cxx:907
 TFITS.cxx:908
 TFITS.cxx:909
 TFITS.cxx:910
 TFITS.cxx:911
 TFITS.cxx:912
 TFITS.cxx:913
 TFITS.cxx:914
 TFITS.cxx:915
 TFITS.cxx:916
 TFITS.cxx:917
 TFITS.cxx:918
 TFITS.cxx:919
 TFITS.cxx:920
 TFITS.cxx:921
 TFITS.cxx:922
 TFITS.cxx:923
 TFITS.cxx:924
 TFITS.cxx:925
 TFITS.cxx:926
 TFITS.cxx:927
 TFITS.cxx:928
 TFITS.cxx:929
 TFITS.cxx:930
 TFITS.cxx:931
 TFITS.cxx:932
 TFITS.cxx:933
 TFITS.cxx:934
 TFITS.cxx:935
 TFITS.cxx:936
 TFITS.cxx:937
 TFITS.cxx:938
 TFITS.cxx:939
 TFITS.cxx:940
 TFITS.cxx:941
 TFITS.cxx:942
 TFITS.cxx:943
 TFITS.cxx:944
 TFITS.cxx:945
 TFITS.cxx:946
 TFITS.cxx:947
 TFITS.cxx:948
 TFITS.cxx:949
 TFITS.cxx:950
 TFITS.cxx:951
 TFITS.cxx:952
 TFITS.cxx:953
 TFITS.cxx:954
 TFITS.cxx:955
 TFITS.cxx:956
 TFITS.cxx:957
 TFITS.cxx:958
 TFITS.cxx:959
 TFITS.cxx:960
 TFITS.cxx:961
 TFITS.cxx:962
 TFITS.cxx:963
 TFITS.cxx:964
 TFITS.cxx:965
 TFITS.cxx:966
 TFITS.cxx:967
 TFITS.cxx:968
 TFITS.cxx:969
 TFITS.cxx:970
 TFITS.cxx:971
 TFITS.cxx:972
 TFITS.cxx:973
 TFITS.cxx:974
 TFITS.cxx:975
 TFITS.cxx:976
 TFITS.cxx:977
 TFITS.cxx:978
 TFITS.cxx:979
 TFITS.cxx:980
 TFITS.cxx:981
 TFITS.cxx:982
 TFITS.cxx:983
 TFITS.cxx:984
 TFITS.cxx:985
 TFITS.cxx:986
 TFITS.cxx:987
 TFITS.cxx:988
 TFITS.cxx:989
 TFITS.cxx:990
 TFITS.cxx:991
 TFITS.cxx:992
 TFITS.cxx:993
 TFITS.cxx:994
 TFITS.cxx:995
 TFITS.cxx:996
 TFITS.cxx:997
 TFITS.cxx:998
 TFITS.cxx:999
 TFITS.cxx:1000
 TFITS.cxx:1001
 TFITS.cxx:1002
 TFITS.cxx:1003
 TFITS.cxx:1004
 TFITS.cxx:1005
 TFITS.cxx:1006
 TFITS.cxx:1007
 TFITS.cxx:1008
 TFITS.cxx:1009
 TFITS.cxx:1010
 TFITS.cxx:1011
 TFITS.cxx:1012
 TFITS.cxx:1013
 TFITS.cxx:1014
 TFITS.cxx:1015
 TFITS.cxx:1016
 TFITS.cxx:1017
 TFITS.cxx:1018
 TFITS.cxx:1019
 TFITS.cxx:1020
 TFITS.cxx:1021
 TFITS.cxx:1022
 TFITS.cxx:1023
 TFITS.cxx:1024
 TFITS.cxx:1025
 TFITS.cxx:1026
 TFITS.cxx:1027
 TFITS.cxx:1028
 TFITS.cxx:1029
 TFITS.cxx:1030
 TFITS.cxx:1031
 TFITS.cxx:1032
 TFITS.cxx:1033
 TFITS.cxx:1034
 TFITS.cxx:1035
 TFITS.cxx:1036
 TFITS.cxx:1037
 TFITS.cxx:1038
 TFITS.cxx:1039
 TFITS.cxx:1040
 TFITS.cxx:1041
 TFITS.cxx:1042
 TFITS.cxx:1043
 TFITS.cxx:1044
 TFITS.cxx:1045
 TFITS.cxx:1046
 TFITS.cxx:1047
 TFITS.cxx:1048
 TFITS.cxx:1049
 TFITS.cxx:1050
 TFITS.cxx:1051
 TFITS.cxx:1052
 TFITS.cxx:1053
 TFITS.cxx:1054
 TFITS.cxx:1055
 TFITS.cxx:1056
 TFITS.cxx:1057
 TFITS.cxx:1058
 TFITS.cxx:1059
 TFITS.cxx:1060
 TFITS.cxx:1061
 TFITS.cxx:1062
 TFITS.cxx:1063
 TFITS.cxx:1064
 TFITS.cxx:1065
 TFITS.cxx:1066
 TFITS.cxx:1067
 TFITS.cxx:1068
 TFITS.cxx:1069
 TFITS.cxx:1070
 TFITS.cxx:1071
 TFITS.cxx:1072
 TFITS.cxx:1073
 TFITS.cxx:1074
 TFITS.cxx:1075
 TFITS.cxx:1076
 TFITS.cxx:1077
 TFITS.cxx:1078
 TFITS.cxx:1079
 TFITS.cxx:1080
 TFITS.cxx:1081
 TFITS.cxx:1082
 TFITS.cxx:1083
 TFITS.cxx:1084
 TFITS.cxx:1085
 TFITS.cxx:1086
 TFITS.cxx:1087
 TFITS.cxx:1088
 TFITS.cxx:1089
 TFITS.cxx:1090
 TFITS.cxx:1091
 TFITS.cxx:1092
 TFITS.cxx:1093
 TFITS.cxx:1094
 TFITS.cxx:1095
 TFITS.cxx:1096
 TFITS.cxx:1097
 TFITS.cxx:1098
 TFITS.cxx:1099
 TFITS.cxx:1100
 TFITS.cxx:1101
 TFITS.cxx:1102
 TFITS.cxx:1103
 TFITS.cxx:1104
 TFITS.cxx:1105
 TFITS.cxx:1106
 TFITS.cxx:1107
 TFITS.cxx:1108
 TFITS.cxx:1109
 TFITS.cxx:1110
 TFITS.cxx:1111
 TFITS.cxx:1112
 TFITS.cxx:1113
 TFITS.cxx:1114
 TFITS.cxx:1115
 TFITS.cxx:1116
 TFITS.cxx:1117
 TFITS.cxx:1118
 TFITS.cxx:1119
 TFITS.cxx:1120
 TFITS.cxx:1121
 TFITS.cxx:1122
 TFITS.cxx:1123
 TFITS.cxx:1124
 TFITS.cxx:1125
 TFITS.cxx:1126
 TFITS.cxx:1127
 TFITS.cxx:1128
 TFITS.cxx:1129
 TFITS.cxx:1130
 TFITS.cxx:1131
 TFITS.cxx:1132
 TFITS.cxx:1133
 TFITS.cxx:1134
 TFITS.cxx:1135
 TFITS.cxx:1136
 TFITS.cxx:1137
 TFITS.cxx:1138
 TFITS.cxx:1139
 TFITS.cxx:1140
 TFITS.cxx:1141
 TFITS.cxx:1142
 TFITS.cxx:1143
 TFITS.cxx:1144
 TFITS.cxx:1145
 TFITS.cxx:1146
 TFITS.cxx:1147
 TFITS.cxx:1148
 TFITS.cxx:1149
 TFITS.cxx:1150
 TFITS.cxx:1151
 TFITS.cxx:1152
 TFITS.cxx:1153
 TFITS.cxx:1154
 TFITS.cxx:1155
 TFITS.cxx:1156
 TFITS.cxx:1157
 TFITS.cxx:1158
 TFITS.cxx:1159
 TFITS.cxx:1160
 TFITS.cxx:1161
 TFITS.cxx:1162
 TFITS.cxx:1163
 TFITS.cxx:1164
 TFITS.cxx:1165
 TFITS.cxx:1166
 TFITS.cxx:1167
 TFITS.cxx:1168
 TFITS.cxx:1169
 TFITS.cxx:1170
 TFITS.cxx:1171
 TFITS.cxx:1172
 TFITS.cxx:1173
 TFITS.cxx:1174
 TFITS.cxx:1175
 TFITS.cxx:1176
 TFITS.cxx:1177
 TFITS.cxx:1178
 TFITS.cxx:1179
 TFITS.cxx:1180
 TFITS.cxx:1181
 TFITS.cxx:1182
 TFITS.cxx:1183
 TFITS.cxx:1184
 TFITS.cxx:1185
 TFITS.cxx:1186
 TFITS.cxx:1187
 TFITS.cxx:1188
 TFITS.cxx:1189
 TFITS.cxx:1190
 TFITS.cxx:1191
 TFITS.cxx:1192
 TFITS.cxx:1193
 TFITS.cxx:1194
 TFITS.cxx:1195
 TFITS.cxx:1196
 TFITS.cxx:1197
 TFITS.cxx:1198
 TFITS.cxx:1199
 TFITS.cxx:1200
 TFITS.cxx:1201
 TFITS.cxx:1202
 TFITS.cxx:1203
 TFITS.cxx:1204
 TFITS.cxx:1205
 TFITS.cxx:1206
 TFITS.cxx:1207
 TFITS.cxx:1208
 TFITS.cxx:1209
 TFITS.cxx:1210
 TFITS.cxx:1211
 TFITS.cxx:1212
 TFITS.cxx:1213
 TFITS.cxx:1214
 TFITS.cxx:1215
 TFITS.cxx:1216
 TFITS.cxx:1217
 TFITS.cxx:1218
 TFITS.cxx:1219
 TFITS.cxx:1220
 TFITS.cxx:1221
 TFITS.cxx:1222
 TFITS.cxx:1223
 TFITS.cxx:1224
 TFITS.cxx:1225
 TFITS.cxx:1226
 TFITS.cxx:1227
 TFITS.cxx:1228
 TFITS.cxx:1229
 TFITS.cxx:1230
 TFITS.cxx:1231
 TFITS.cxx:1232
 TFITS.cxx:1233
 TFITS.cxx:1234
 TFITS.cxx:1235
 TFITS.cxx:1236
 TFITS.cxx:1237
 TFITS.cxx:1238
 TFITS.cxx:1239
 TFITS.cxx:1240
 TFITS.cxx:1241
 TFITS.cxx:1242
 TFITS.cxx:1243
 TFITS.cxx:1244
 TFITS.cxx:1245
 TFITS.cxx:1246
 TFITS.cxx:1247
 TFITS.cxx:1248
 TFITS.cxx:1249
 TFITS.cxx:1250
 TFITS.cxx:1251
 TFITS.cxx:1252
 TFITS.cxx:1253
 TFITS.cxx:1254
 TFITS.cxx:1255
 TFITS.cxx:1256
 TFITS.cxx:1257
 TFITS.cxx:1258
 TFITS.cxx:1259
 TFITS.cxx:1260
 TFITS.cxx:1261
 TFITS.cxx:1262
 TFITS.cxx:1263
 TFITS.cxx:1264
 TFITS.cxx:1265
 TFITS.cxx:1266
 TFITS.cxx:1267
 TFITS.cxx:1268
 TFITS.cxx:1269
 TFITS.cxx:1270
 TFITS.cxx:1271
 TFITS.cxx:1272
 TFITS.cxx:1273
 TFITS.cxx:1274
 TFITS.cxx:1275
 TFITS.cxx:1276
 TFITS.cxx:1277
 TFITS.cxx:1278
 TFITS.cxx:1279
 TFITS.cxx:1280
 TFITS.cxx:1281
 TFITS.cxx:1282
 TFITS.cxx:1283
 TFITS.cxx:1284
 TFITS.cxx:1285
 TFITS.cxx:1286
 TFITS.cxx:1287
 TFITS.cxx:1288
 TFITS.cxx:1289
 TFITS.cxx:1290
 TFITS.cxx:1291
 TFITS.cxx:1292
 TFITS.cxx:1293
 TFITS.cxx:1294
 TFITS.cxx:1295
 TFITS.cxx:1296
 TFITS.cxx:1297
 TFITS.cxx:1298
 TFITS.cxx:1299
 TFITS.cxx:1300
 TFITS.cxx:1301
 TFITS.cxx:1302
 TFITS.cxx:1303
 TFITS.cxx:1304
 TFITS.cxx:1305
 TFITS.cxx:1306
 TFITS.cxx:1307
 TFITS.cxx:1308
 TFITS.cxx:1309
 TFITS.cxx:1310
 TFITS.cxx:1311
 TFITS.cxx:1312
 TFITS.cxx:1313
 TFITS.cxx:1314
 TFITS.cxx:1315
 TFITS.cxx:1316
 TFITS.cxx:1317
 TFITS.cxx:1318
 TFITS.cxx:1319
 TFITS.cxx:1320
 TFITS.cxx:1321
 TFITS.cxx:1322
 TFITS.cxx:1323
 TFITS.cxx:1324
 TFITS.cxx:1325
 TFITS.cxx:1326
 TFITS.cxx:1327
 TFITS.cxx:1328
 TFITS.cxx:1329
 TFITS.cxx:1330
 TFITS.cxx:1331
 TFITS.cxx:1332
 TFITS.cxx:1333
 TFITS.cxx:1334
 TFITS.cxx:1335
 TFITS.cxx:1336
 TFITS.cxx:1337
 TFITS.cxx:1338
 TFITS.cxx:1339
 TFITS.cxx:1340
 TFITS.cxx:1341
 TFITS.cxx:1342
 TFITS.cxx:1343
 TFITS.cxx:1344
 TFITS.cxx:1345
 TFITS.cxx:1346
 TFITS.cxx:1347
 TFITS.cxx:1348
 TFITS.cxx:1349
 TFITS.cxx:1350
 TFITS.cxx:1351
 TFITS.cxx:1352
 TFITS.cxx:1353
 TFITS.cxx:1354
 TFITS.cxx:1355
 TFITS.cxx:1356
 TFITS.cxx:1357
 TFITS.cxx:1358
 TFITS.cxx:1359
 TFITS.cxx:1360
 TFITS.cxx:1361
 TFITS.cxx:1362
 TFITS.cxx:1363
 TFITS.cxx:1364
 TFITS.cxx:1365
 TFITS.cxx:1366
 TFITS.cxx:1367
 TFITS.cxx:1368
 TFITS.cxx:1369
 TFITS.cxx:1370
 TFITS.cxx:1371
 TFITS.cxx:1372
 TFITS.cxx:1373
 TFITS.cxx:1374
 TFITS.cxx:1375
 TFITS.cxx:1376
 TFITS.cxx:1377
 TFITS.cxx:1378
 TFITS.cxx:1379
 TFITS.cxx:1380
 TFITS.cxx:1381
 TFITS.cxx:1382
 TFITS.cxx:1383
 TFITS.cxx:1384
 TFITS.cxx:1385
 TFITS.cxx:1386
 TFITS.cxx:1387
 TFITS.cxx:1388
 TFITS.cxx:1389
 TFITS.cxx:1390
 TFITS.cxx:1391
 TFITS.cxx:1392
 TFITS.cxx:1393
 TFITS.cxx:1394
 TFITS.cxx:1395
 TFITS.cxx:1396
 TFITS.cxx:1397
 TFITS.cxx:1398
 TFITS.cxx:1399
 TFITS.cxx:1400
 TFITS.cxx:1401
 TFITS.cxx:1402
 TFITS.cxx:1403
 TFITS.cxx:1404
 TFITS.cxx:1405
 TFITS.cxx:1406
 TFITS.cxx:1407
 TFITS.cxx:1408
 TFITS.cxx:1409
 TFITS.cxx:1410
 TFITS.cxx:1411
 TFITS.cxx:1412
 TFITS.cxx:1413
 TFITS.cxx:1414
 TFITS.cxx:1415
 TFITS.cxx:1416
 TFITS.cxx:1417
 TFITS.cxx:1418
 TFITS.cxx:1419
 TFITS.cxx:1420
 TFITS.cxx:1421
 TFITS.cxx:1422
 TFITS.cxx:1423
 TFITS.cxx:1424
 TFITS.cxx:1425
 TFITS.cxx:1426
 TFITS.cxx:1427
 TFITS.cxx:1428
 TFITS.cxx:1429
 TFITS.cxx:1430
 TFITS.cxx:1431