Measure the image moments of source using adaptive Gaussian weights. More...

#include <SdssShape.h>

Inheritance diagram for lsst::meas::base::SdssShapeAlgorithm:

Public Types
typedef SdssShapeControl	Control

typedef SdssShapeResult	Result

typedef SdssShapeResultKey	ResultKey

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

virtual void	measure (afw::table::SourceRecord &measRecord, afw::image::Exposure< float > const &exposure) const
	Called to measure a single child source in an image.

virtual void	fail (afw::table::SourceRecord &measRecord, MeasurementError *error=nullptr) const
	Handle an exception thrown by the current algorithm by setting flags in the given record.

template<typename ImageT >
SdssShapeResult	computeAdaptiveMoments (ImageT const &image, geom::Point2D const &center, bool negative, Control const &control)

void	measureForced (afw::table::SourceRecord &measRecord, afw::image::Exposure< float > const &exposure, afw::table::SourceRecord const &refRecord, afw::geom::SkyWcs const &refWcs) const override
	Called to measure a single child source in an image.

void	measureNForced (afw::table::SourceCatalog const &measCat, afw::image::Exposure< float > const &exposure, afw::table::SourceCatalog const &refRecord, afw::geom::SkyWcs const &refWcs) const override
	Called to simultaneously measure all children in a deblend family, in a single image.

virtual void	measureN (afw::table::SourceCatalog const &measCat, afw::image::Exposure< float > const &exposure) const
	Called to simultaneously measure all children in a deblend family, in a single image.

std::string	getLogName () const

Static Public Member Functions
static FlagDefinitionList const &	getFlagDefinitions ()

template<typename ImageT >
static Result	computeAdaptiveMoments (ImageT const &image, geom::Point2D const &position, bool negative=false, Control const &ctrl=Control())
	Compute the adaptive Gaussian-weighted moments of an image.

template<typename ImageT >
static FluxResult	computeFixedMomentsFlux (ImageT const &image, afw::geom::ellipses::Quadrupole const &shape, geom::Point2D const &position)
	Compute the instFlux within a fixed Gaussian aperture.

Static Public Attributes
static unsigned int const	N_FLAGS = 6

static FlagDefinition const	FAILURE = flagDefinitions.addFailureFlag()

static FlagDefinition const	UNWEIGHTED_BAD

static FlagDefinition const	UNWEIGHTED

static FlagDefinition const	SHIFT

static FlagDefinition const	MAXITER

static FlagDefinition const	PSF_SHAPE_BAD

Protected Attributes
std::string	_logName

Detailed Description

Measure the image moments of source using adaptive Gaussian weights.

This algorithm measures the weighted second moments of an image using a Gaussian weight function, which is iteratively updated to match the current weights. If this iteration does not converge, it can fall back to using unweighted moments, which can be significantly noisier.

See Bernstein & Jarvis, 2002, for more information on this type of algorithm. Note that the code here makes no attempt to correct for the PSF; for PSF corrected ellipticities using weighted moments please use the shapeHSM package.

Definition at line 150 of file SdssShape.h.

Member Typedef Documentation

◆ Control

typedef SdssShapeControl lsst::meas::base::SdssShapeAlgorithm::Control

Definition at line 162 of file SdssShape.h.

◆ Result

typedef SdssShapeResult lsst::meas::base::SdssShapeAlgorithm::Result

Definition at line 163 of file SdssShape.h.

◆ ResultKey

typedef SdssShapeResultKey lsst::meas::base::SdssShapeAlgorithm::ResultKey

Definition at line 164 of file SdssShape.h.

Constructor & Destructor Documentation

◆ SdssShapeAlgorithm()

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

Definition at line 759 of file SdssShape.cc.

        : _ctrl(ctrl),
          _resultKey(ResultKey::addFields(schema, name, ctrl.doMeasurePsf)),
          _centroidExtractor(schema, name) {}

Member Function Documentation

◆ computeAdaptiveMoments() [1/2]

template<typename ImageT >

SdssShapeResult lsst::meas::base::SdssShapeAlgorithm::computeAdaptiveMoments	(	ImageT const &	image,
		geom::Point2D const &	center,
		bool	negative,
		Control const &	control )

Definition at line 766 of file SdssShape.cc.

                                                                                                  {
    double xcen = center.getX();  // object's column position
    double ycen = center.getY();  // object's row position
 
    xcen -= image.getX0();  // work in image Pixel coordinates
    ycen -= image.getY0();
 
    float shiftmax = control.maxShift;  // Max allowed centroid shift
    if (shiftmax < 2) {
        shiftmax = 2;
    } else if (shiftmax > 10) {
        shiftmax = 10;
    }
 
    SdssShapeResult result;
    try {
        result.flags[FAILURE.number] =
                !getAdaptiveMoments(image, control.background, xcen, ycen, shiftmax, &result, control.maxIter,
                                    control.tol1, control.tol2, negative);
    } catch (pex::exceptions::Exception &err) {
        result.flags[FAILURE.number] = true;
    }
    if (result.flags[UNWEIGHTED.number] || result.flags[SHIFT.number]) {
        // These are also considered fatal errors in terms of the quality of the results,
        // even though they do produce some results.
        result.flags[FAILURE.number] = true;
    }
    double IxxIyy = result.getQuadrupole().getIxx() * result.getQuadrupole().getIyy();
    double Ixy_sq = result.getQuadrupole().getIxy();
    Ixy_sq *= Ixy_sq;
    double epsilon = 1.0e-6;
    if (IxxIyy < (1.0 + epsilon) * Ixy_sq)
    // We are checking that Ixx*Iyy > (1 + epsilon)*Ixy*Ixy where epsilon is suitably small. The
    // value of epsilon used here is a magic number subject to future review (DM-5801 was marked won't fix).
    {
        if (!result.flags[FAILURE.number]) {
            throw LSST_EXCEPT(pex::exceptions::LogicError,
                              (boost::format("computeAdaptiveMoments IxxIxy %d < (1 + eps=%d)*(Ixy^2=%d);"
                                             " implying singular moments without any flag set")
                              % IxxIyy % epsilon % Ixy_sq).str());
        }
    }
 
    // getAdaptiveMoments() just computes the zeroth moment in result.instFlux (and its error in
    // result.instFluxErr, result.instFlux_xx_Cov, etc.)  That's related to the instFlux by some geometric
    // factors, which we apply here.
    // This scale factor is just the inverse of the normalization constant in a bivariate normal distribution:
    // 2*pi*sigma_x*sigma_y*sqrt(1-rho^2), where rho is the correlation.
    // This happens to be twice the ellipse area pi*sigma_maj*sigma_min, which is identically equal to:
    // pi*sqrt(I_xx*I_yy - I_xy^2); i.e. pi*sqrt(determinant(I))
    double instFluxScale = geom::TWOPI * std::sqrt(IxxIyy - Ixy_sq);
 
    result.instFlux *= instFluxScale;
    result.instFluxErr *= instFluxScale;
    result.x += image.getX0();
    result.y += image.getY0();
 
    if (ImageAdaptor<ImageT>::hasVariance) {
        result.instFlux_xx_Cov *= instFluxScale;
        result.instFlux_yy_Cov *= instFluxScale;
        result.instFlux_xy_Cov *= instFluxScale;
    }
 
    return result;
}

◆ computeAdaptiveMoments() [2/2]

template<typename ImageT >

static Result lsst::meas::base::SdssShapeAlgorithm::computeAdaptiveMoments	(	ImageT const &	image,
		geom::Point2D const &	position,
		bool	negative = false,
		Control const &	ctrl = Control() )

static

Compute the adaptive Gaussian-weighted moments of an image.

Parameters

[in]	image	An Image or MaskedImage instance with int, float, or double pixels. This need not be a small postage stamp (the pixel region actually used in the fit will be a subset of this image determined automatically).
[in]	position	Center position of the object to be measured, in the image's PARENT coordinates.
[in]	negative	Boolean, specify if the source is in negative instFlux space
[in]	ctrl	Control object specifying the details of how the object is to be measured.

◆ computeFixedMomentsFlux()

template<typename ImageT >

template FluxResult lsst::meas::base::SdssShapeAlgorithm::computeFixedMomentsFlux	(	ImageT const &	image,
		afw::geom::ellipses::Quadrupole const &	shape,
		geom::Point2D const &	position )

static

Compute the instFlux within a fixed Gaussian aperture.

Parameters

[in]	image	An Image or MaskedImage instance with int, float, or double pixels. This need not be a small postage stamp (the pixel region actually used in the fit will be a subset of this image determined automatically).
[in]	shape	Ellipse object specifying the 1-sigma contour of the Gaussian.
[in]	position	Center position of the object to be measured, in the image's PARENT coordinates.

Definition at line 834 of file SdssShape.cc.

                                                                                  {
    // while arguments to computeFixedMomentsFlux are in PARENT coordinates, the implementation is LOCAL.
    geom::Point2D localCenter = center - geom::Extent2D(image.getXY0());
 
    geom::BoxI const bbox = computeAdaptiveMomentsBBox(image.getBBox(afw::image::LOCAL), localCenter,
                                                       shape.getIxx(), shape.getIxy(), shape.getIyy());
 
    std::tuple<std::pair<bool, double>, double, double, double> weights =
            getWeights(shape.getIxx(), shape.getIxy(), shape.getIyy());
 
    FluxResult result;
 
    if (!std::get<0>(weights).first) {
        throw pex::exceptions::InvalidParameterError("Input shape is singular");
    }
 
    double const w11 = std::get<1>(weights);
    double const w12 = std::get<2>(weights);
    double const w22 = std::get<3>(weights);
    bool const interp = shouldInterp(shape.getIxx(), shape.getIyy(), std::get<0>(weights).second);
 
    double sum0 = 0;  //  sum of pixel values weighted by a Gaussian
    calcmom(ImageAdaptor<ImageT>().getImage(image), localCenter.getX(), localCenter.getY(), bbox,
                      0.0, interp, w11, w12, w22, sum0);
 
    result.instFlux = sum0 * 2.0;
 
    if (ImageAdaptor<ImageT>::hasVariance) {
        int ix = static_cast<int>(center.getX() - image.getX0());
        int iy = static_cast<int>(center.getY() - image.getY0());
        if (!image.getBBox(afw::image::LOCAL).contains(geom::Point2I(ix, iy))) {
            throw LSST_EXCEPT(pex::exceptions::RuntimeError,
                              (boost::format("Center (%d,%d) not in image (%dx%d)") % ix % iy %
                               image.getWidth() % image.getHeight())
                                      .str());
        }
        double var = ImageAdaptor<ImageT>().getVariance(image, ix, iy);
        // 0th moment (i0) error = sqrt(var / wArea); instFlux (error) = 2 * wArea * i0 (error)
        double const wArea = geom::PI * std::sqrt(shape.getDeterminant());
        result.instFluxErr = 2 * std::sqrt(var * wArea);
    }
 
    return result;
}

◆ fail()

void lsst::meas::base::SdssShapeAlgorithm::fail	(	afw::table::SourceRecord &	measRecord,
		MeasurementError *	error = nullptr ) const

virtual

Handle an exception thrown by the current algorithm by setting flags in the given record.

fail() is called by the measurement framework when an exception is allowed to propagate out of one the algorithm's measure() methods. It should generally set both a general failure flag for the algorithm as well as a specific flag indicating the error condition, if possible. To aid in this, if the exception was an instance of MeasurementError, it will be passed in, carrying information about what flag to set.

An algorithm can also to chose to set flags within its own measure() methods, and then just return, rather than throw an exception. However, fail() should be implemented even when all known failure modes do not throw exceptions, to ensure that unexpected exceptions thrown in lower-level code are properly handled.

Implements lsst::meas::base::BaseAlgorithm.

Definition at line 916 of file SdssShape.cc.

                                                                                             {
    _resultKey.getFlagHandler().handleFailure(measRecord, error);
}

◆ getFlagDefinitions()

FlagDefinitionList const & lsst::meas::base::SdssShapeAlgorithm::getFlagDefinitions ( )

static

Definition at line 58 of file SdssShape.cc.

58{ return flagDefinitions; }

◆ getLogName()

std::string lsst::meas::base::BaseAlgorithm::getLogName ( ) const

inlineinherited

Definition at line 66 of file Algorithm.h.

66{ return _logName; }

lsst::meas::base::BaseAlgorithm::_logName

std::string _logName

Definition Algorithm.h:69

◆ measure()

void lsst::meas::base::SdssShapeAlgorithm::measure	(	afw::table::SourceRecord &	measRecord,
		afw::image::Exposure< float > const &	exposure ) const

virtual

Called to measure a single child source in an image.

Before this method is called, all neighbors will be replaced with noise, using the outputs of the deblender. Outputs should be saved in the given SourceRecord, which can also be used to obtain centroid (see SafeCentroidExtractor) and shape (see SafeShapeExtractor) information.

Implements lsst::meas::base::SingleFrameAlgorithm.

Definition at line 881 of file SdssShape.cc.

                                                                                {
    bool negative = false;
 
    try {
        negative = measRecord.get(measRecord.getSchema().find<afw::table::Flag>("flags_negative").key);
    } catch (pexExcept::Exception &e) {
    }
    SdssShapeResult result = computeAdaptiveMoments(
            exposure.getMaskedImage(), _centroidExtractor(measRecord, _resultKey.getFlagHandler()), negative,
            _ctrl);
 
    if (_ctrl.doMeasurePsf) {
        // Compute moments of Psf model.  In the interest of implementing this quickly, we're just
        // calling Psf::computeShape(), which delegates to SdssShapeResult::computeAdaptiveMoments
        // for all nontrivial Psf classes.  But this could in theory save the results of a shape
        // computed some other way as part of base_SdssShape, which might be confusing.  We should
        // fix this eventually either by making Psf shape measurement not part of base_SdssShape, or
        // by making the measurements stored with shape.sdss always computed via the
        // SdssShapeAlgorithm instead of delegating to the Psf class.
        try {
            std::shared_ptr<afw::detection::Psf const> psf = exposure.getPsf();
            if (!psf) {
                result.flags[PSF_SHAPE_BAD.number] = true;
            } else {
                _resultKey.setPsfShape(measRecord, psf->computeShape(geom::Point2D(result.x, result.y)));
            }
        } catch (pex::exceptions::Exception &err) {
            result.flags[PSF_SHAPE_BAD.number] = true;
        }
    }
 
    measRecord.set(_resultKey, result);
}

◆ measureForced()

void lsst::meas::base::SimpleAlgorithm::measureForced	(	afw::table::SourceRecord &	measRecord,
		afw::image::Exposure< float > const &	exposure,
		afw::table::SourceRecord const &	refRecord,
		afw::geom::SkyWcs const &	refWcs ) const

inlineoverridevirtualinherited

Called to measure a single child source in an image.

Before this method is called, all neighbors will be replaced with noise, using the outputs of the deblender. Outputs should be saved in the given SourceRecord, which can also be used to obtain centroid (see SafeCentroidExtractor) and shape (see SafeShapeExtractor) information.

Implements lsst::meas::base::ForcedAlgorithm.

Reimplemented in lsst::meas::extensions::photometryKron::KronFluxAlgorithm.

Definition at line 172 of file Algorithm.h.

                                                                           {
        measure(measRecord, exposure);
    }

◆ measureN()

void lsst::meas::base::SingleFrameAlgorithm::measureN	(	afw::table::SourceCatalog const &	measCat,
		afw::image::Exposure< float > const &	exposure ) const

virtualinherited

Called to simultaneously measure all children in a deblend family, in a single image.

Outputs should be saved in the given SourceCatalog, which can also be used to obtain centroid (see SafeCentroidExtractor) and shape (see SafeShapeExtractor) information.

The default implementation simply throws an exception, indicating that simultaneous measurement is not supported.

Definition at line 31 of file Algorithm.cc.

                                                                                   {
    throw LSST_EXCEPT(pex::exceptions::LogicError, "measureN not implemented for this algorithm");
}

◆ measureNForced()

void lsst::meas::base::SimpleAlgorithm::measureNForced	(	afw::table::SourceCatalog const &	measCat,
		afw::image::Exposure< float > const &	exposure,
		afw::table::SourceCatalog const &	refRecord,
		afw::geom::SkyWcs const &	refWcs ) const

inlineoverridevirtualinherited

Called to simultaneously measure all children in a deblend family, in a single image.

Outputs should be saved in the given SourceCatalog, which can also be used to obtain centroid (see SafeCentroidExtractor) and shape (see SafeShapeExtractor) information.

The default implementation simply throws an exception, indicating that simultaneous measurement is not supported.

Reimplemented from lsst::meas::base::ForcedAlgorithm.

Definition at line 179 of file Algorithm.h.

                                                                            {
        measureN(measCat, exposure);
    }

Member Data Documentation

◆ _logName

std::string lsst::meas::base::BaseAlgorithm::_logName

protectedinherited

Definition at line 69 of file Algorithm.h.

◆ FAILURE

FlagDefinition const lsst::meas::base::SdssShapeAlgorithm::FAILURE = flagDefinitions.addFailureFlag()

static

Definition at line 155 of file SdssShape.h.

◆ MAXITER

FlagDefinition const lsst::meas::base::SdssShapeAlgorithm::MAXITER

static

Initial value:

=

flagDefinitions.add("flag_maxIter", "Too many iterations in adaptive moments")

Definition at line 159 of file SdssShape.h.

◆ N_FLAGS

unsigned int const lsst::meas::base::SdssShapeAlgorithm::N_FLAGS = 6

static

Definition at line 154 of file SdssShape.h.

◆ PSF_SHAPE_BAD

FlagDefinition const lsst::meas::base::SdssShapeAlgorithm::PSF_SHAPE_BAD

static

Initial value:

=

flagDefinitions.add("flag_psf", "Failure in measuring PSF model shape")

Definition at line 160 of file SdssShape.h.

◆ SHIFT

FlagDefinition const lsst::meas::base::SdssShapeAlgorithm::SHIFT

static

Initial value:

=

flagDefinitions.add("flag_shift", "centroid shifted by more than the maximum allowed amount")

Definition at line 158 of file SdssShape.h.

◆ UNWEIGHTED

FlagDefinition const lsst::meas::base::SdssShapeAlgorithm::UNWEIGHTED

static

Initial value:

= flagDefinitions.add(

"flag_unweighted", "Weighted moments converged to an invalid value; using unweighted moments")

Definition at line 157 of file SdssShape.h.

◆ UNWEIGHTED_BAD

FlagDefinition const lsst::meas::base::SdssShapeAlgorithm::UNWEIGHTED_BAD

static

Initial value:

=

flagDefinitions.add("flag_unweightedBad", "Both weighted and unweighted moments were invalid")

Definition at line 156 of file SdssShape.h.

The documentation for this class was generated from the following files:

/j/snowflake/release/lsstsw/stack/lsst-scipipe-9.0.0/Linux64/meas_base/g217e2c1bcf+937a289c59/include/lsst/meas/base/SdssShape.h
/j/snowflake/release/lsstsw/stack/lsst-scipipe-9.0.0/Linux64/meas_base/g217e2c1bcf+937a289c59/src/SdssShape.cc

Public Types

Public Member Functions

Static Public Member Functions

Static Public Attributes

Protected Attributes

Detailed Description

Member Typedef Documentation

◆ Control

◆ Result

◆ ResultKey

Constructor & Destructor Documentation

◆ SdssShapeAlgorithm()

Member Function Documentation

◆ computeAdaptiveMoments() [1/2]

◆ computeAdaptiveMoments() [2/2]

◆ computeFixedMomentsFlux()

◆ fail()

◆ getFlagDefinitions()

◆ getLogName()

◆ measure()

◆ measureForced()

◆ measureN()

◆ measureNForced()

Member Data Documentation

◆ _logName

◆ FAILURE

◆ MAXITER

◆ N_FLAGS

◆ PSF_SHAPE_BAD

◆ SHIFT

◆ UNWEIGHTED

◆ UNWEIGHTED_BAD