source: trunk/src/FitsIO/subsection.cc @ 326

// -----------------------------------------------------------------------
// subsection.cc: Make sure subsections are valid for a particular
//                FITS file.
// -----------------------------------------------------------------------
// Copyright (C) 2006, Matthew Whiting, ATNF
//
// This program is free software; you can redistribute it and/or modify it
// under the terms of the GNU General Public License as published by the
// Free Software Foundation; either version 2 of the License, or (at your
// option) any later version.
//
// Duchamp is distributed in the hope that it will be useful, but WITHOUT
// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
// FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
// for more details.
//
// You should have received a copy of the GNU General Public License
// along with Duchamp; if not, write to the Free Software Foundation,
// Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307, USA
//
// Correspondence concerning Duchamp may be directed to:
//    Internet email: Matthew.Whiting [at] atnf.csiro.au
//    Postal address: Dr. Matthew Whiting
//                    Australia Telescope National Facility, CSIRO
//                    PO Box 76
//                    Epping NSW 1710
//                    AUSTRALIA
// -----------------------------------------------------------------------
#include <iostream>
#include <sstream>
#include <string>
#include <vector>
#include <math.h>
#include <wcs.h>
#define WCSLIB_GETWCSTAB // define this so that we don't try and redefine 
                         //  wtbarr (this is a problem when using gcc v.4+
#include <fitsio.h>
#include <param.hh>
#include <duchamp.hh>
#include <Utils/Section.hh>

void Param::setOffsets(struct wcsprm *wcs)
{
  /** 
   * If there is a subsection being used, set the offset values
   * according to the correct dimensions given by the WCS struct.  
   * If not, set the offsets to zero.
   * \param wcs The WCSLIB wcsprm struct that defines which axis is which.
   */
  if(this->flagSubsection){ // if so, then the offsets array is defined.
    this->xSubOffset = this->pixelSec.getStart(wcs->lng);
    this->ySubOffset = this->pixelSec.getStart(wcs->lat);
    this->zSubOffset = this->pixelSec.getStart(wcs->spec);
  }
  else{// else they should be 0
    this->xSubOffset = this->ySubOffset = this->zSubOffset = 0;
  }
}

void Param::setOffsets(struct wcsprm &wcs)
{
  /** 
   * If there is a subsection being used, set the offset values
   * according to the correct dimensions given by the WCS struct.  
   * If not, set the offsets to zero.
   * \param wcs The WCSLIB wcsprm struct that defines which axis is which.
   */
  if(this->flagSubsection){ // if so, then the offsets array is defined.
    this->xSubOffset = this->pixelSec.getStart(wcs.lng);
    this->ySubOffset = this->pixelSec.getStart(wcs.lat);
    this->zSubOffset = this->pixelSec.getStart(wcs.spec);
  }
  else{// else they should be 0
    this->xSubOffset = this->ySubOffset = this->zSubOffset = 0;
  }
}

int Param::verifySubsection()
{
  /**
   * Checks that the subsection strings (the pixel and stats
   * subsections) are in the appropriate format, with the correct
   * number of entries (one for each axis).
   *
   * This reads the dimensional information from the FITS file, and
   * uses this with the Section::parse() function to make sure each
   * section is OK.
   *
   * \return SUCCESS/FAILURE depending on outcome of the
   * Section::parse() calls. Also FAILURE if something goes wrong with
   * the FITS access.
   */

  if(!this->flagSubsection && !this->flagStatSec){
    // if we get here, we are using neither subsection
    return SUCCESS;
  }
  else{
    // otherwise, at least one of the subsections is being used, and
    // so we need to check them

    // First open the requested FITS file and check its existence and 
    //  number of axes.
    int numAxes,status = 0;  /* MUST initialize status */
    fitsfile *fptr;         

    // Make sure the FITS file exists
    int exists;
    fits_file_exists(this->imageFile.c_str(),&exists,&status);
    if(exists<=0){
      fits_report_error(stderr, status);
      duchampWarning("Cube Reader", "Requested image does not exist!\n");
      return FAILURE;
    }
    // Open the FITS file
    if( fits_open_file(&fptr,this->imageFile.c_str(),READONLY,&status) ){
      fits_report_error(stderr, status);
      return FAILURE;
    }
    // Read the size information -- number of axes and their sizes
    status = 0;
    if(fits_get_img_dim(fptr, &numAxes, &status)){
      fits_report_error(stderr, status);
      return FAILURE;
    }
    long *dimAxes = new long[numAxes];
    for(int i=0;i<numAxes;i++) dimAxes[i]=1;
    status = 0;
    if(fits_get_img_size(fptr, numAxes, dimAxes, &status)){
      fits_report_error(stderr, status);
      return FAILURE;
    }
    // Close the FITS file.
    status = 0;
    fits_close_file(fptr, &status);
    if (status){
      duchampWarning("Cube Reader","Error closing file: ");
      fits_report_error(stderr, status);
    }

    ///////////////////
    // Now parse the subsection and make sure all that works.
  
    std::vector<long> dim(numAxes);
    for(int i=0;i<numAxes;i++) dim[i] = dimAxes[i];
    delete [] dimAxes;
  
    if(this->flagSubsection)
      if(this->pixelSec.parse(dim)==FAILURE) return FAILURE;

    if(this->flagStatSec)
      if(this->statSec.parse(dim)==FAILURE)  return FAILURE;
  
    return SUCCESS;
  }

}