SdssShape.h

Go to the documentation of this file.

// -*- lsst-c++ -*-
/*
 * LSST Data Management System
 * Copyright 2008-2016 AURA/LSST.
 *
 * 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 <https://www.lsstcorp.org/LegalNotices/>.
 */

#ifndef LSST_MEAS_BASE_SdssShape_h_INCLUDED
#define LSST_MEAS_BASE_SdssShape_h_INCLUDED

#include <bitset>

#include "lsst/pex/config.h"
#include "lsst/afw/image/Exposure.h"
#include "lsst/afw/geom/ellipses/Quadrupole.h"
#include "lsst/afw/table/aggregates.h"
#include "lsst/meas/base/FluxUtilities.h"
#include "lsst/meas/base/CentroidUtilities.h"
#include "lsst/meas/base/ShapeUtilities.h"
#include "lsst/meas/base/InputUtilities.h"
#include "lsst/meas/base/Algorithm.h"

namespace lsst {
namespace meas {
namespace base {

class SdssShapeResult;

class SdssShapeControl {
public:
    LSST_CONTROL_FIELD(background, double, "Additional value to add to background");
    LSST_CONTROL_FIELD(maxIter, int, "Maximum number of iterations");
    LSST_CONTROL_FIELD(maxShift, double, "Maximum centroid shift, limited to 2-10");
    LSST_CONTROL_FIELD(tol1, float, "Convergence tolerance for e1,e2");
    LSST_CONTROL_FIELD(tol2, float, "Convergence tolerance for FWHM");
    LSST_CONTROL_FIELD(doMeasurePsf, bool, "Whether to also compute the shape of the PSF model");

    SdssShapeControl() : background(0.0), maxIter(100), maxShift(), tol1(1E-5), tol2(1E-4),
                         doMeasurePsf(true) {}
};

class SdssShapeResultKey : public afw::table::FunctorKey<SdssShapeResult> {
public:

    static SdssShapeResultKey addFields(
        afw::table::Schema & schema,
        std::string const & name,
        bool doMeasurePsf
    );

    SdssShapeResultKey() {}

    SdssShapeResultKey(afw::table::SubSchema const & s);

    virtual SdssShapeResult get(afw::table::BaseRecord const & record) const;

    virtual void set(afw::table::BaseRecord & record, SdssShapeResult const & value) const;

    virtual afw::geom::ellipses::Quadrupole getPsfShape(afw::table::BaseRecord const & record) const;

    virtual void setPsfShape(afw::table::BaseRecord & record,
                             afw::geom::ellipses::Quadrupole const & value) const;

    bool operator==(SdssShapeResultKey const & other) const;
    bool operator!=(SdssShapeResultKey const & other) const { return !(*this == other); }

    bool isValid() const;

    FlagHandler const & getFlagHandler() const { return _flagHandler; }

private:
    bool _includePsf;
    ShapeResultKey _shapeResult;
    CentroidResultKey _centroidResult;
    FluxResultKey _fluxResult;
    afw::table::QuadrupoleKey _psfShapeResult;
    afw::table::Key<ErrElement> _flux_xx_Cov;
    afw::table::Key<ErrElement> _flux_yy_Cov;
    afw::table::Key<ErrElement> _flux_xy_Cov;
    FlagHandler _flagHandler;
};

class SdssShapeAlgorithm : public SimpleAlgorithm {
public:

    typedef SdssShapeControl Control;
    typedef SdssShapeResult Result;
    typedef SdssShapeResultKey ResultKey;

    // NOTE: In order to accommodate the optional setting of additional fields when running with
    //       doMeasurePsf = true (do set extra fields) or false (do NOT set extra fields), all of
    //       the code in SdssShape assumes that PSF_SHAPE_BAD is the last entry in the enum list.
    //       If new flags are added, be sure to add them above the PSF_SHAPE_BAD entry.
    enum {
        FAILURE=FlagHandler::FAILURE,
        UNWEIGHTED_BAD,
        UNWEIGHTED,
        SHIFT,
        MAXITER,
        PSF_SHAPE_BAD,  // NOTE: PSF_SHAPE_BAD must be the last entry in the enum list
        N_FLAGS
    };

    SdssShapeAlgorithm(Control const & ctrl, std::string const & name, afw::table::Schema & schema);

    template <typename ImageT>
    static Result computeAdaptiveMoments(
        ImageT const & image,
        afw::geom::Point2D const & position,
        bool negative=false,
        Control const & ctrl=Control()
    );

    template <typename ImageT>
    static FluxResult computeFixedMomentsFlux(
        ImageT const & image,
        afw::geom::ellipses::Quadrupole const & shape,
        afw::geom::Point2D const & position
    );

    virtual void measure(
        afw::table::SourceRecord & measRecord,
        afw::image::Exposure<float> const & exposure
    ) const;

    virtual void fail(
        afw::table::SourceRecord & measRecord,
        MeasurementError * error=NULL
    ) const;

private:
    Control _ctrl;
    ResultKey _resultKey;
    SafeCentroidExtractor _centroidExtractor;
};

class SdssShapeResult : public ShapeResult, public CentroidResult, public FluxResult {
public:
    ErrElement flux_xx_Cov; 
    ErrElement flux_yy_Cov; 
    ErrElement flux_xy_Cov; 

#ifndef SWIG
    std::bitset<SdssShapeAlgorithm::N_FLAGS> flags; 
#endif

    bool getFlag(int index) const { return flags[index]; }

    SdssShapeResult(); 

};

class SdssShapeTransform : public BaseTransform {
public:
    typedef SdssShapeControl Control;

    SdssShapeTransform(Control const & ctrl, std::string const & name, afw::table::SchemaMapper & mapper);

    /*
     * @brief Perform transformation from inputCatalog to outputCatalog.
     *
     * @param[in]     inputCatalog   Source of data to be transformed
     * @param[in,out] outputCatalog  Container for transformed results
     * @param[in]     wcs            World coordinate system under which transformation will take place
     * @param[in]     calib          Photometric calibration under which transformation will take place
     * @throws        LengthError    Catalog sizes do not match
     */
    virtual void operator()(afw::table::SourceCatalog const & inputCatalog,
                            afw::table::BaseCatalog & outputCatalog,
                            afw::image::Wcs const & wcs,
                            afw::image::Calib const & calib) const;
private:
    FluxTransform _fluxTransform;
    CentroidTransform _centroidTransform;
    ShapeResultKey _outShapeKey;
    afw::table::QuadrupoleKey _outPsfShapeKey;
    bool _transformPsf;
};

}}} // namespace lsst::meas::base

#endif // !LSST_MEAS_BASE_SdssShape_h_INCLUDED