fits.h

Go to the documentation of this file.

// -*- lsst-c++ -*-
#ifndef LSST_AFW_fits_h_INCLUDED
#define LSST_AFW_fits_h_INCLUDED

#include <string>

#include <boost/format.hpp>

#include "lsst/base.h"
#include "lsst/pex/exceptions.h"
#include "lsst/daf/base.h"
#include "ndarray.h"

namespace lsst { namespace afw { namespace fits {

LSST_EXCEPTION_TYPE(FitsError, lsst::pex::exceptions::IoError, lsst::afw::fits::FitsError)


LSST_EXCEPTION_TYPE(FitsTypeError, lsst::afw::fits::FitsError, lsst::afw::fits::FitsTypeError)

#ifndef SWIG // only want SWIG to see the exceptions; everything else is too low-level for Python.

class HeaderIterationFunctor {
public:

    virtual void operator()(
        std::string const & key,
        std::string const & value,
        std::string const & comment
    ) = 0;

    virtual ~HeaderIterationFunctor() {}

};

std::string makeErrorMessage(std::string const & fileName="", int status=0, std::string const & msg="");
inline std::string makeErrorMessage(std::string const & fileName, int status, boost::format const & msg) {
    return makeErrorMessage(fileName, status, msg.str());
}

std::string makeErrorMessage(void * fptr, int status=0, std::string const & msg="");
inline std::string makeErrorMessage(void * fptr, int status, boost::format const & msg) {
    return makeErrorMessage(fptr, status, msg.str());
}

#define LSST_FITS_EXCEPT(type, fitsObj, ...) \
    type(LSST_EXCEPT_HERE, lsst::afw::fits::makeErrorMessage((fitsObj).fptr, (fitsObj).status, __VA_ARGS__))

#define LSST_FITS_CHECK_STATUS(fitsObj, ...)                            \
    if ((fitsObj).status != 0) throw LSST_FITS_EXCEPT(lsst::afw::fits::FitsError, fitsObj, __VA_ARGS__)

template <typename T> int getBitPix();

class MemFileManager {
public:

    MemFileManager() : _ptr(0), _len(0), _managed(true) {}

    explicit MemFileManager(std::size_t len) : _ptr(0), _len(0), _managed(true) { reset(len); }

    MemFileManager(void * ptr, std::size_t len) : _ptr(ptr), _len(len), _managed(false) {}

    void reset();

    void reset(std::size_t len);

    void reset(void * ptr, std::size_t len) { reset(); _ptr = ptr; _len = len; _managed = false; }

    ~MemFileManager() { reset(); }

    // No copying
    MemFileManager (const MemFileManager&) = delete;
    MemFileManager& operator=(const MemFileManager&) = delete;

    // No moving
    MemFileManager (MemFileManager&&) = delete;
    MemFileManager& operator=(MemFileManager&&) = delete;

    void* getData() const { return _ptr; }

    std::size_t getLength() const { return _len; }

private:

    friend class Fits;

    void * _ptr;
    std::size_t _len;
    bool _managed;
};

class Fits {
    template <typename T> void createImageImpl(int nAxis, long * nAxes);
    template <typename T> void writeImageImpl(T const * data, int nElements);
    template <typename T> void readImageImpl(int nAxis, T * data, long * begin, long * end, long * increment);
    void getImageShapeImpl(int maxDim, long * nAxes);
public:

    enum BehaviorFlags {
        AUTO_CLOSE = 0x01, // Close files when the Fits object goes out of scope if fptr != NULL
        AUTO_CHECK = 0x02  // Call LSST_FITS_CHECK_STATUS after every cfitsio call
    };

    std::string getFileName() const;

    int getHdu();

    void setHdu(int hdu, bool relative=false);

    int countHdus();

    template <typename T>
    void updateKey(std::string const & key, T const & value, std::string const & comment);
    void updateKey(std::string const & key, char const * value, std::string const & comment) {
        updateKey(key, std::string(value), comment);
    }
    template <typename T>
    void updateKey(std::string const & key, T const & value);
    void updateKey(std::string const & key, char const * value) {
        updateKey(key, std::string(value));
    }


    template <typename T>
    void writeKey(std::string const & key, T const & value, std::string const & comment);
    void writeKey(std::string const & key, char const * value, std::string const & comment) {
        writeKey(key, std::string(value), comment);
    }
    template <typename T>
    void writeKey(std::string const & key, T const & value);
    void writeKey(std::string const & key, char const * value) {
        writeKey(key, std::string(value));
    }

    template <typename T>
    void updateColumnKey(std::string const & prefix, int n, T const & value, std::string const & comment);
    void updateColumnKey(std::string const & prefix, int n, char const * value, std::string const & comment) {
        updateColumnKey(prefix, n, std::string(value), comment);
    }
    template <typename T>
    void updateColumnKey(std::string const & prefix, int n, T const & value);
    void updateColumnKey(std::string const & prefix, int n, char const * value) {
        updateColumnKey(prefix, n, std::string(value));
    }

    template <typename T>
    void writeColumnKey(std::string const & prefix, int n, T const & value, std::string const & comment);
    void writeColumnKey(std::string const & prefix, int n, char const * value, std::string const & comment) {
        writeColumnKey(prefix, n, std::string(value), comment);
    }
    template <typename T>
    void writeColumnKey(std::string const & prefix, int n, T const & value);
    void writeColumnKey(std::string const & prefix, int n, char const * value) {
        writeColumnKey(prefix, n, std::string(value));
    }

    void writeMetadata(daf::base::PropertySet const & metadata);

    void readMetadata(daf::base::PropertySet & metadata, bool strip=false);

    template <typename T>
    void readKey(std::string const & key, T & value);

    void forEachKey(HeaderIterationFunctor & functor);

    void createEmpty();

    template <typename PixelT, int N>
    void createImage(ndarray::Vector<ndarray::Size,N> const & shape) {
        ndarray::Vector<long,N> nAxes(shape.reverse());
        createImageImpl<PixelT>(N, nAxes.elems);
    }

    template <typename PixelT>
    void createImage(long x, long y) {
        long naxes[2] = { x, y };
        createImageImpl<PixelT>(2, naxes);
    }

    template <typename T, int N, int C>
    void writeImage(ndarray::Array<T const,N,C> const & array) {
        ndarray::Array<T const,N,N> contiguous = ndarray::dynamic_dimension_cast<2>(array);
        if (contiguous.empty()) contiguous = ndarray::copy(array);
        writeImageImpl(contiguous.getData(), contiguous.getNumElements());
    }

    int getImageDim();

    template <int N>
    ndarray::Vector<ndarray::Size,N> getImageShape() {
        ndarray::Vector<long,N> nAxes(1);
        getImageShapeImpl(N, nAxes.elems);
        ndarray::Vector<ndarray::Size,N> shape;
        for (int i = 0; i < N; ++i) shape[i] = nAxes[N-i-1];
        return shape;
    }

    template <typename T>
    bool checkImageType();

    template <typename T, int N>
    void readImage(
        ndarray::Array<T,N,N> const & array,
        ndarray::Vector<int,N> const & offset
    ) {
        ndarray::Vector<long,N> begin(offset.reverse());
        ndarray::Vector<long,N> end(begin);
        end += array.getShape().reverse();
        ndarray::Vector<long,N> increment(1);
        begin += increment;  // first FITS pixel is 1, not 0
        readImageImpl(N, array.getData(), begin.elems, end.elems, increment.elems);
    }

    void createTable();

    template <typename T>
    int addColumn(std::string const & ttype, int size, std::string const & comment);

    template <typename T>
    int addColumn(std::string const & ttype, int size);

    std::size_t addRows(std::size_t nRows);

    std::size_t countRows();

    template <typename T>
    void writeTableArray(std::size_t row, int col, int nElements, T const * value);

    template <typename T>
    void writeTableScalar(std::size_t row, int col, T value) { writeTableArray(row, col, 1, &value); }

    void writeTableScalar(std::size_t row, int col, std::string const & value);

    template <typename T>
    void readTableArray(std::size_t row, int col, int nElements, T * value);

    template <typename T>
    void readTableScalar(std::size_t row, int col, T & value) { readTableArray(row, col, 1, &value); }

    void readTableScalar(std::size_t row, int col, std::string & value);

    long getTableArraySize(int col);

    long getTableArraySize(std::size_t row, int col);

    Fits() : fptr(0), status(0), behavior(0) {}

    Fits(std::string const & filename, std::string const & mode, int behavior);

    Fits(MemFileManager & manager, std::string const & mode, int behavior);

    void closeFile();

    ~Fits() { if ((fptr) && (behavior & AUTO_CLOSE)) closeFile(); }

    // No copying
    Fits (const Fits&) = delete;
    Fits& operator=(const Fits&) = delete;

    // No moving
    Fits (Fits&&) = delete;
    Fits& operator=(Fits&&) = delete;

    void * fptr;  // the actual cfitsio fitsfile pointer; void to avoid including fitsio.h here.
    int status;   // the cfitsio status indicator that gets passed to every cfitsio call.
    int behavior; // bitwise OR of BehaviorFlags
};

#endif // !SWIG

PTR(daf::base::PropertyList) readMetadata(std::string const & fileName, int hdu=0, bool strip=false);
PTR(daf::base::PropertyList) readMetadata(fits::MemFileManager & manager, int hdu=0, bool strip=false);
PTR(daf::base::PropertyList) readMetadata(fits::Fits & fitsfile, bool strip=false);

}}} 

#endif // !LSST_AFW_fits_h_INCLUDED