GridTransform.cc
Go to the documentation of this file.
1 // -*- lsst-c++ -*-
2
3 /*
4  * LSST Data Management System
5  * Copyright 2008, 2009, 2010 LSST Corporation.
6  *
7  * This product includes software developed by the
8  * LSST Project (http://www.lsst.org/).
9  *
10  * This program is free software: you can redistribute it and/or modify
12  * the Free Software Foundation, either version 3 of the License, or
13  * (at your option) any later version.
14  *
15  * This program is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  * GNU General Public License for more details.
19  *
20  * You should have received a copy of the LSST License Statement and
21  * the GNU General Public License along with this program. If not,
22  * see <http://www.lsstcorp.org/LegalNotices/>.
23  */
29
30 namespace lsst {
31 namespace afw {
32 namespace geom {
33 namespace ellipses {
34
37
39  return _eig.operatorInverseSqrt();
40 }
41
42 BaseCore::GridTransform::operator lsst::geom::LinearTransform() const {
43  return lsst::geom::LinearTransform(_eig.operatorInverseSqrt());
44 }
45
47  /*
48  Grid transform is easiest to differentiate in the ReducedShear/DeterminantRadius parametrization.
49  But we actually differentiate the inverse of the transform, and then use
50  $dM^{-1}/dt = -M^{-1} dM/dt M^{-1} to compute the derivative of the inverse. 51 52 The inverse of the grid transform in ReducedShear/DeterminantRadius is: 53$\frac{r}{\sqrt{1-g^2}}(\sigma_x + g_1 \sigma_z + g2 \sigma_y)$, where$\sigma_i\$ are the
54  Pauli spin matrices.
55  */
57  C core;
58  Jacobian rhs = core.dAssign(_input);
59  double g1 = core.getE1();
60  double g2 = core.getE2();
61  double g = core.getEllipticity().getE();
63  double beta = 1.0 - g * g;
64  double alpha = r / std::sqrt(beta);
65
66  Eigen::Matrix2d sigma_z, sigma_y;
67  sigma_z << 1.0, 0.0, 0.0, -1.0;
68  sigma_y << 0.0, 1.0, 1.0, 0.0;
69  Eigen::Matrix2d t = _eig.operatorSqrt();
70  Eigen::Matrix2d tInv = _eig.operatorInverseSqrt();
71  Eigen::Matrix2d dt_dg1 = t * g1 / beta + alpha * sigma_z;
72  Eigen::Matrix2d dt_dg2 = t * g2 / beta + alpha * sigma_y;
73  Eigen::Matrix2d dt_dr = t * (1.0 / r);
74  Eigen::Matrix2d dtInv_dg1 = -tInv * dt_dg1 * tInv;
75  Eigen::Matrix2d dtInv_dg2 = -tInv * dt_dg2 * tInv;
76  Eigen::Matrix2d dtInv_dr = -tInv * dt_dr * tInv;
77
79  mid(lsst::geom::LinearTransform::XX, C::E1) = dtInv_dg1(0, 0);
81  dtInv_dg1(0, 1);
82  mid(lsst::geom::LinearTransform::YY, C::E1) = dtInv_dg1(1, 1);
83  mid(lsst::geom::LinearTransform::XX, C::E2) = dtInv_dg2(0, 0);
85  dtInv_dg2(0, 1);
86  mid(lsst::geom::LinearTransform::YY, C::E2) = dtInv_dg2(1, 1);
87  mid(lsst::geom::LinearTransform::XX, C::RADIUS) = dtInv_dr(0, 0);
89  dtInv_dr(0, 1);
90  mid(lsst::geom::LinearTransform::YY, C::RADIUS) = dtInv_dr(1, 1);
91  return mid * rhs;
92 }
93
94 double BaseCore::GridTransform::getDeterminant() const { return sqrt(1.0 / _eig.eigenvalues().prod()); }
95
97  return lsst::geom::LinearTransform(_eig.operatorSqrt());
98 }
99
100 Ellipse::GridTransform::GridTransform(Ellipse const& input) : _input(input), _coreGt(input.getCore()) {}
101
103  lsst::geom::AffineTransform::Matrix r = lsst::geom::AffineTransform::Matrix::Zero();
104  r.block<2, 2>(0, 0) = _coreGt.getMatrix();
105  r.block<2, 1>(0, 2) = -r.block<2, 2>(0, 0) * _input.getCenter().asEigen();
106  r(2, 2) = 1.0;
107  return r;
108 }
109
111  DerivativeMatrix r = DerivativeMatrix::Zero();
112  lsst::geom::LinearTransform linear = _coreGt;
113  r.block<4, 3>(0, 0) = _coreGt.d();
114  double x = -_input.getCenter().getX();
115  double y = -_input.getCenter().getY();
132  return r;
133 }
134
135 double Ellipse::GridTransform::getDeterminant() const { return _coreGt.getDeterminant(); }
136
137 Ellipse::GridTransform::operator lsst::geom::AffineTransform() const {
138  lsst::geom::LinearTransform linear = _coreGt;
139  return lsst::geom::AffineTransform(linear, linear(lsst::geom::Point2D() - _input.getCenter()));
140 }
141
143  return lsst::geom::AffineTransform(_coreGt.inverted(), lsst::geom::Extent2D(_input.getCenter()));
144 }
145 } // namespace ellipses
146 } // namespace geom
147 } // namespace afw
148 } // namespace lsst
double x
int y
Definition: SpanSet.cc:48
double getDeterminant() const
Return the determinant of the lsst::geom::LinearTransform.
lsst::geom::LinearTransform inverted() const
Return the inverse of the lsst::geom::LinearTransform;.
DerivativeMatrix d() const
Return the derivative of the transform with respect to input core.
GridTransform(BaseCore const &input)
Standard constructor.
Eigen::Matrix< double, 4, 3 > DerivativeMatrix
Matrix type for derivative with respect to ellipse parameters.
Definition: GridTransform.h:51
lsst::geom::LinearTransform::Matrix getMatrix() const
Return the transform matrix as an Eigen object.
A base class for parametrizations of the "core" of an ellipse - the ellipticity and size.
Definition: BaseCore.h:55
Jacobian dAssign(BaseCore const &other)
Assign other to this and return the derivative of the conversion, d(this)/d(other).
Definition: BaseCore.cc:168
Eigen::Matrix3d Jacobian
Parameter Jacobian matrix type.
Definition: BaseCore.h:64
double getDeterminant() const
Return the determinant of the lsst::geom::AffineTransform.
DerivativeMatrix d() const
Return the derivative of transform with respect to input ellipse.
lsst::geom::AffineTransform::Matrix getMatrix() const
Return the transform matrix as an Eigen object.
Eigen::Matrix< double, 6, 5 > DerivativeMatrix
Matrix type for derivative with respect to input ellipse parameters.
Definition: GridTransform.h:85
lsst::geom::AffineTransform inverted() const
Return the inverse of the AffineTransform.
GridTransform(Ellipse const &input)
Standard constructor.
An ellipse defined by an arbitrary BaseCore and a center point.
Definition: Ellipse.h:51
An ellipse core with quadrupole moments as parameters.