LSST Applications  21.0.0-172-gfb10e10a+18fedfabac,22.0.0+297cba6710,22.0.0+80564b0ff1,22.0.0+8d77f4f51a,22.0.0+a28f4c53b1,22.0.0+dcf3732eb2,22.0.1-1-g7d6de66+2a20fdde0d,22.0.1-1-g8e32f31+297cba6710,22.0.1-1-geca5380+7fa3b7d9b6,22.0.1-12-g44dc1dc+2a20fdde0d,22.0.1-15-g6a90155+515f58c32b,22.0.1-16-g9282f48+790f5f2caa,22.0.1-2-g92698f7+dcf3732eb2,22.0.1-2-ga9b0f51+7fa3b7d9b6,22.0.1-2-gd1925c9+bf4f0e694f,22.0.1-24-g1ad7a390+a9625a72a8,22.0.1-25-g5bf6245+3ad8ecd50b,22.0.1-25-gb120d7b+8b5510f75f,22.0.1-27-g97737f7+2a20fdde0d,22.0.1-32-gf62ce7b1+aa4237961e,22.0.1-4-g0b3f228+2a20fdde0d,22.0.1-4-g243d05b+871c1b8305,22.0.1-4-g3a563be+32dcf1063f,22.0.1-4-g44f2e3d+9e4ab0f4fa,22.0.1-42-gca6935d93+ba5e5ca3eb,22.0.1-5-g15c806e+85460ae5f3,22.0.1-5-g58711c4+611d128589,22.0.1-5-g75bb458+99c117b92f,22.0.1-6-g1c63a23+7fa3b7d9b6,22.0.1-6-g50866e6+84ff5a128b,22.0.1-6-g8d3140d+720564cf76,22.0.1-6-gd805d02+cc5644f571,22.0.1-8-ge5750ce+85460ae5f3,master-g6e05de7fdc+babf819c66,master-g99da0e417a+8d77f4f51a,w.2021.48
LSST Data Management Base Package
Namespaces | Classes | Typedefs | Functions
lsst::meas::astrom Namespace Reference

Namespaces

 approximateWcs
 
 astrometry
 
 denormalizeMatches
 
 detail
 
 directMatch
 
 display
 
 fitAffineWcs
 
 fitSipDistortion
 
 fitTanSipWcs
 
 matchOptimisticB
 
 matchOptimisticBTask
 
 matchPessimisticB
 
 pessimistic_pattern_matcher_b_3D
 
 ref_match
 
 setMatchDistance
 
 sip
 
 verifyWcs
 
 version
 

Classes

struct  RecordProxy
 A wrapper around a SimpleRecord or SourceRecord that allows us to record a pixel position in a way that is independent of the record type. More...
 
struct  ProxyPair
 
struct  MatchOptimisticBControl
 
class  PolynomialTransform
 A 2-d coordinate transform represented by a pair of standard polynomials (one for each coordinate). More...
 
class  ScaledPolynomialTransform
 A 2-d coordinate transform represented by a lazy composition of an AffineTransform, a PolynomialTransform, and another AffineTransform. More...
 
class  OutlierRejectionControl
 Control object for outlier rejection in ScaledPolynomialTransformFitter. More...
 
class  ScaledPolynomialTransformFitter
 A fitter class for scaled polynomial transforms. More...
 
class  SipTransformBase
 Base class for SIP transform objects. More...
 
class  SipForwardTransform
 A transform that maps pixel coordinates to intermediate world coordinates according to the SIP convention. More...
 
class  SipReverseTransform
 A transform that maps intermediate world coordinates to pixel coordinates according to the SIP convention. More...
 

Typedefs

typedef std::vector< RecordProxyProxyVector
 

Functions

template<typename MatchT >
afw::math::Statistics makeMatchStatistics (std::vector< MatchT > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl=afw::math::StatisticsControl())
 Compute statistics of the distance field of a match list. More...
 
template<typename MatchT >
afw::math::Statistics makeMatchStatisticsInPixels (afw::geom::SkyWcs const &wcs, std::vector< MatchT > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl=afw::math::StatisticsControl())
 Compute statistics of on-detector radial separation for a match list, in pixels. More...
 
template<typename MatchT >
afw::math::Statistics makeMatchStatisticsInRadians (afw::geom::SkyWcs const &wcs, std::vector< MatchT > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl=afw::math::StatisticsControl())
 Compute statistics of on-sky radial separation for a match list, in radians. More...
 
ProxyVector makeProxies (afw::table::SourceCatalog const &sourceCat, afw::geom::SkyWcs const &distortedWcs, afw::geom::SkyWcs const &tanWcs)
 
ProxyVector makeProxies (afw::table::SimpleCatalog const &posRefCat, afw::geom::SkyWcs const &tanWcs)
 
afw::table::ReferenceMatchVector matchOptimisticB (afw::table::SimpleCatalog const &posRefCat, afw::table::SourceCatalog const &sourceCat, MatchOptimisticBControl const &control, afw::geom::SkyWcs const &wcs, int posRefBegInd=0, bool verbose=false)
 Match sources to stars in a position reference catalog using optimistic pattern matching B. More...
 
std::pair< size_t, size_t > find_candidate_reference_pair_range (float src_dist, ndarray::Array< float, 1, 1 > const &ref_dist_array, float max_dist_rad)
 Find the range of reference spokes within a spoke distance tolerance of our source spoke. More...
 
std::vector< std::pair< size_t, size_t > > create_pattern_spokes (ndarray::Array< double, 1, 1 > const &src_ctr, ndarray::Array< double, 2, 1 > const &src_delta_array, ndarray::Array< double, 1, 1 > const &src_dist_array, ndarray::Array< double, 1, 1 > const &ref_ctr, ndarray::Array< double, 1, 1 > const &proj_ref_ctr_delta, ndarray::Array< float, 1, 1 > const &ref_dist_array, ndarray::Array< uint16_t, 1, 1 > const &ref_id_array, ndarray::Array< double, 2, 1 > const &reference_array, double max_dist_rad, size_t n_match)
 Create the individual spokes that make up the pattern now that the shift and rotation are within tolerance. More...
 
int check_spoke (double cos_theta_src, double sin_theta_src, ndarray::Array< double, 1, 1 > const &ref_ctr, ndarray::Array< double, 1, 1 > const &proj_ref_ctr_delta, double proj_ref_ctr_dist_sq, std::pair< size_t, size_t > const &candidate_range, ndarray::Array< uint16_t, 1, 1 > const &ref_id_array, ndarray::Array< double, 2, 1 > const &reference_array, double src_sin_tol)
 Check the opening angle between the first spoke of our pattern for the source object against the reference object. More...
 
PolynomialTransform compose (geom::AffineTransform const &t1, PolynomialTransform const &t2)
 Return a PolynomialTransform that is equivalent to the composition t1(t2()) More...
 
PolynomialTransform compose (PolynomialTransform const &t1, geom::AffineTransform const &t2)
 Return a PolynomialTransform that is equivalent to the composition t1(t2()) More...
 
std::shared_ptr< afw::geom::SkyWcsmakeWcs (SipForwardTransform const &sipForward, SipReverseTransform const &sipReverse, geom::SpherePoint const &skyOrigin)
 Create a new TAN SIP Wcs from a pair of SIP transforms and the sky origin. More...
 
std::shared_ptr< afw::geom::SkyWcstransformWcsPixels (afw::geom::SkyWcs const &wcs, geom::AffineTransform const &s)
 Create a new SkyWcs whose pixel coordinate system has been transformed via an affine transform. More...
 
std::shared_ptr< afw::geom::SkyWcsrotateWcsPixelsBy90 (afw::geom::SkyWcs const &wcs, int nQuarter, geom::Extent2I const &dimensions)
 Return a new SkyWcs that represents a rotation of the image it corresponds to about the image's center. More...
 
 PYBIND11_MODULE (makeMatchStatistics, mod)
 
 PYBIND11_MODULE (matchOptimisticB, mod)
 
 PYBIND11_MODULE (pessimisticPatternMatcherUtils, mod)
 
 PYBIND11_MODULE (polynomialTransform, mod)
 
 PYBIND11_MODULE (scaledPolynomialTransformFitter, mod)
 
 PYBIND11_MODULE (sipTransform, mod)
 
template afw::math::Statistics makeMatchStatistics< afw::table::ReferenceMatch > (std::vector< afw::table::ReferenceMatch > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl)
 
template afw::math::Statistics makeMatchStatisticsInPixels< afw::table::ReferenceMatch > (afw::geom::SkyWcs const &wcs, std::vector< afw::table::ReferenceMatch > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl)
 
template afw::math::Statistics makeMatchStatisticsInRadians< afw::table::ReferenceMatch > (afw::geom::SkyWcs const &wcs, std::vector< afw::table::ReferenceMatch > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl)
 
template afw::math::Statistics makeMatchStatistics< afw::table::SourceMatch > (std::vector< afw::table::SourceMatch > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl)
 
template afw::math::Statistics makeMatchStatisticsInPixels< afw::table::SourceMatch > (afw::geom::SkyWcs const &wcs, std::vector< afw::table::SourceMatch > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl)
 
template afw::math::Statistics makeMatchStatisticsInRadians< afw::table::SourceMatch > (afw::geom::SkyWcs const &wcs, std::vector< afw::table::SourceMatch > const &matchList, int const flags, afw::math::StatisticsControl const &sctrl)
 

Typedef Documentation

◆ ProxyVector

Definition at line 52 of file matchOptimisticB.h.

Function Documentation

◆ check_spoke()

int lsst::meas::astrom::check_spoke ( double  cos_theta_src,
double  sin_theta_src,
ndarray::Array< double, 1, 1 > const &  ref_ctr,
ndarray::Array< double, 1, 1 > const &  proj_ref_ctr_delta,
double  proj_ref_ctr_dist_sq,
std::pair< size_t, size_t > const &  candidate_range,
ndarray::Array< uint16_t, 1, 1 > const &  ref_id_array,
ndarray::Array< double, 2, 1 > const &  reference_array,
double  src_sin_tol 
)

Check the opening angle between the first spoke of our pattern for the source object against the reference object.

This method makes heavy use of the small angle approximation to perform the comparison.

Parameters
[in]cos_theta_srcCosine of the angle between the current candidate source spoke and the first spoke.
[in]sin_theta_srcSine of the angle between the current candidate source spoke and the first spoke.
[in]ref_ctr3 vector of the candidate reference center
[in]proj_ref_ctr_deltaPlane projected first spoke in the reference pattern using the pattern center as normal.
[in]proj_ref_ctr_dist_sqSquared length of the projected vector.
[in]candidate_rangeMin and max index locations in ref_id_array that have pair lengths within the tolerance range.
[in]ref_id_arrayArray of id lookups into the master reference array that our center id object is paired with. uint16 type is to limit the amount of memory used and is set in the pessimistic_pattern_matcher_3d python class with reference catalogs trimmed by the matchPessimisticB runner class.
[in]reference_arrayArray of three vectors representing the locations of all reference objects.
[in]src_sin_tolSine of tolerance allowed between source and reference spoke opening angles.
Returns
ID of the candidate reference object successfully matched or -1 if no match is found.

Definition at line 137 of file pessimisticPatternMatcherUtils.cc.

141  {
142  // Loop over our candidate reference objects. candidate_range is the min
143  // and max of for pair candidates and are view into ref_id_array. Here we
144  // start from the midpoint of min and max values and step outward.
145  size_t midpoint = (candidate_range.first + candidate_range.second) / 2;
146  for (size_t idx = 0; idx < candidate_range.second - candidate_range.first; idx++) {
147  if (idx % 2 == 0) {
148  midpoint = midpoint + idx;
149  } else {
150  midpoint = midpoint - idx;
151  }
152  // Compute the delta vector from the pattern center.
153  uint16_t ref_id = ref_id_array[midpoint];
154  ndarray::Array<double, 1, 1> ref_delta = copy(reference_array[ref_id] - ref_ctr);
155 
156  double ref_dot = ndarray::asEigenMatrix(ref_delta).dot(ndarray::asEigenMatrix(ref_ctr));
157  ndarray::Array<double, 1, 1> proj_ref_delta = copy(ref_delta - ref_dot * ref_ctr);
158  // Compute the cos between our "center" reference vector and the
159  // current reference candidate.
160  auto proj_ref_delta_eigen = ndarray::asEigenMatrix(proj_ref_delta);
161  auto proj_ref_ctr_delta_eigen = ndarray::asEigenMatrix(proj_ref_ctr_delta);
162  double proj_delta_dist_sq = proj_ref_delta_eigen.dot(proj_ref_delta_eigen);
163  double geom_dist_ref = sqrt(proj_ref_ctr_dist_sq * proj_delta_dist_sq);
164  double cos_theta_ref = proj_ref_delta_eigen.dot(proj_ref_ctr_delta_eigen) / geom_dist_ref;
165 
166  // Make sure we can safely make the comparison in case
167  // our "center" and candidate vectors are mostly aligned.
168  double cos_sq_comparison;
169  if (cos_theta_ref * cos_theta_ref < (1 - src_sin_tol * src_sin_tol)) {
170  cos_sq_comparison = (cos_theta_src - cos_theta_ref) * (cos_theta_src - cos_theta_ref) /
171  (1 - cos_theta_ref * cos_theta_ref);
172  } else {
173  cos_sq_comparison = (cos_theta_src - cos_theta_ref) * (cos_theta_src - cos_theta_ref) /
174  (src_sin_tol * src_sin_tol);
175  }
176  // Test the difference of the cosine of the reference angle against
177  // the source angle. Assumes that the delta between the two is
178  // small.
179  if (cos_sq_comparison > src_sin_tol * src_sin_tol) {
180  continue;
181  }
182  // The cosine tests the magnitude of the angle but not
183  // its direction. To do that we need to know the sine as well.
184  // This cross product calculation does that.
185  Eigen::Vector3d cross_ref =
186  proj_ref_delta_eigen.head<3>().cross(proj_ref_ctr_delta_eigen.head<3>()) / geom_dist_ref;
187  double sin_theta_ref = cross_ref.dot(ndarray::asEigenMatrix(ref_ctr));
188  // Check the value of the cos again to make sure that it is not
189  // near zero.
190  double sin_comparison;
191  if (abs(cos_theta_src) < src_sin_tol) {
192  sin_comparison = (sin_theta_src - sin_theta_ref) / src_sin_tol;
193  } else {
194  sin_comparison = (sin_theta_src - sin_theta_ref) / cos_theta_ref;
195  }
196 
197  // Return the correct id of the candidate we found.
198  if (abs(sin_comparison) < src_sin_tol) {
199  return ref_id;
200  }
201  }
202  return -1;
203 }
T copy(T... args)
Angle abs(Angle const &a)
Definition: Angle.h:106
T sqrt(T... args)

◆ compose() [1/2]

PolynomialTransform lsst::meas::astrom::compose ( geom::AffineTransform const &  t1,
PolynomialTransform const &  t2 
)

Return a PolynomialTransform that is equivalent to the composition t1(t2())

The returned composition would be exact in ideal arithmetic, but may suffer from significant round-off error for high-order polynomials.

Definition at line 214 of file PolynomialTransform.cc.

214  {
215  typedef geom::AffineTransform AT;
216  PolynomialTransform result(t2.getOrder());
217 
218  result._xCoeffs.deep() = t2._xCoeffs * t1[AT::XX] + t2._yCoeffs * t1[AT::XY];
219  result._yCoeffs.deep() = t2._xCoeffs * t1[AT::YX] + t2._yCoeffs * t1[AT::YY];
220  result._xCoeffs(0, 0) += t1[AT::X];
221  result._yCoeffs(0, 0) += t1[AT::Y];
222  return result;
223 }
py::object result
Definition: _schema.cc:429
An affine coordinate transformation consisting of a linear transformation and an offset.

◆ compose() [2/2]

PolynomialTransform lsst::meas::astrom::compose ( PolynomialTransform const &  t1,
geom::AffineTransform const &  t2 
)

Return a PolynomialTransform that is equivalent to the composition t1(t2())

The returned composition would be exact in ideal arithmetic, but may suffer from significant round-off error for high-order polynomials.

Definition at line 225 of file PolynomialTransform.cc.

225  {
226  typedef geom::AffineTransform AT;
227  int const order = t1.getOrder();
228  if (order < 1) {
229  PolynomialTransform t1a(1);
230  t1a._xCoeffs(0, 0) = t1._xCoeffs(0, 0);
231  t1a._yCoeffs(0, 0) = t1._yCoeffs(0, 0);
232  return compose(t1a, t2);
233  }
234  detail::BinomialMatrix binomial(order);
235  // For each of these, (e.g.) a[n] == pow(a, n)
236  auto const t2u = detail::computePowers(t2[AT::X], order);
237  auto const t2v = detail::computePowers(t2[AT::Y], order);
238  auto const t2uu = detail::computePowers(t2[AT::XX], order);
239  auto const t2uv = detail::computePowers(t2[AT::XY], order);
240  auto const t2vu = detail::computePowers(t2[AT::YX], order);
241  auto const t2vv = detail::computePowers(t2[AT::YY], order);
242  PolynomialTransform result(order);
243  for (int p = 0; p <= order; ++p) {
244  for (int m = 0; m <= p; ++m) {
245  for (int j = 0; j <= m; ++j) {
246  for (int q = 0; p + q <= order; ++q) {
247  for (int n = 0; n <= q; ++n) {
248  for (int k = 0; k <= n; ++k) {
249  double z = binomial(p, m) * t2u[p - m] * binomial(m, j) * t2uu[j] * t2uv[m - j] *
250  binomial(q, n) * t2v[q - n] * binomial(n, k) * t2vu[k] * t2vv[n - k];
251  result._xCoeffs(j + k, m + n - j - k) += t1._xCoeffs(p, q) * z;
252  result._yCoeffs(j + k, m + n - j - k) += t1._yCoeffs(p, q) * z;
253  } // k
254  } // n
255  } // q
256  } // j
257  } // m
258  } // p
259  return result;
260 }
double z
Definition: Match.cc:44
int m
Definition: SpanSet.cc:48
table::Key< int > order

◆ create_pattern_spokes()

std::vector< std::pair< size_t, size_t > > lsst::meas::astrom::create_pattern_spokes ( ndarray::Array< double, 1, 1 > const &  src_ctr,
ndarray::Array< double, 2, 1 > const &  src_delta_array,
ndarray::Array< double, 1, 1 > const &  src_dist_array,
ndarray::Array< double, 1, 1 > const &  ref_ctr,
ndarray::Array< double, 1, 1 > const &  proj_ref_ctr_delta,
ndarray::Array< float, 1, 1 > const &  ref_dist_array,
ndarray::Array< uint16_t, 1, 1 > const &  ref_id_array,
ndarray::Array< double, 2, 1 > const &  reference_array,
double  max_dist_rad,
size_t  n_match 
)

Create the individual spokes that make up the pattern now that the shift and rotation are within tolerance.

If we can't create a valid pattern we exit early with a partial result.

Parameters
[in]src_ctr3 vector of the source pinwheel center
[in]src_delta_arrayArray of 3 vector deltas between the source center and the pairs that make up the remaining spokes of the pinwheel.
[in]src_dist_arrayArray of the distances of each src_delta in the pinwheel.
[in]ref_ctr3 vector of the candidate reference center.
[in]proj_ref_ctr_deltaPlane projected 3 vector formed from the center point of the candidate pin-wheel and the second point in the pattern to create the first spoke pair. This is the candidate pair that was matched in the main _construct_pattern_and_shift_rot_matrix loop.
[in]ref_dist_arrayArray of vector distances for each of the reference pairs
[in]ref_id_arrayArray of id lookups into the master reference array that our center id object is paired with.
[in]reference_arrayFull 3 vector data for the reference catalog.
[in]max_dist_radMaximum search radius for distances.
[in]n_matchNumber of source deltas that must be matched into the reference deltas in order to consider this a successful pattern match.
Returns
Return pairs of reference ids and their matched src ids.

Append the successful indices to our list. The src_idx needs an extra iteration to skip the first and second source objects.

Append the successful indices to our list. The src_idx needs an extra iteration to skip the first and second source objects.

Definition at line 47 of file pessimisticPatternMatcherUtils.cc.

52  {
53  // Struct where we will be putting our results.
55 
56  // Counter for number of spokes we failed to find a reference
57  // candidate for. We break the loop if we haven't found enough.
58  size_t n_fail = 0;
59 
60  // Plane project the center/first spoke of the source pattern using
61  // the center vector of the pattern as normal.
62  auto src_ctr_eigen = ndarray::asEigenMatrix(src_ctr);
63  double src_delta_ctr_dot = ndarray::asEigenMatrix(src_delta_array[0]).dot(src_ctr_eigen);
64 
65  ndarray::Array<double, 1, 1> proj_src_ctr_delta = copy(src_delta_array[0] - src_delta_ctr_dot * src_ctr);
66  auto proj_src_ctr_delta_eigen = ndarray::asEigenMatrix(proj_src_ctr_delta);
67  double proj_src_ctr_dist_sq = proj_src_ctr_delta_eigen.dot(proj_src_ctr_delta_eigen);
68 
69  // Pre - compute the squared length of the projected reference vector.
70  auto proj_ref_ctr_delta_eigen = ndarray::asEigenMatrix(proj_ref_ctr_delta);
71  double proj_ref_ctr_dist_sq = proj_ref_ctr_delta_eigen.dot(proj_ref_ctr_delta_eigen);
72 
73  // Loop over the source pairs.
74  // Value of sin where sin(theta) ~= theta to within 0.1%.
75  double max_sin_tol = 0.0447;
76  for (size_t src_idx = 1; src_idx < src_dist_array.size(); src_idx++) {
77  if (n_fail > src_dist_array.size() - (n_match - 1)) {
78  break;
79  }
80  // Find the reference pairs that include our candidate pattern center
81  // and sort them in increasing delta. Check this first so we don't
82  // compute anything else if no candidates exist.
83  std::pair<size_t, size_t> candidate_range =
84  find_candidate_reference_pair_range(src_dist_array[src_idx], ref_dist_array, max_dist_rad);
85  if (candidate_range.first == candidate_range.second) {
86  n_fail++;
87  continue;
88  }
89 
90  // Given our length tolerance we can use it to compute a tolerance
91  // on the angle between our spoke.
92  double src_sin_tol = max_dist_rad / (src_dist_array[src_idx] + max_dist_rad);
93 
94  // Test if the small angle approximation will still hold. This is
95  // defined as when sin(theta) ~= theta to within 0.1% of each
96  // other. If the implied opening angle is too large we set it to
97  // the 0.1% threshold.
98  if (src_sin_tol > max_sin_tol) {
99  src_sin_tol = max_sin_tol;
100  }
101 
102  // Plane project the candidate source spoke and compute the cosine
103  // and sine of the opening angle.
104  double proj_src_delta_dot = ndarray::asEigenMatrix(src_delta_array[src_idx]).dot(src_ctr_eigen);
105 
106  ndarray::Array<double, 1, 1> proj_src_delta =
107  copy(src_delta_array[src_idx] - proj_src_delta_dot * src_ctr);
108  auto proj_src_delta_eigen = ndarray::asEigenMatrix(proj_src_delta);
109  double geom_dist_src = sqrt(proj_src_delta_eigen.dot(proj_src_delta_eigen) * proj_src_ctr_dist_sq);
110 
111  // Compute cosine and sine of the delta vector opening angle.
112  double cos_theta_src = proj_src_delta_eigen.dot(proj_src_ctr_delta_eigen) / geom_dist_src;
113  Eigen::Vector3d cross_src =
114  proj_src_delta_eigen.head<3>().cross(proj_src_ctr_delta_eigen.head<3>()) / geom_dist_src;
115  double sin_theta_src = cross_src.dot(src_ctr_eigen);
116 
117  // Test the spokes and return the id of the reference object.
118  // Return -1 if no match is found.
119  int ref_id =
120  check_spoke(cos_theta_src, sin_theta_src, ref_ctr, proj_ref_ctr_delta, proj_ref_ctr_dist_sq,
121  candidate_range, ref_id_array, reference_array, src_sin_tol);
122  if (ref_id < 0) {
123  n_fail++;
124  continue;
125  }
126 
129  output_spokes.push_back(std::make_pair(static_cast<size_t>(ref_id), src_idx + 1));
130  if (output_spokes.size() >= n_match - 2) {
131  break;
132  }
133  }
134  return output_spokes;
135 }
T make_pair(T... args)
int check_spoke(double cos_theta_src, double sin_theta_src, ndarray::Array< double, 1, 1 > const &ref_ctr, ndarray::Array< double, 1, 1 > const &proj_ref_ctr_delta, double proj_ref_ctr_dist_sq, std::pair< size_t, size_t > const &candidate_range, ndarray::Array< uint16_t, 1, 1 > const &ref_id_array, ndarray::Array< double, 2, 1 > const &reference_array, double src_sin_tol)
Check the opening angle between the first spoke of our pattern for the source object against the refe...
std::pair< size_t, size_t > find_candidate_reference_pair_range(float src_dist, ndarray::Array< float, 1, 1 > const &ref_dist_array, float max_dist_rad)
Find the range of reference spokes within a spoke distance tolerance of our source spoke.
T push_back(T... args)
T size(T... args)

◆ find_candidate_reference_pair_range()

std::pair< size_t, size_t > lsst::meas::astrom::find_candidate_reference_pair_range ( float  src_dist,
ndarray::Array< float, 1, 1 > const &  ref_dist_array,
float  max_dist_rad 
)

Find the range of reference spokes within a spoke distance tolerance of our source spoke.

Returns an the min and max index spanning src_dist +/- max_dist_rad.

Parameters
[in]src_distfloat value of the distance we would like to search for in the reference array in radians.
[in]ref_dist_arraysorted array of distances in radians.
[in]max_dist_radmaximum plus/minus search to find in the reference array in radians.
Returns
pair of indices for the min and max range of indices to search.

Definition at line 35 of file pessimisticPatternMatcherUtils.cc.

36  {
37  auto itr =
38  std::upper_bound(ref_dist_array.begin(), ref_dist_array.end(), src_dist - max_dist_rad - 1e-16);
39  auto itrEnd =
40  std::upper_bound(ref_dist_array.begin(), ref_dist_array.end(), src_dist + max_dist_rad + 1e-16);
41 
42  size_t startIdx = itr - ref_dist_array.begin();
43  size_t endIdx = itrEnd - ref_dist_array.begin();
44  return std::make_pair(startIdx, endIdx);
45 }
T upper_bound(T... args)

◆ makeMatchStatistics()

template<typename MatchT >
afw::math::Statistics lsst::meas::astrom::makeMatchStatistics ( std::vector< MatchT > const &  matchList,
int const  flags,
afw::math::StatisticsControl const &  sctrl = afw::math::StatisticsControl() 
)

Compute statistics of the distance field of a match list.

Parameters
[in]matchListlist of matchList between reference objects and sources; fields read:
  • distance: distance between source and reference object, in arbitrary units; the resulting statistics have the same units as distance
[in]flagswhat to calculate; OR constants such as lsst::afw::math::MEAN, MEANCLIP, STDDEV, MEDIAN, defined in lsst/afw/math/Statitics.h's Property enum
[in]sctrlstatistics configuration

Definition at line 33 of file makeMatchStatistics.cc.

34  {
35  if (matchList.empty()) {
36  throw LSST_EXCEPT(pexExcept::RuntimeError, "matchList is empty");
37  }
39  val.reserve(matchList.size());
40 
41  for (auto const& match : matchList) {
42  val.push_back(match.distance);
43  }
44  return afw::math::makeStatistics(val, flags, sctrl);
45 }
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition: Exception.h:48
Reports errors that are due to events beyond the control of the program.
Definition: Runtime.h:104
T empty(T... args)
Statistics makeStatistics(lsst::afw::image::Image< Pixel > const &img, lsst::afw::image::Mask< image::MaskPixel > const &msk, int const flags, StatisticsControl const &sctrl=StatisticsControl())
Handle a watered-down front-end to the constructor (no variance)
Definition: Statistics.h:359
ImageT val
Definition: CR.cc:146

◆ makeMatchStatistics< afw::table::ReferenceMatch >()

◆ makeMatchStatistics< afw::table::SourceMatch >()

◆ makeMatchStatisticsInPixels()

template<typename MatchT >
afw::math::Statistics lsst::meas::astrom::makeMatchStatisticsInPixels ( afw::geom::SkyWcs const &  wcs,
std::vector< MatchT > const &  matchList,
int const  flags,
afw::math::StatisticsControl const &  sctrl = afw::math::StatisticsControl() 
)

Compute statistics of on-detector radial separation for a match list, in pixels.

Parameters
[in]wcsWCS describing pixel to sky transformation
[in]matchListlist of matchList between reference objects and sources; fields read:
  • first: reference object; only the coord is read
  • second: source; only the centroid is read
[in]flagswhat to calculate; OR constants such as lsst::afw::math::MEAN, MEANCLIP, STDDEV, MEDIAN, defined in lsst/afw/math/Statitics.h's Property enum
[in]sctrlstatistics configuration

Definition at line 48 of file makeMatchStatistics.cc.

50  {
51  if (matchList.empty()) {
52  throw LSST_EXCEPT(pexExcept::RuntimeError, "matchList is empty");
53  }
55  val.reserve(matchList.size());
56 
57  for (auto const& match : matchList) {
58  auto refPtr = match.first;
59  auto srcPtr = match.second;
60  auto srcX = srcPtr->getX();
61  auto srcY = srcPtr->getY();
62  auto refPos = wcs.skyToPixel(refPtr->getCoord());
63  auto refX = refPos[0];
64  auto refY = refPos[1];
65  val.push_back(::hypot(srcX - refX, srcY - refY));
66  }
67  return afw::math::makeStatistics(val, flags, sctrl);
68 }
table::Key< table::Array< std::uint8_t > > wcs
Definition: SkyWcs.cc:66
T hypot(T... args)

◆ makeMatchStatisticsInPixels< afw::table::ReferenceMatch >()

◆ makeMatchStatisticsInPixels< afw::table::SourceMatch >()

◆ makeMatchStatisticsInRadians()

template<typename MatchT >
afw::math::Statistics lsst::meas::astrom::makeMatchStatisticsInRadians ( afw::geom::SkyWcs const &  wcs,
std::vector< MatchT > const &  matchList,
int const  flags,
afw::math::StatisticsControl const &  sctrl = afw::math::StatisticsControl() 
)

Compute statistics of on-sky radial separation for a match list, in radians.

Parameters
[in]wcsWCS describing pixel to sky transformation
[in]matchListlist of matchList between reference objects and sources; fields read:
  • first: reference object; only the coord is read
  • second: source; only the centroid is read
[in]flagswhat to calculate; OR constants such as lsst::afw::math::MEAN, MEANCLIP, STDDEV, MEDIAN, defined in lsst/afw/math/Statitics.h's Property enum
[in]sctrlstatistics configuration

Definition at line 71 of file makeMatchStatistics.cc.

73  {
74  if (matchList.empty()) {
75  throw LSST_EXCEPT(pexExcept::RuntimeError, "matchList is empty");
76  }
78  val.reserve(matchList.size());
79 
80  for (auto const& match : matchList) {
81  auto refPtr = match.first;
82  auto srcPtr = match.second;
83  auto refCoord = refPtr->getCoord();
84  auto srcCoord = wcs.pixelToSky(srcPtr->getCentroid());
85  auto angSep = refCoord.separation(srcCoord);
86  val.push_back(angSep.asRadians());
87  }
88  return afw::math::makeStatistics(val, flags, sctrl);
89 }

◆ makeMatchStatisticsInRadians< afw::table::ReferenceMatch >()

◆ makeMatchStatisticsInRadians< afw::table::SourceMatch >()

◆ makeProxies() [1/2]

ProxyVector lsst::meas::astrom::makeProxies ( afw::table::SimpleCatalog const &  posRefCat,
afw::geom::SkyWcs const &  tanWcs 
)

Definition at line 448 of file matchOptimisticB.cc.

448  {
449  auto coordKey = afwTable::CoordKey(posRefCat.getSchema()["coord"]);
450  ProxyVector r;
451  r.reserve(posRefCat.size());
452  for (afwTable::SimpleCatalog::const_iterator posRefPtr = posRefCat.begin(); posRefPtr != posRefCat.end();
453  ++posRefPtr) {
454  r.push_back(RecordProxy(posRefPtr, tanWcs.skyToPixel(posRefPtr->get(coordKey))));
455  }
456  return r;
457 }
A FunctorKey used to get or set celestial coordinates from a pair of lsst::geom::Angle keys.
Definition: aggregates.h:210
typename Base::const_iterator const_iterator
Definition: SortedCatalog.h:50
T reserve(T... args)
A wrapper around a SimpleRecord or SourceRecord that allows us to record a pixel position in a way th...

◆ makeProxies() [2/2]

ProxyVector lsst::meas::astrom::makeProxies ( afw::table::SourceCatalog const &  sourceCat,
afw::geom::SkyWcs const &  distortedWcs,
afw::geom::SkyWcs const &  tanWcs 
)

Definition at line 436 of file matchOptimisticB.cc.

437  {
438  ProxyVector r;
439  r.reserve(sourceCat.size());
440  for (afwTable::SourceCatalog::const_iterator sourcePtr = sourceCat.begin(); sourcePtr != sourceCat.end();
441  ++sourcePtr) {
442  r.push_back(
443  RecordProxy(sourcePtr, tanWcs.skyToPixel(distortedWcs.pixelToSky(sourcePtr->getCentroid()))));
444  }
445  return r;
446 }

◆ makeWcs()

std::shared_ptr< afw::geom::SkyWcs > lsst::meas::astrom::makeWcs ( SipForwardTransform const &  sipForward,
SipReverseTransform const &  sipReverse,
geom::SpherePoint const &  skyOrigin 
)

Create a new TAN SIP Wcs from a pair of SIP transforms and the sky origin.

Parameters
[in]sipForwardMapping from pixel coordinates to intermediate world coordinates.
[in]sipReverseMapping from intermediate world coordinates to pixel coordinates.
[in]skyOriginICRS position of the gnomonic projection that maps sky coordinates to intermediate world coordinates (CRVAL).
Exceptions
pex::exceptions::InvalidParameterErrorif the forward and reverse SIP transforms have different CRPIX values or CD matrices.

Definition at line 148 of file SipTransform.cc.

150  {
151  if (!sipForward.getPixelOrigin().asEigen().isApprox(sipReverse.getPixelOrigin().asEigen())) {
152  std::ostringstream oss;
153  oss << "SIP forward and reverse transforms have inconsistent CRPIX: " << sipForward.getPixelOrigin()
154  << " != " << sipReverse.getPixelOrigin();
155  throw LSST_EXCEPT(pex::exceptions::InvalidParameterError, oss.str());
156  }
157  if (!sipForward.getCdMatrix().getMatrix().isApprox(sipReverse.getCdMatrix().getMatrix())) {
158  std::ostringstream oss;
159  oss << "SIP forward and reverse transforms have inconsistent CD matrix: " << sipForward.getCdMatrix()
160  << "\n!=\n"
161  << sipReverse.getCdMatrix();
162  throw LSST_EXCEPT(pex::exceptions::InvalidParameterError, oss.str());
163  }
164  Eigen::MatrixXd sipA(ndarray::asEigenMatrix(sipForward.getPoly().getXCoeffs()));
165  Eigen::MatrixXd sipB(ndarray::asEigenMatrix(sipForward.getPoly().getYCoeffs()));
166  Eigen::MatrixXd sipAP(ndarray::asEigenMatrix(sipReverse.getPoly().getXCoeffs()));
167  Eigen::MatrixXd sipBP(ndarray::asEigenMatrix(sipReverse.getPoly().getYCoeffs()));
168 
169  return afw::geom::makeTanSipWcs(sipForward.getPixelOrigin(), skyOrigin,
170  sipForward.getCdMatrix().getMatrix(), sipA, sipB, sipAP, sipBP);
171 }
std::shared_ptr< SkyWcs > makeTanSipWcs(lsst::geom::Point2D const &crpix, lsst::geom::SpherePoint const &crval, Eigen::Matrix2d const &cdMatrix, Eigen::MatrixXd const &sipA, Eigen::MatrixXd const &sipB)
Construct a TAN-SIP SkyWcs with forward SIP distortion terms and an iterative inverse.
Definition: SkyWcs.cc:538
T str(T... args)

◆ matchOptimisticB()

afwTable::ReferenceMatchVector lsst::meas::astrom::matchOptimisticB ( afw::table::SimpleCatalog const &  posRefCat,
afw::table::SourceCatalog const &  sourceCat,
MatchOptimisticBControl const &  control,
afw::geom::SkyWcs const &  wcs,
int  posRefBegInd = 0,
bool  verbose = false 
)

Match sources to stars in a position reference catalog using optimistic pattern matching B.

Optimistic Pattern Matching is described in V. Tabur 2007, PASA, 24, 189-198 "Fast Algorithms for Matching CCD Images to a Stellar Catalogue"

Parameters
[in]posRefCatcatalog of position reference stars; fields read:
  • "coord"
  • control.refFluxField
[in]sourceCatcatalog of detected sources; fields read:
  • "Centroid_x"
  • "Centroid_y"
  • control.refFluxField
[in]wcsestimated WCS
[in]controlcontrol object
[in]posRefBegIndindex of first start to use in posRefCat
[in]verbosetrue to print diagnostic information to std::cout
Returns
a list of matches; the d field may be set, but you should not rely on it

Definition at line 459 of file matchOptimisticB.cc.

463  {
464  control.validate();
465  if (posRefCat.empty()) {
466  throw LSST_EXCEPT(pexExcept::InvalidParameterError, "no entries in posRefCat");
467  }
468  if (sourceCat.empty()) {
469  throw LSST_EXCEPT(pexExcept::InvalidParameterError, "no entries in sourceCat");
470  }
471  if (posRefBegInd < 0) {
472  throw LSST_EXCEPT(pexExcept::InvalidParameterError, "posRefBegInd < 0");
473  }
474  if (static_cast<size_t>(posRefBegInd) >= posRefCat.size()) {
475  throw LSST_EXCEPT(pexExcept::InvalidParameterError, "posRefBegInd too big");
476  }
477  double const maxRotationRad = geom::degToRad(control.maxRotationDeg);
478 
479  // Create an undistorted Wcs to project everything with
480  // We'll anchor it at the center of the area.
481  geom::Extent2D srcCenter(0, 0);
482  for (auto iter = sourceCat.begin(); iter != sourceCat.end(); ++iter) {
483  srcCenter += geom::Extent2D(iter->getCentroid());
484  }
485  srcCenter /= sourceCat.size();
486 
487  sphgeom::Vector3d refCenter(0, 0, 0);
488  for (auto iter = posRefCat.begin(); iter != posRefCat.end(); ++iter) {
489  refCenter += iter->getCoord().getVector();
490  }
491  refCenter /= posRefCat.size();
492 
493  auto tanWcs =
494  afw::geom::makeSkyWcs(geom::Point2D(srcCenter), geom::SpherePoint(refCenter), wcs.getCdMatrix());
495 
496  ProxyVector posRefProxyCat = makeProxies(posRefCat, *tanWcs);
497  ProxyVector sourceProxyCat = makeProxies(sourceCat, wcs, *tanWcs);
498 
499  // sourceSubCat contains at most the numBrightStars brightest sources, sorted by decreasing flux
500  ProxyVector sourceSubCat =
501  selectPoint(sourceProxyCat, sourceCat.getSchema().find<double>(control.sourceFluxField).key,
502  control.numBrightStars);
503 
504  // posRefSubCat skips the initial posRefBegInd brightest reference objects and contains
505  // at most the next len(sourceSubCat) + 25 brightest reference objects, sorted by decreasing flux
506  ProxyVector posRefSubCat =
507  selectPoint(posRefProxyCat, posRefCat.getSchema().find<double>(control.refFluxField).key,
508  sourceSubCat.size() + 25, posRefBegInd);
509  if (verbose) {
510  std::cout << "Catalog sizes: " << sourceSubCat.size() << " " << posRefSubCat.size() << std::endl;
511  }
512 
513  // Construct a list of pairs of position reference stars sorted by increasing separation
514  std::vector<ProxyPair> posRefPairList;
515  size_t const posRefCatSubSize = posRefSubCat.size();
516  for (size_t i = 0; i < posRefCatSubSize - 1; i++) {
517  for (size_t j = i + 1; j < posRefCatSubSize; j++) {
518  posRefPairList.push_back(ProxyPair(posRefSubCat[i], posRefSubCat[j]));
519  }
520  }
521  std::sort(posRefPairList.begin(), posRefPairList.end(), cmpPair);
522 
523  // Construct a list of pairs of sources sorted by increasing separation
524  std::vector<ProxyPair> sourcePairList;
525  size_t const sourceSubCatSize = sourceSubCat.size();
526  for (size_t i = 0; i < sourceSubCatSize - 1; i++) {
527  for (size_t j = i + 1; j < sourceSubCatSize; j++) {
528  sourcePairList.push_back(ProxyPair(sourceSubCat[i], sourceSubCat[j]));
529  }
530  }
531  std::sort(sourcePairList.begin(), sourcePairList.end(), cmpPair);
532 
534  afwTable::ReferenceMatchVector matPairSave;
536 
537  size_t const fullShapeSize = control.numPointsForShape - 1; // Max size of shape array
538  for (size_t ii = 0; ii < sourcePairList.size(); ii++) {
539  ProxyPair p = sourcePairList[ii];
540 
542  searchPair(posRefPairList, p, control.matchingAllowancePix, maxRotationRad);
543 
544  // If candidate pairs are found
545  if (q.size() != 0) {
546  std::vector<ProxyPair> srcMatPair;
547  std::vector<ProxyPair> catMatPair;
548 
549  // Go through candidate pairs
550  for (size_t l = 0; l < q.size(); l++) {
551  // sign matters, so don't use deltaAng; no need to wrap because
552  // the result is used with deltaAng later
553  double dpa = p.pa - q[l].pa;
554 
555  srcMatPair.clear();
556  catMatPair.clear();
557 
558  srcMatPair.push_back(p);
559  catMatPair.push_back(q[l]);
560 
561  if (verbose) {
562  std::cout << "p dist: " << p.distance << " pa: " << p.pa << std::endl;
563  std::cout << "q dist: " << q[l].distance << " pa: " << q[l].pa << std::endl;
564  }
565 
566  for (size_t k = 0; k < sourceSubCat.size(); k++) {
567  if (p.first == sourceSubCat[k] || p.second == sourceSubCat[k]) continue;
568 
569  ProxyPair pp(p.first, sourceSubCat[k]);
570 
571  std::vector<ProxyPair>::iterator r = searchPair3(
572  posRefPairList, pp, q[l], control.matchingAllowancePix, dpa, maxRotationRad);
573  if (r != posRefPairList.end()) {
574  srcMatPair.push_back(pp);
575  catMatPair.push_back(*r);
576  if (verbose) {
577  std::cout << " p dist: " << pp.distance << " pa: " << pp.pa << std::endl;
578  std::cout << " r dist: " << (*r).distance << " pa: " << (*r).pa << std::endl;
579  }
580  if (srcMatPair.size() == fullShapeSize) {
581  break;
582  }
583  }
584  }
585 
586  bool goodMatch = false;
587  if (srcMatPair.size() == fullShapeSize) {
588  goodMatch = true;
589  for (size_t k = 1; k < catMatPair.size(); k++) {
590  if (catMatPair[0].first != catMatPair[k].first) {
591  goodMatch = false;
592  break;
593  }
594  }
595  }
596 
597  if (goodMatch && srcMatPair.size() == fullShapeSize) {
598  ProxyVector srcMat;
599  ProxyVector catMat;
600 
601  srcMat.push_back(srcMatPair[0].first);
602  catMat.push_back(catMatPair[0].first);
603  for (size_t k = 0; k < srcMatPair.size(); k++) {
604  srcMat.push_back(srcMatPair[k].second);
605  catMat.push_back(catMatPair[k].second);
606  }
607 
608  boost::shared_array<double> coeff = polyfit(1, srcMat, catMat);
609 
610  if (verbose) {
611  for (size_t k = 0; k < srcMat.size(); k++) {
612  std::cout << "circle(" << srcMat[k].getX() << "," << srcMat[k].getY()
613  << ",10) # color=green" << std::endl;
614  std::cout << "circle(" << catMat[k].getX() << "," << catMat[k].getY()
615  << ",10) # color=red" << std::endl;
616  std::cout << "line(" << srcMat[0].getX() << "," << srcMat[0].getY() << ","
617  << srcMat[k].getX() << "," << srcMat[k].getY()
618  << ") # line=0 0 color=green" << std::endl;
619  std::cout << "line(" << catMat[0].getX() << "," << catMat[0].getY() << ","
620  << catMat[k].getX() << "," << catMat[k].getY()
621  << ") # line=0 0 color=red" << std::endl;
622  }
623  }
624 
625  double a = coeff[1];
626  double b = coeff[2];
627  double c = coeff[4];
628  double d = coeff[5];
629  geom::Angle const theta(std::acos((a * b + c * d) /
630  (std::sqrt(a * a + c * c) * std::sqrt(b * b + d * d))),
631  geom::radians);
632  if (verbose) {
633  std::cout << "Linear fit from match:" << std::endl;
634  std::cout << coeff[0] << " " << coeff[1] << " " << coeff[2] << std::endl;
635  std::cout << coeff[3] << " " << coeff[4] << " " << coeff[5] << std::endl;
636  std::cout << "Determinant (max " << control.maxDeterminant << "): ";
637  std::cout << coeff[1] * coeff[5] - coeff[2] * coeff[4] - 1. << std::endl;
638  std::cout << "Angle between axes (deg; allowed 90 +/- ";
639  std::cout << control.allowedNonperpDeg << "): ";
640  std::cout << theta.asDegrees() << std::endl;
641  }
642  if (std::fabs(coeff[1] * coeff[5] - coeff[2] * coeff[4] - 1.) > control.maxDeterminant ||
643  std::fabs(theta.asDegrees() - 90) > control.allowedNonperpDeg ||
644  std::fabs(coeff[0]) > control.maxOffsetPix ||
645  std::fabs(coeff[3]) > control.maxOffsetPix) {
646  if (verbose) {
647  std::cout << "Bad; continuing" << std::endl;
648  }
649  continue;
650  } else {
651  double x0, y0, x1, y1;
652  int num = 0;
653  srcMat.clear();
654  catMat.clear();
655  for (size_t i = 0; i < sourceSubCat.size(); i++) {
656  x0 = sourceSubCat[i].getX();
657  y0 = sourceSubCat[i].getY();
658  transform(1, coeff, x0, y0, &x1, &y1);
659  auto refObjDist =
660  searchNearestPoint(posRefSubCat, x1, y1, control.matchingAllowancePix);
661  if (refObjDist.first != posRefSubCat.end()) {
662  num++;
663  srcMat.push_back(sourceSubCat[i]);
664  catMat.push_back(*refObjDist.first);
665  if (verbose) {
666  std::cout << "Match: " << x0 << "," << y0 << " --> " << x1 << "," << y1
667  << " <==> " << refObjDist.first->getX() << ","
668  << refObjDist.first->getY() << std::endl;
669  }
670  }
671  }
672  if (num <= control.numPointsForShape) {
673  // Can get matrix = 0,0,0,0; everything matches a single catalog object
674  if (verbose) {
675  std::cout << "Insufficient initial matches; continuing" << std::endl;
676  }
677  continue;
678  }
679  coeff = polyfit(1, srcMat, catMat);
680  if (verbose) {
681  std::cout << "Coefficients from initial matching:" << std::endl;
682  for (size_t i = 0; i < 6; ++i) {
683  std::cout << coeff[i] << " ";
684  }
685  std::cout << std::endl;
686  }
687 
688  matPair = FinalVerify(coeff, posRefProxyCat, sourceProxyCat,
689  control.matchingAllowancePix, verbose);
690  if (verbose) {
691  std::cout << "Number of matches: " << matPair.size() << " vs "
692  << control.minMatchedPairs << std::endl;
693  }
694  if (matPair.size() <= static_cast<std::size_t>(control.minMatchedPairs)) {
695  if (verbose) {
696  std::cout << "Insufficient final matches; continuing" << std::endl;
697  }
698  if (matPair.size() > matPairSave.size()) {
699  matPairSave = matPair;
700  }
701  continue;
702  } else {
703  if (verbose) {
704  std::cout << "Finish" << std::endl;
705  }
706  matPairCand.push_back(matPair);
707  if (matPairCand.size() == 3) {
708  goto END;
709  }
710  }
711  }
712  }
713  }
714  }
715  }
716 
717 END:
718  if (matPairCand.size() == 0) {
719  return matPairSave;
720  } else {
721  size_t nmatch = matPairCand[0].size();
722  afwTable::ReferenceMatchVector matPairRet = matPairCand[0];
723  for (size_t i = 1; i < matPairCand.size(); i++) {
724  if (matPairCand[i].size() > nmatch) {
725  nmatch = matPairCand[i].size();
726  matPairRet = matPairCand[i];
727  }
728  }
729  return matPairRet;
730  }
731 }
table::Key< int > transform
table::Key< int > b
table::Key< int > a
T acos(T... args)
T begin(T... args)
A class representing an angle.
Definition: Angle.h:127
Point in an unspecified spherical coordinate system.
Definition: SpherePoint.h:57
Reports invalid arguments.
Definition: Runtime.h:66
T clear(T... args)
T end(T... args)
T endl(T... args)
T fabs(T... args)
std::shared_ptr< SkyWcs > makeSkyWcs(daf::base::PropertySet &metadata, bool strip=false)
Construct a SkyWcs from FITS keywords.
Definition: SkyWcs.cc:521
Extent< double, 2 > Extent2D
Definition: Extent.h:400
constexpr double degToRad(double x) noexcept
Definition: Angle.h:51
constexpr AngleUnit radians
constant with units of radians
Definition: Angle.h:108
ProxyVector makeProxies(afw::table::SourceCatalog const &sourceCat, afw::geom::SkyWcs const &distortedWcs, afw::geom::SkyWcs const &tanWcs)
T sort(T... args)
table::Key< table::Array< double > > coeff
Definition: PsfexPsf.cc:362

◆ PYBIND11_MODULE() [1/6]

lsst::meas::astrom::PYBIND11_MODULE ( makeMatchStatistics  ,
mod   
)

Definition at line 47 of file makeMatchStatistics.cc.

47  {
48  py::module::import("lsst.afw.math");
49 
50 
51  declareMakeMatchStatistics<afw::table::ReferenceMatch>(mod);
52  declareMakeMatchStatistics<afw::table::SourceMatch>(mod);
53 }

◆ PYBIND11_MODULE() [2/6]

lsst::meas::astrom::PYBIND11_MODULE ( matchOptimisticB  ,
mod   
)

Definition at line 91 of file matchOptimisticB.cc.

91  {
92  declareRecordProxy(mod);
93  declareProxyPair(mod);
94  declareMatchOptimisticBControl(mod);
95 
96  mod.def("makeProxies",
97  (ProxyVector(*)(afw::table::SourceCatalog const &, afw::geom::SkyWcs const &,
98  afw::geom::SkyWcs const &)) &
100  "sourceCat"_a, "distortedWcs"_a, "tanWcs"_a);
101  mod.def("makeProxies",
102  (ProxyVector(*)(afw::table::SimpleCatalog const &, afw::geom::SkyWcs const &)) & makeProxies,
103  "posRefCat"_a, "tanWcs"_a);
104 
105  mod.def("matchOptimisticB", &matchOptimisticB, "posRefCat"_a, "sourceCat"_a, "control"_a, "wcs"_a,
106  "posRefBegInd"_a = 0, "verbose"_a = false);
107 }
SortedCatalogT< SourceRecord > SourceCatalog
Definition: fwd.h:85
SortedCatalogT< SimpleRecord > SimpleCatalog
Definition: fwd.h:79
std::vector< RecordProxy > ProxyVector
afw::table::ReferenceMatchVector matchOptimisticB(afw::table::SimpleCatalog const &posRefCat, afw::table::SourceCatalog const &sourceCat, MatchOptimisticBControl const &control, afw::geom::SkyWcs const &wcs, int posRefBegInd=0, bool verbose=false)
Match sources to stars in a position reference catalog using optimistic pattern matching B.

◆ PYBIND11_MODULE() [3/6]

lsst::meas::astrom::PYBIND11_MODULE ( pessimisticPatternMatcherUtils  ,
mod   
)

Definition at line 40 of file pessimisticPatternMatcherUtils.cc.

40  {
41  mod.def("find_candidate_reference_pair_range", &find_candidate_reference_pair_range, "src_dist"_a,
42  "ref_dist_array"_a, "max_dist_rad"_a);
43  mod.def("create_pattern_spokes", &create_pattern_spokes, "src_ctr"_a, "src_delta_array"_a,
44  "src_dist_array"_a, "ref_ctr"_a, "proj_ref_ctr_delta"_a, "ref_dist_array"_a, "ref_id_array"_a,
45  "reference_array"_a, "max_dist_rad"_a, "n_match"_a);
46  mod.def("check_spoke", &check_spoke, "cos_theta_src"_a, "sin_theta_src"_a, "ref_ctr"_a,
47  "proj_ref_ctr_delta"_a, "proj_ref_ctr_dist_sq"_a, "candidate_range"_a, "ref_id_array"_a,
48  "reference_array"_a, "src_sin_tol"_a);
49 }
std::vector< std::pair< size_t, size_t > > create_pattern_spokes(ndarray::Array< double, 1, 1 > const &src_ctr, ndarray::Array< double, 2, 1 > const &src_delta_array, ndarray::Array< double, 1, 1 > const &src_dist_array, ndarray::Array< double, 1, 1 > const &ref_ctr, ndarray::Array< double, 1, 1 > const &proj_ref_ctr_delta, ndarray::Array< float, 1, 1 > const &ref_dist_array, ndarray::Array< uint16_t, 1, 1 > const &ref_id_array, ndarray::Array< double, 2, 1 > const &reference_array, double max_dist_rad, size_t n_match)
Create the individual spokes that make up the pattern now that the shift and rotation are within tole...

◆ PYBIND11_MODULE() [4/6]

lsst::meas::astrom::PYBIND11_MODULE ( polynomialTransform  ,
mod   
)

Definition at line 102 of file polynomialTransform.cc.

102  {
103  declarePolynomialTransform(mod);
104  declareScaledPolynomialTransform(mod);
105 
106  mod.def("compose",
107  (PolynomialTransform(*)(geom::AffineTransform const &, PolynomialTransform const &)) & compose,
108  "t1"_a, "t2"_a);
109  mod.def("compose",
110  (PolynomialTransform(*)(PolynomialTransform const &, geom::AffineTransform const &)) & compose,
111  "t1"_a, "t2"_a);
112 }

◆ PYBIND11_MODULE() [5/6]

lsst::meas::astrom::PYBIND11_MODULE ( scaledPolynomialTransformFitter  ,
mod   
)

Definition at line 71 of file scaledPolynomialTransformFitter.cc.

71  {
72  declareOutlierRejectionControl(mod);
73  declareScaledPolynomialTransformFitter(mod);
74 }

◆ PYBIND11_MODULE() [6/6]

lsst::meas::astrom::PYBIND11_MODULE ( sipTransform  ,
mod   
)

Definition at line 104 of file sipTransform.cc.

104  {
105  declareSipTransformBase(mod);
106  declareSipForwardTransform(mod);
107  declareSipReverseTransform(mod);
108 
109  mod.def("makeWcs", makeWcs, "sipForward"_a, "sipReverse"_a, "skyOrigin"_a);
110  mod.def("transformWcsPixels", transformWcsPixels, "wcs"_a, "s"_a);
111  mod.def("rotateWcsPixelsBy90", rotateWcsPixelsBy90, "wcs"_a, "nQuarter"_a, "dimensions"_a);
112 }
std::shared_ptr< afw::geom::SkyWcs > transformWcsPixels(afw::geom::SkyWcs const &wcs, geom::AffineTransform const &s)
Create a new SkyWcs whose pixel coordinate system has been transformed via an affine transform.
std::shared_ptr< afw::geom::SkyWcs > rotateWcsPixelsBy90(afw::geom::SkyWcs const &wcs, int nQuarter, geom::Extent2I const &dimensions)
Return a new SkyWcs that represents a rotation of the image it corresponds to about the image's cente...
std::shared_ptr< afw::geom::SkyWcs > makeWcs(SipForwardTransform const &sipForward, SipReverseTransform const &sipReverse, geom::SpherePoint const &skyOrigin)
Create a new TAN SIP Wcs from a pair of SIP transforms and the sky origin.

◆ rotateWcsPixelsBy90()

std::shared_ptr< afw::geom::SkyWcs > lsst::meas::astrom::rotateWcsPixelsBy90 ( afw::geom::SkyWcs const &  wcs,
int  nQuarter,
geom::Extent2I const &  dimensions 
)

Return a new SkyWcs that represents a rotation of the image it corresponds to about the image's center.

Parameters
[in]wcsOriginal SkyWcs to be rotated.
[in]nQuarterNumber of 90 degree rotations (positive is counterclockwise).
[in]dimensionsWidth and height of the image.

Definition at line 179 of file SipTransform.cc.

180  {
181  geom::Extent2D offset;
182  switch (nQuarter % 4) {
183  case 0:
184  offset = geom::Extent2D(0, 0);
185  break;
186  case 1:
187  offset = geom::Extent2D(dimensions.getY() - 1, 0);
188  break;
189  case 2:
190  offset = geom::Extent2D(dimensions - geom::Extent2I(1, 1));
191  break;
192  case 3:
193  offset = geom::Extent2D(0, dimensions.getX() - 1);
194  break;
195  }
196  auto rot = geom::LinearTransform::makeRotation(nQuarter * 90.0 * geom::degrees);
197  return transformWcsPixels(wcs, geom::AffineTransform(rot, offset));
198 }
afw::table::PointKey< int > dimensions
Definition: GaussianPsf.cc:48
static LinearTransform makeRotation(Angle t) noexcept
constexpr AngleUnit degrees
constant with units of degrees
Definition: Angle.h:109

◆ transformWcsPixels()

std::shared_ptr< afw::geom::SkyWcs > lsst::meas::astrom::transformWcsPixels ( afw::geom::SkyWcs const &  wcs,
geom::AffineTransform const &  s 
)

Create a new SkyWcs whose pixel coordinate system has been transformed via an affine transform.

Parameters
[in]wcsOriginal SkyWcs object.
[in]sAffineTransform to apply to the pixel coordinate system.
Returns
a new Wcs that satisfies the following:
newWcs = transformWcsPixels(wcs, s);
assert(newWcs.skyToPixel(sky), s(wcs.skyToPixel(sky)));
assert(newWcs.pixelToSky(pixel), wcs.pixelToSky(s.inverted()(pixel)));
table::PointKey< int > pixel
for all sky coordinates sky and pixel coordinates pixel.

Definition at line 173 of file SipTransform.cc.

174  {
175  auto affineTransform22 = afw::geom::makeTransform(s);
176  return afw::geom::makeModifiedWcs(*affineTransform22->inverted(), wcs, true);
177 }
std::shared_ptr< TransformPoint2ToPoint2 > makeTransform(lsst::geom::AffineTransform const &affine)
Wrap an lsst::geom::AffineTransform as a Transform.
std::shared_ptr< SkyWcs > makeModifiedWcs(TransformPoint2ToPoint2 const &pixelTransform, SkyWcs const &wcs, bool modifyActualPixels)
Create a new SkyWcs whose pixels are transformed by pixelTransform, as described below.
Definition: SkyWcs.cc:486