Skip to content

Commit

Permalink
Initial try at new affine fitter with comments.
Browse files Browse the repository at this point in the history
  • Loading branch information
morriscb committed Jul 8, 2019
1 parent f193c2a commit 562905a
Showing 1 changed file with 19 additions and 15 deletions.
34 changes: 19 additions & 15 deletions python/lsst/meas/astrom/fitAffineWcs.py
Original file line number Diff line number Diff line change
Expand Up @@ -57,18 +57,24 @@ def _chi_func(x, ref_points, src_pixels, wcs_maker):
wcs = wcs_maker.makeWcs(x[:2], x[2:].reshape((2, 2)))

output_separations = []
# Fit both sky to pixel and pixel to sky to avoid any non-invertible
# affine matrixes.
for ref, src in zip(ref_points, src_pixels):
sky_sep = ref.getTangentPlaneOffset(wcs.pixelToSky(src))
output_separations.append(sky_sep[0].asArcseconds())
output_separations.append(sky_sep[1].asArcseconds())
xy_sep = src - wcs.skyToPixel(ref)
# Convert the pixel separations to units, arcseconds to match units
# of sky separation.
output_separations.append(xy_sep[0] * wcs.getPixelScale(src).asArcseconds())
output_separations.append(xy_sep[1] * wcs.getPixelScale(src).asArcseconds())

# return np.sum(np.array(output_separations) ** 2)
return output_separations


# Keeping this around for now in case any of the fit parameters need to be
# configurable. Likely the maximum allowed shift magnitude (parameter 2 in the
# fit.)
class FitAffineWcsConfig(pexConfig.Config):
"""Config for FitTanSipWcsTask."""
pass
Expand Down Expand Up @@ -137,13 +143,18 @@ def fitWcs(self,
- ``scatterOnSky`` : median on-sky separation between reference
objects and sources in "matches" (`lsst.afw.geom.Angle`)
"""
# Create a data-structure that decomposes the input Wcs frames and
# appends the new transform.
wcs_maker = TransformedSkyWcsMaker(initWcs)

ref_points = []
src_pixels = []
offset_dir = 0
offset_dist = 0
vect = sphgeom.Vector3d(0, 0, 0)
# Grab reference coordinates and source centroids. Compute the average
# direction and separation between the reference and the sources.
# I'm not sure if bearingTo should be computed from the src to ref
# or ref to source.
for match in matches:
ref_coord = match.first.getCoord()
ref_points.append(ref_coord)
Expand All @@ -157,32 +168,25 @@ def fitWcs(self,
if offset_dir > 180:
offset_dir = offset_dir - 360

# fit = minimize(_chi_func,
# x0=[offset_dir, offset_dist, 1., 0., 0., 1.],
# args=(ref_points,
# src_pixels,
# wcs_maker),
# bounds=((-180, 180),
# (0, None),
# (-2, 2),
# (-2, 2),
# (-2, 2),
# (-2, 2)),
# tol=1e-16)
# Best performing fitter in scipy tried so far (vs. default settings in
# minimize). Fits all current test cases with a scatter of a most 0.15
# arcseconds. exits early because of the xTol value which cannot be
# disabled in scipy1.2.1.
fit = least_squares(_chi_func,
x0=[offset_dir, offset_dist, 1., 0., 0., 1.],
args=(ref_points,
src_pixels,
wcs_maker),
method='dogbox',
bounds=[[-180, 0, -np.inf, -np.inf, -np.inf, -np.inf],
[180, np.inf, np.inf, np.inf, np.inf, np.inf]],
[180, np.inf, np.inf, np.inf, np.inf, np.inf]],
ftol=2.3e-16,
gtol=2.31e-16,
xtol=2.3e-16)

wcs = wcs_maker.makeWcs(fit.x[:2], fit.x[2:].reshape((2, 2)))

# Copied from other fit*WcsTasks.
if refCat is not None:
self.log.debug("Updating centroids in refCat")
lsst.afw.table.updateRefCentroids(wcs, refList=refCat)
Expand Down

0 comments on commit 562905a

Please sign in to comment.