File SoftenedLinearPrior.h

namespace lsst

Class for a simple mapping implementing a generic AstrometryTransform.

Remove all non-astronomical counts from the Chunk Exposure’s pixels.

Forward declarations for lsst::utils::Cache

For details on the Cache class, see the Cache.h file.

It uses a template rather than a pointer so that the derived classes can use the specifics of the transform. The class simplePolyMapping overloads a few routines.

A base class for image defects

Numeric constants used by the Integrate.h integrator routines.

Compute Image Statistics

Note

Gauss-Kronrod-Patterson quadrature coefficients for use in quadpack routine qng. These coefficients were calculated with 101 decimal digit arithmetic by L. W. Fullerton, Bell Labs, Nov 1981.

Note

The Statistics class itself can only handle lsst::afw::image::MaskedImage() types. The philosophy has been to handle other types by making them look like lsst::afw::image::MaskedImage() and reusing that code. Users should have no need to instantiate a Statistics object directly, but should use the overloaded makeStatistics() factory functions.

namespace meas
namespace modelfit
class SoftenedLinearPrior : public lsst::meas::modelfit::Prior
#include <SoftenedLinearPrior.h>

A prior that’s linear in radius and flat in ellipticity, with a cubic roll-off at the edges.

Public Types

typedef SoftenedLinearPriorControl Control

Public Functions

SoftenedLinearPrior(Control const &ctrl = Control())
Scalar evaluate(ndarray::Array<Scalar const, 1, 1> const &nonlinear, ndarray::Array<Scalar const, 1, 1> const &amplitudes) const

Evaluate the prior at the given point in nonlinear and amplitude space.

Parameters
  • [in] nonlinear: Vector of nonlinear parameters

  • [in] amplitudes: Vector of linear parameters

void evaluateDerivatives(ndarray::Array<Scalar const, 1, 1> const &nonlinear, ndarray::Array<Scalar const, 1, 1> const &amplitudes, ndarray::Array<Scalar, 1, 1> const &nonlinearGradient, ndarray::Array<Scalar, 1, 1> const &amplitudeGradient, ndarray::Array<Scalar, 2, 1> const &nonlinearHessian, ndarray::Array<Scalar, 2, 1> const &amplitudeHessian, ndarray::Array<Scalar, 2, 1> const &crossHessian) const

Evaluate the derivatives of the prior at the given point in nonlinear and amplitude space.

Note that while the model is linear in the amplitudes, the prior is not necessarily linear in the amplitudes, so we do care about second derivatives w.r.t. amplitudes.

Parameters
  • [in] nonlinear: Vector of nonlinear parameters

  • [in] amplitudes: Vector of linear parameters

  • [in] nonlinearGradient: First derivative w.r.t. nonlinear parameters

  • [in] amplitudeGradient: First derivative w.r.t. linear parameters parameters

  • [in] nonlinearHessian: Second derivative w.r.t. nonlinear parameters

  • [in] amplitudeHessian: Second derivative w.r.t. linear parameters parameters

  • [in] crossHessian: Second derivative cross term of d(nonlinear)d(amplitudes); shape is [nonlinearDim, amplitudeDim].

Scalar marginalize(Vector const &gradient, Matrix const &hessian, ndarray::Array<Scalar const, 1, 1> const &nonlinear) const

Return the -log amplitude integral of the prior*likelihood product.

If \(\alpha\) are the amplitudes, \(\theta\) are the nonlinear parameters, and \(D\) is the data, then this method represents \(P(\alpha,\theta)\) by computing

\[ -\ln\left[\int\!P(D|\alpha,\theta)\,P(\alpha,\theta)\,d\alpha\right] \]
at fixed \(\theta\). Because \(\alpha\) are linear parameters, \(P(D|\alpha,\theta)\) is Gaussian in \(\alpha\), and because \(\theta\) is fixed, it’s usually convenient to think of the integral as:
\[ -ln\left[P(\theta)\int\!P(D|\alpha,\theta)\,P(\alpha|\theta)\,d\alpha\right] \]
Thus, we marginalize the likelihood in \(\alpha\) at fixed \(\theta\), and then multiply by the prior on \(\theta\).

We also assume the likelihood \(P(D|\alpha,\theta)\) is Gaussian in \(\alpha\), which is generally true because \(\alpha\) defined such that the model is linear in them, and the noise on the data is generally Gaussian. In detail, we represent the likelihood at fixed \(\theta\) as

\[ P(D|\alpha,\theta) = A e^{-g^T\alpha - \frac{1}{2}\alpha^T H \alpha} \]
The normalization \(A\) can be brought outside the integral as a constant to be added to the return value, so it is not passed as an argument to this function.

Parameters
  • [in] gradient: Gradient of the -log likelihood in \(\alpha\) at fixed \(\theta\); the vector \(g\) in the equation above.

  • [in] hessian: Second derivatives of of the -log likelihood in \(\alpha\) at fixed \(\theta\); the matrix \(H\) in the equation above.

  • [in] nonlinear: The nonlinear parameters \(\theta\).

Scalar maximize(Vector const &gradient, Matrix const &hessian, ndarray::Array<Scalar const, 1, 1> const &nonlinear, ndarray::Array<Scalar, 1, 1> const &amplitudes) const

Compute the amplitude vector that maximizes the prior x likelihood product.

Return

The -log(posterior) at the computed amplitude point.

Parameters
  • [in] gradient: Gradient of the -log likelihood in \(\alpha\) at fixed \(\theta\); the vector \(g\) in the equation above.

  • [in] hessian: Second derivatives of of the -log likelihood in \(\alpha\) at fixed \(\theta\); the matrix \(H\) in the equation above.

  • [in] nonlinear: The nonlinear parameters \(\theta\).

  • [out] amplitudes: The posterior-maximum amplitude parameters \(\alpha\).

void drawAmplitudes(Vector const &gradient, Matrix const &fisher, ndarray::Array<Scalar const, 1, 1> const &nonlinear, afw::math::Random &rng, ndarray::Array<Scalar, 2, 1> const &amplitudes, ndarray::Array<Scalar, 1, 1> const &weights, bool multiplyWeights = false) const

Draw a set of Monte Carlo amplitude vectors.

This provides a Monte Carlo approach to extracting the conditional amplitude distribution that is integrated by the marginalize() method.

Parameters
  • [in] gradient: Gradient of the -log likelihood in \(\alpha\) at fixed \(\theta\).

  • [in] hessian: Second derivatives of of the -log likelihood in \(\alpha\) at fixed \(\theta\).

  • [in] nonlinear: The nonlinear parameters \(\theta\) at which we are evaluating the conditional distribution \(P(\alpha|\theta)\).

  • [inout] rng: Random number generator.

  • [out] amplitudes: The Monte Carlo sample of amplitude parameters \(\alpha\). The number of rows sets the number of samples, while the number of columns must match the dimensionality of \(\alpha\).

  • [out] weights: The weights of the Monte Carlo samples; should asymptotically average to one.

  • [in] multiplyWeights: If true, multiply weight vector instead of overwriting it.

Control const &getControl() const

Private Functions

Scalar _evaluate(ndarray::Array<Scalar const, 1, 1> const &nonlinear) const

Private Members

Control _ctrl
double _logRadiusP1
double _logRadiusSlope
double _logRadiusMinRampFraction
double _logRadiusMaxRampFraction
double _ellipticityMaxRampFraction
Eigen::Matrix<double, 4, 1, Eigen::DontAlign> _logRadiusPoly1
Eigen::Matrix<double, 4, 1, Eigen::DontAlign> _logRadiusPoly2
Eigen::Matrix<double, 4, 1, Eigen::DontAlign> _ellipticityPoly
struct SoftenedLinearPriorControl

Public Functions

lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(ellipticityMaxOuter, double, "Maximum ellipticity magnitude (conformal shear units)")
lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(ellipticityMaxInner, double, "Ellipticity magnitude (conformal shear units) at which the softened cutoff begins")
lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(logRadiusMinOuter, double, "Minimum ln(radius)")
lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(logRadiusMinInner, double, "ln(radius) at which the softened cutoff begins towards the minimum")
lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(logRadiusMaxOuter, double, "Maximum ln(radius)")
lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(logRadiusMaxInner, double, "ln(radius) at which the softened cutoff begins towards the maximum")
lsst::meas::modelfit::SoftenedLinearPriorControl::LSST_CONTROL_FIELD(logRadiusMinMaxRatio, double, "The ratio P(logRadiusMinInner)/P(logRadiusMaxInner)")
SoftenedLinearPriorControl()