File LeastSqFitter2d.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
astrom
-
namespace
sip
-
template<class
FittingFunc
>
classLeastSqFitter2d
- #include <LeastSqFitter2d.h>
Fit an lsst::afw::math::Function1 object to a set of data points in two dimensions.
The class is templated over the kind of object to fit. Note that we fit the 1d object in each dimension, not the 2d one.
Input is a list of x and y ordinates for a set of points, the z coordinate, and the uncertainties, s. order is order of the polynomial to fit (e.g if the templated function is lsst::afw::math::PolynomialFunction1, then order=3 => fit a function of the form \(ax^2+bx+c\)
- See
- Template Parameters
FittingFunc
: The 1d function to fit in both dimensions. Must inherit from lsst::afw::math::Function1
- Parameters
x
: Ordinate of points to fity
: Ordinate of points to fitz
: Co-ordinate of pionts to fits
: 1 \(\sigma\) uncertainties in zorder
: Polynomial order to fit
Public Functions
-
LeastSqFitter2d
(const std::vector<double> &x, const std::vector<double> &y, const std::vector<double> &z, const std::vector<double> &s, int order) Fit a 2d polynomial to a set of data points z(x, y)
- Template Parameters
FittingFunc
: The type of function to fit in each dimension. This function should extend the base class of lsst::afw::math::Function1. Note that this is a 1d function
- Parameters
x
: vector of x positions of datay
: vector of y positions of dataz
: Value of data for a given x,y. \(z = z_i = z_i(x_i, y_i)\)s
: Vector of measured uncertainties in the values of zorder
: Order of 2d function to fit
-
Eigen::MatrixXd
getParams
() Build up a triangular matrix of the parameters. The shape of the matrix is such that the values correspond to the coefficients of the following polynomials
1 y y^2 y^3
x xy xy^2 0
x^2 x^2y 0 0
x^3 0 0 0
(order==4)
where row x column < order
-
Eigen::MatrixXd
getErrors
() Companion function to getParams(). Returns uncertainties in the parameters as a matrix.
-
double
valueAt
(double x, double y) Return the value of the best fit function at a given position (x,y)
-
std::vector<double>
residuals
() Return a vector of residuals of the fit (i.e the difference between the input z values, and the value of the fitting function at that point.
-
double
getChiSq
() Return a measure of the goodness of fit.
\[ \chi^2 = \sum \left( \frac{z_i - f(x_i, y_i)}{s_i} \right)^2 \].
-
double
getReducedChiSq
() Return a measure of the goodness of fit.
\[ \chi_r^2 = \sum \left( \frac{z_i - f(x_i, y_i)}{s} \right)^2 \div (N-p) \]Where \( N \) is the number of data points, and \( p \) is the number of parameters in the fit.
-
template<class
-
namespace
-
namespace