117 def fitWcs(self, matches, initWcs, bbox=None, refCat=None, sourceCat=None, exposure=None):
118 """Fit a TAN-SIP WCS from a list of reference object/source matches.
119
120 Parameters
121 ----------
122 matches : `list` of `lsst.afw.table.ReferenceMatch`
123 A sequence of reference object/source matches.
124 The following fields are read:
125 - match.first (reference object) coord
126 - match.second (source) centroid
127
128 The following fields are written:
129 - match.first (reference object) centroid
130 - match.second (source) centroid
131 - match.distance (on sky separation, in radians)
132
133 initWcs : `lsst.afw.geom.SkyWcs`
134 An initial WCS whose CD matrix is used as the final CD matrix.
135 bbox : `lsst.geom.Box2I`
136 The region over which the WCS will be valid (PARENT pixel coordinates);
137 if `None` or an empty box then computed from matches
138 refCat : `lsst.afw.table.SimpleCatalog`
139 Reference object catalog, or `None`.
140 If provided then all centroids are updated with the new WCS,
141 otherwise only the centroids for ref objects in matches are updated.
142 Required fields are "centroid_x", "centroid_y", "coord_ra", and "coord_dec".
143 sourceCat : `lsst.afw.table.SourceCatalog`
144 Source catalog, or `None`.
145 If provided then coords are updated with the new WCS;
146 otherwise only the coords for sources in matches are updated.
147 Required input fields are "slot_Centroid_x", "slot_Centroid_y",
148 "slot_Centroid_xErr", "slot_Centroid_yErr", and optionally
149 "slot_Centroid_x_y_Cov". The "coord_ra" and "coord_dec" fields
150 will be updated but are not used as input.
151 exposure : `lsst.afw.image.Exposure`
152 An Exposure or other displayable image on which matches can be
153 overplotted. Ignored (and may be `None`) if display-based debugging
154 is not enabled via lsstDebug.
155
156 Returns
157 -------
158 An lsst.pipe.base.Struct with the following fields:
159 - wcs : `lsst.afw.geom.SkyWcs`
160 The best-fit WCS.
161 - scatterOnSky : `lsst.geom.Angle`
162 The median on-sky separation between reference objects and
163 sources in "matches", as an `lsst.geom.Angle`
164 """
165 import lsstDebug
169
170 if bbox is None:
172 for match in matches:
173 bbox.include(match.second.getCentroid())
175
176 wcs = self.makeInitialWcs(matches, initWcs)
178
179
180
181
182
183
184
185 revFitter = ScaledPolynomialTransformFitter.fromMatches(self.config.order, matches, wcs,
186 self.config.refUncertainty)
187 revFitter.fit()
188 for nIter in range(self.config.numRejIter):
189 revFitter.updateModel()
190 intrinsicScatter = revFitter.updateIntrinsicScatter()
191 clippedSigma, nRejected = revFitter.rejectOutliers(self.outlierRejectionCtrl)
192 self.log.debug(
193 "Iteration %s: intrinsic scatter is %4.3f pixels, "
194 "rejected %d outliers at %3.2f sigma.",
195 nIter+1, intrinsicScatter, nRejected, clippedSigma
196 )
197 if display:
198 displayFrame = self.display(revFitter, exposure=exposure, bbox=bbox,
199 frame=displayFrame, displayPause=displayPause)
200 revFitter.fit()
201 revScaledPoly = revFitter.getTransform()
202
203
204
205 sipReverse = SipReverseTransform.convert(revScaledPoly, wcs.getPixelOrigin(), cdMatrix)
206
207
208
209
210
211
213 gridBBoxPix.grow(self.config.gridBorder)
214
215
216
217
219 for point in gridBBoxPix.getCorners():
221 gridBBoxIwc.include(cdMatrix(point))
222 fwdFitter = ScaledPolynomialTransformFitter.fromGrid(self.config.order, gridBBoxIwc,
223 self.config.nGridX, self.config.nGridY,
224 revScaledPoly)
225 fwdFitter.fit()
226
227 fwdScaledPoly = fwdFitter.getTransform()
228 sipForward = SipForwardTransform.convert(fwdScaledPoly, wcs.getPixelOrigin(), cdMatrix)
229
230
231
232 wcs = makeWcs(sipForward, sipReverse, wcs.getSkyOrigin())
233
234 if refCat is not None:
235 self.log.debug("Updating centroids in refCat")
237 else:
238 self.log.warning("Updating reference object centroids in match list; refCat is None")
240
241 if sourceCat is not None:
242 self.log.debug("Updating coords in sourceCat")
244 else:
245 self.log.warning("Updating source coords in match list; sourceCat is None")
247
248 self.log.debug("Updating distance in match list")
249 setMatchDistance(matches)
250
251 stats = makeMatchStatisticsInRadians(wcs, matches, lsst.afw.math.MEDIAN)
252 scatterOnSky = stats.getValue()*lsst.geom.radians
253
254 if scatterOnSky.asArcseconds() > self.config.maxScatterArcsec:
255 raise lsst.pipe.base.TaskError(
256 "Fit failed: median scatter on sky = %0.3f arcsec > %0.3f config.maxScatterArcsec" %
257 (scatterOnSky.asArcseconds(), self.config.maxScatterArcsec))
258
259 return lsst.pipe.base.Struct(
260 wcs=wcs,
261 scatterOnSky=scatterOnSky,
262 )
263
A floating-point coordinate rectangle geometry.
An integer coordinate rectangle.
void updateRefCentroids(geom::SkyWcs const &wcs, ReferenceCollection &refList)
Update centroids in a collection of reference objects.
void updateSourceCoords(geom::SkyWcs const &wcs, SourceCollection &sourceList, bool include_covariance=true)
Update sky coordinates in a collection of source objects.