LSST Applications  21.0.0-147-g0e635eb1+1acddb5be5,22.0.0+052faf71bd,22.0.0+1ea9a8b2b2,22.0.0+6312710a6c,22.0.0+729191ecac,22.0.0+7589c3a021,22.0.0+9f079a9461,22.0.1-1-g7d6de66+b8044ec9de,22.0.1-1-g87000a6+536b1ee016,22.0.1-1-g8e32f31+6312710a6c,22.0.1-10-gd060f87+016f7cdc03,22.0.1-12-g9c3108e+df145f6f68,22.0.1-16-g314fa6d+c825727ab8,22.0.1-19-g93a5c75+d23f2fb6d8,22.0.1-19-gb93eaa13+aab3ef7709,22.0.1-2-g8ef0a89+b8044ec9de,22.0.1-2-g92698f7+9f079a9461,22.0.1-2-ga9b0f51+052faf71bd,22.0.1-2-gac51dbf+052faf71bd,22.0.1-2-gb66926d+6312710a6c,22.0.1-2-gcb770ba+09e3807989,22.0.1-20-g32debb5+b8044ec9de,22.0.1-23-gc2439a9a+fb0756638e,22.0.1-3-g496fd5d+09117f784f,22.0.1-3-g59f966b+1e6ba2c031,22.0.1-3-g849a1b8+f8b568069f,22.0.1-3-gaaec9c0+c5c846a8b1,22.0.1-32-g5ddfab5d3+60ce4897b0,22.0.1-4-g037fbe1+64e601228d,22.0.1-4-g8623105+b8044ec9de,22.0.1-5-g096abc9+d18c45d440,22.0.1-5-g15c806e+57f5c03693,22.0.1-7-gba73697+57f5c03693,master-g6e05de7fdc+c1283a92b8,master-g72cdda8301+729191ecac,w.2021.39
LSST Data Management Base Package
approximateWcs.py
Go to the documentation of this file.
1 #
2 # LSST Data Management System
3 # Copyright 2008-2017 LSST Corporation.
4 #
5 # This product includes software developed by the
6 # LSST Project (http://www.lsst.org/).
7 #
8 # This program is free software: you can redistribute it and/or modify
9 # it under the terms of the GNU General Public License as published by
10 # the Free Software Foundation, either version 3 of the License, or
11 # (at your option) any later version.
12 #
13 # This program is distributed in the hope that it will be useful,
14 # but WITHOUT ANY WARRANTY; without even the implied warranty of
15 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 # GNU General Public License for more details.
17 #
18 # You should have received a copy of the LSST License Statement and
19 # the GNU General Public License along with this program. If not,
20 # see <http://www.lsstcorp.org/LegalNotices/>.
21 #
22 
23 __all__ = ["approximateWcs"]
24 
25 import numpy as np
26 
27 import lsst.geom
28 import lsst.afw.table as afwTable
29 import lsst.afw.geom as afwGeom
30 from lsst.meas.base import SingleFrameMeasurementTask
31 from lsst.meas.astrom.sip import makeCreateWcsWithSip
32 from lsst.afw.geom.utils import assertWcsAlmostEqualOverBBox
33 
34 
36  """A fake unit test case class that will enable us to call
37  assertWcsAlmostEqualOverBBox from the method approximateWcs"""
38 
39  def fail(self, msgStr):
40  raise UserWarning("WCS fitting failed " + msgStr)
41 
42 
43 def approximateWcs(wcs, bbox, order=3, nx=20, ny=20, iterations=3,
44  skyTolerance=0.001*lsst.geom.arcseconds, pixelTolerance=0.02, useTanWcs=False):
45  """Approximate an existing WCS as a TAN-SIP WCS
46 
47  The fit is performed by evaluating the WCS at a uniform grid of points
48  within a bounding box.
49 
50  Parameters
51  ----------
52  wcs : `lsst.afw.geom.SkyWcs`
53  wcs to approximate
54  bbox : `lsst.geom.Box2I`
55  the region over which the WCS will be fit
56  order : `int`
57  order of SIP fit
58  nx : `int`
59  number of grid points along x
60  ny : `int`
61  number of grid points along y
62  iterations : `int`
63  number of times to iterate over fitting
64  skyTolerance : `lsst.geom.Angle`
65  maximum allowed difference in world coordinates between
66  input wcs and approximate wcs (default is 0.001 arcsec)
67  pixelTolerance : `float`
68  maximum allowed difference in pixel coordinates between
69  input wcs and approximate wcs (default is 0.02 pixels)
70  useTanWcs : `bool`
71  send a TAN version of wcs to the fitter? It is documented to require that,
72  but I don't think the fitter actually cares
73 
74  Returns
75  -------
76  fitWcs : `lsst.afw.geom.SkyWcs`
77  the fit TAN-SIP WCS
78  """
79  if useTanWcs:
80  crpix = wcs.getPixelOrigin()
81  crval = wcs.getSkyOrigin()
82  cdMatrix = wcs.getCdMatrix(crpix)
83  tanWcs = afwGeom.makeSkyWcs(crpix=crpix, crval=crval, cdMatrix=cdMatrix)
84  else:
85  tanWcs = wcs
86 
87  # create a matchList consisting of a grid of points covering the bbox
88  refSchema = afwTable.SimpleTable.makeMinimalSchema()
89  refCoordKey = afwTable.CoordKey(refSchema["coord"])
90  refCat = afwTable.SimpleCatalog(refSchema)
91 
92  sourceSchema = afwTable.SourceTable.makeMinimalSchema()
93  SingleFrameMeasurementTask(schema=sourceSchema) # expand the schema
94  sourceCentroidKey = afwTable.Point2DKey(sourceSchema["slot_Centroid"])
95 
96  sourceCat = afwTable.SourceCatalog(sourceSchema)
97 
98  matchList = []
99 
100  bboxd = lsst.geom.Box2D(bbox)
101  for x in np.linspace(bboxd.getMinX(), bboxd.getMaxX(), nx):
102  for y in np.linspace(bboxd.getMinY(), bboxd.getMaxY(), ny):
103  pixelPos = lsst.geom.Point2D(x, y)
104  skyCoord = wcs.pixelToSky(pixelPos)
105 
106  refObj = refCat.addNew()
107  refObj.set(refCoordKey, skyCoord)
108 
109  source = sourceCat.addNew()
110  source.set(sourceCentroidKey, pixelPos)
111 
112  matchList.append(afwTable.ReferenceMatch(refObj, source, 0.0))
113 
114  # The TAN-SIP fitter is fitting x and y separately, so we have to iterate to make it converge
115  for indx in range(iterations):
116  sipObject = makeCreateWcsWithSip(matchList, tanWcs, order, bbox)
117  tanWcs = sipObject.getNewWcs()
118  fitWcs = sipObject.getNewWcs()
119 
120  mockTest = _MockTestCase()
121  assertWcsAlmostEqualOverBBox(mockTest, wcs, fitWcs, bbox, maxDiffSky=skyTolerance,
122  maxDiffPix=pixelTolerance)
123 
124  return fitWcs
A FunctorKey used to get or set celestial coordinates from a pair of lsst::geom::Angle keys.
Definition: aggregates.h:210
Custom catalog class for record/table subclasses that are guaranteed to have an ID,...
Definition: SortedCatalog.h:42
A floating-point coordinate rectangle geometry.
Definition: Box.h:413
def assertWcsAlmostEqualOverBBox(testCase, wcs0, wcs1, bbox, maxDiffSky=0.01 *lsst.geom.arcseconds, maxDiffPix=0.01, nx=5, ny=5, msg="WCSs differ")
Definition: utils.py:156
std::shared_ptr< SkyWcs > makeSkyWcs(daf::base::PropertySet &metadata, bool strip=false)
Construct a SkyWcs from FITS keywords.
Definition: SkyWcs.cc:521
def approximateWcs(wcs, bbox, order=3, nx=20, ny=20, iterations=3, skyTolerance=0.001 *lsst.geom.arcseconds, pixelTolerance=0.02, useTanWcs=False)
CreateWcsWithSip< MatchT > makeCreateWcsWithSip(std::vector< MatchT > const &matches, afw::geom::SkyWcs const &linearWcs, int const order, geom::Box2I const &bbox=geom::Box2I(), int const ngrid=0)
Factory function for CreateWcsWithSip.
Lightweight representation of a geometric match between two records.
Definition: Match.h:67