CreateWcsWithSip.h

Go to the documentation of this file.

// -*- LSST-C++ -*-

/* 
 * LSST Data Management System
 * Copyright 2008, 2009, 2010 LSST Corporation.
 * 
 * This product includes software developed by the
 * LSST Project (http://www.lsst.org/).
 *
 * 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 3 of the License, or
 * (at your option) any later version.
 * 
 * This program 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 LSST License Statement and 
 * the GNU General Public License along with this program.  If not, 
 * see <http://www.lsstcorp.org/LegalNotices/>.
 */
 

#ifndef CREATE_WCS_WITH_SIP
#define CREATE_WCS_WITH_SIP

#include <vector>

#include "boost/shared_ptr.hpp"
#include "lsst/base.h"
#include "Eigen/Core"

#include "lsst/afw/table/Match.h"
#include "lsst/afw/geom/Box.h"
#include "lsst/afw/geom/Point.h"
#include "lsst/afw/geom/Angle.h"
#include "lsst/pex/logging/Log.h"

namespace lsst { 
    namespace afw {
        namespace image {
            class Wcs;
            class TanWcs;
        }
    }

namespace meas { 
namespace astrom { 
namespace sip {


template<class MatchT>
class CreateWcsWithSip {
public:

    typedef boost::shared_ptr<CreateWcsWithSip> Ptr;
    typedef boost::shared_ptr<CreateWcsWithSip const> ConstPtr;

    CreateWcsWithSip(
        std::vector<MatchT> const & matches,
        afw::image::Wcs const & linearWcs,
        int const order,
        afw::geom::Box2I const& bbox = afw::geom::Box2I(),
        int const ngrid=0
    );

    PTR(afw::image::TanWcs) getNewWcs() { return _newWcs; }

    double getScatterInPixels() const;

    afw::geom::Angle getScatterOnSky() const;

    double getLinearScatterInPixels() const;

    afw::geom::Angle getLinearScatterOnSky() const;

    int getOrder() const { return  _sipA.rows(); }
    int getNPoints() const { return _matches.size(); }
    int getNGrid() const { return _ngrid; }

private:

    lsst::pex::logging::Log _log;
    
    std::vector<MatchT> const _matches;
    afw::geom::Box2I mutable _bbox;
    int _ngrid;                         // grid size to calculate inverse SIP coefficients (1-D)
    CONST_PTR(afw::image::Wcs) _linearWcs;
    // _sipOrder is polynomial order for forward transform.
    // _reverseSipOrder is order for reverse transform, not necessarily the same.
    int const _sipOrder, _reverseSipOrder;      

    Eigen::MatrixXd _sipA, _sipB;
    Eigen::MatrixXd _sipAp, _sipBp;

    PTR(afw::image::TanWcs) _newWcs;

    void _calculateForwardMatrices();
    void _calculateReverseMatrices();
    
    afw::geom::Point2D _getCrvalAsGeomPoint() const;
};    

template<class MatchT>
CreateWcsWithSip<MatchT> makeCreateWcsWithSip(
    std::vector<MatchT> const & matches,
    afw::image::Wcs const& linearWcs,
    int const order,
    afw::geom::Box2I const& bbox = afw::geom::Box2I(),
    int const ngrid=0
) {
    return CreateWcsWithSip<MatchT>(matches, linearWcs, order, bbox, ngrid);
}

}}}}

#endif