From 4cabcff1f643b864d9dde71f35371bc03046c322 Mon Sep 17 00:00:00 2001 From: Robert Lupton the Good Date: Tue, 16 Nov 2021 14:55:57 -1000 Subject: [PATCH 1/5] Return the matched fiducials and minimum distances --- python/pfs/utils/coordinates/transform.py | 57 ++++++++++++++++++++--- 1 file changed, 50 insertions(+), 7 deletions(-) diff --git a/python/pfs/utils/coordinates/transform.py b/python/pfs/utils/coordinates/transform.py index d866ded..bd52490 100644 --- a/python/pfs/utils/coordinates/transform.py +++ b/python/pfs/utils/coordinates/transform.py @@ -10,16 +10,42 @@ def matchIds(u, v, x, y, fid, matchRadius=2): - """Given a set of points (x, y, fid) which define the ids for the points (x, y) - and a set of points (u, v), return an array of the ids for (u, v) + """Match set of measured points (u, v) to a set of points (x, y) with ids fid + + If there are multiple possible matches from (u, v) to (x, y), only return the best match + + Parameters + ---------- + u : `array` of `float` + x-positions of points to be matched + v : `array` of `float` + y-positions of points to be matched + x : `array` of `float` + known x-positions of points + y : `array` of `float` + known y-positions of points + fid : `array` of `int` + ids for points (x, y) + matchRadius: `float` + maximum allowed match radius + + Returns + ------- + fids: `array` of `int` of length len(u) + matched values of fid or -1 + distance: array of `float` of length len(u) + distance from (u, v) to nearest element of (x, y) """ fid_out = np.full_like(u, -1, dtype=fid.dtype) for i, (up, vp) in enumerate(zip(u, v)): d = np.hypot(x - up, y - vp) - if min(d) < matchRadius: - fid_out[i] = fid[np.argmin(d)] + distance[i] = min(d) + if distance[i] < matchRadius: + matched_fid = fid[np.argmin(d)] + + fid_out[i] = matched_fid - return fid_out + return fid_out, distance def fromCameraName(cameraName, *args, **kwargs): @@ -161,6 +187,12 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=2, nMatchMin=0.75, fi nMatchMin: Minimum number of permissible matches (if <= 1, interpreted as the fraction of the number of fiducials) fig matplotlib figure for displays; or None + + Returns: + fids: array of `int` + array of length mcs_data giving indices into fiducials or -1 + distance: array of `float` + array of length mcs_data giving distance in mm to nearest fiducial """ if nMatchMin <= 1: nMatchMin *= len(fiducials.fiducialId) @@ -178,7 +210,7 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=2, nMatchMin=0.75, fi xd, yd = ptd.mcsToPfi(mcs_x_pix, mcs_y_pix) del ptd - fid = matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) + fid, dmin = matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) nMatch = sum(fid > 0) if fig is not None: @@ -227,6 +259,8 @@ def create_artists(self, legend, orig_handle, res = scipy.optimize.minimize(distortion, distortion.getArgs(), method='Powell') self.mcsDistort.setArgs(res.x) + xd, yd = self.mcsToPfi(mcs_x_pix, mcs_y_pix) + return matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) # # Note that these two routines use instance values of mcs_boresight_[xy]_pix # and assume that the pfi origin is at (0, 0) @@ -352,6 +386,12 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75): matchRadius: Radius to match points and fiducials (mm) nMatchMin: Minimum number of permissible matches (if <= 1, interpreted as the fraction of the number of fiducials) + + Returns: + fids: array of `int` + array of length mcs_data giving indices into fiducials or -1 + distance: array of `float` + array of length mcs_data giving distance in mm to nearest fiducial """ if nMatchMin <= 1: nMatchMin *= len(fiducials.fiducialId) @@ -364,7 +404,7 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75): fiducialId = fiducials.fiducialId.to_numpy() xd, yd = self.mcsToPfi(mcs_x_pix, mcs_y_pix) - fid = matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) + fid, dmin = matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) nMatch = sum(fid > 0) if nMatch < nMatchMin: raise RuntimeError(f"I only matched {nMatch} out of {len(fiducialId)} fiducial fibres") @@ -374,6 +414,9 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75): asrdDistort.setArgs(res.x) self.mcsDistort = asrdDistort + + xd, yd = self.mcsToPfi(mcs_x_pix, mcs_y_pix) + return matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) def mcsToPfi(self, x, y): """transform ASRD 71M camera pixels to pfi mm From 7cf9b83b22f6bec90819ca8ba594d3c90b62736a Mon Sep 17 00:00:00 2001 From: Robert Lupton the Good Date: Tue, 16 Nov 2021 14:56:43 -1000 Subject: [PATCH 2/5] If > 1 measured positions matches a fiducial only keep the best match --- python/pfs/utils/coordinates/transform.py | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/python/pfs/utils/coordinates/transform.py b/python/pfs/utils/coordinates/transform.py index bd52490..55fe3b8 100644 --- a/python/pfs/utils/coordinates/transform.py +++ b/python/pfs/utils/coordinates/transform.py @@ -37,12 +37,23 @@ def matchIds(u, v, x, y, fid, matchRadius=2): distance from (u, v) to nearest element of (x, y) """ fid_out = np.full_like(u, -1, dtype=fid.dtype) + distance = np.empty_like(u) + + matched_distances = {} + for i, (up, vp) in enumerate(zip(u, v)): d = np.hypot(x - up, y - vp) distance[i] = min(d) if distance[i] < matchRadius: matched_fid = fid[np.argmin(d)] + if matched_fid in matched_distances: # we've already found a match + if distance[i] > matched_distances[matched_fid]: # the old one is better + continue + + fid_out[fid_out == matched_fid] = -1 # invalidate old match + + matched_distances[matched_fid] = distance[i] fid_out[i] = matched_fid return fid_out, distance From 3932eb91344e00cd9520aebe82bd04fa16ad389a Mon Sep 17 00:00:00 2001 From: Robert Lupton the Good Date: Tue, 16 Nov 2021 14:57:30 -1000 Subject: [PATCH 3/5] Support plotting for both types of supported camera Moving the code out into a method makes support for any future cameras easy too --- python/pfs/utils/coordinates/transform.py | 72 +++++++++++++---------- 1 file changed, 40 insertions(+), 32 deletions(-) diff --git a/python/pfs/utils/coordinates/transform.py b/python/pfs/utils/coordinates/transform.py index 55fe3b8..1a39154 100644 --- a/python/pfs/utils/coordinates/transform.py +++ b/python/pfs/utils/coordinates/transform.py @@ -224,37 +224,8 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=2, nMatchMin=0.75, fi fid, dmin = matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) nMatch = sum(fid > 0) - if fig is not None: - from matplotlib.legend_handler import HandlerPatch - class HandlerEllipse(HandlerPatch): - def create_artists(self, legend, orig_handle, - xdescent, ydescent, width, height, fontsize, trans): - center = 0.5 * width - 0.5 * xdescent, 0.5 * height - 0.5 * ydescent - p = Circle(center, matchRadius - #height=height + ydescent - ) - self.update_prop(p, orig_handle, legend) - p.set_transform(trans) - return [p] - - - ax = fig.gca() - - for x, y in zip(x_fid_mm, y_fid_mm): - c = Circle((x, y), matchRadius, color='red', alpha=0.5) - ax.add_patch(c) - #plt.plot(x_fid_mm, y_fid_mm, 'o', label="fiducials") - plt.plot(xd, yd, '.', label="detections") - plt.plot(xd[fid > 0], yd[fid > 0], '+', color='black', label="matches") - - handles, labels = ax.get_legend_handles_labels() - if True: - handles += [c] - labels += ["search"] - plt.legend(handles, labels, handler_map={Circle: HandlerEllipse()}) - plt.title(f"Matched {nMatch} points") - plt.gca().set_aspect('equal') - + self._plotMatches(fig, x_fid_mm, y_fid_mm, xd, yd, fid, matchRadius, nMatch) + if nMatch < nMatchMin: raise RuntimeError(f"I only matched {nMatch} out of {len(fiducialId)} fiducial fibres") @@ -370,6 +341,39 @@ def pfiToMcs(self, x, y, niter=5, lam=1.0): return xyout[0], xyout[1] + def _plotMatches(self, fig, x_fid_mm, y_fid_mm, xd, yd, fid, matchRadius, nMatch): + if fig is None: + return + + from matplotlib.legend_handler import HandlerPatch + class HandlerEllipse(HandlerPatch): + def create_artists(self, legend, orig_handle, + xdescent, ydescent, width, height, fontsize, trans): + center = 0.5 * width - 0.5 * xdescent, 0.5 * height - 0.5 * ydescent + p = Circle(center, matchRadius + #height=height + ydescent + ) + self.update_prop(p, orig_handle, legend) + p.set_transform(trans) + return [p] + + + ax = fig.gca() + + for x, y in zip(x_fid_mm, y_fid_mm): + c = Circle((x, y), matchRadius, color='red', alpha=0.5) + ax.add_patch(c) + #plt.plot(x_fid_mm, y_fid_mm, 'o', label="fiducials") + plt.plot(xd, yd, '.', label="detections") + plt.plot(xd[fid > 0], yd[fid > 0], '+', color='black', label="matches") + + handles, labels = ax.get_legend_handles_labels() + if True: + handles += [c] + labels += ["search"] + plt.legend(handles, labels, handler_map={Circle: HandlerEllipse()}) + plt.title(f"Matched {nMatch} points") + plt.gca().set_aspect('equal') class ASRD71MTransform(PfiTransform): def __init__(self, altitude=90, insrot=0, applyDistortion=True): @@ -386,7 +390,7 @@ def setParams(self, altitude=90, insrot=0): self.mcsDistort.setArgs([-3.76006171e+02, -2.68710420e+02, -9.24753269e-01, -5.75212519e-01, -2.25647580e-13]) - def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75): + def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75, fig=None): """Update our estimate of the transform, based on the positions of fiducial fibres mcs_data: Pandas DataFrame containing mcs_center_x_pix, mcs_center_y_pix @@ -397,6 +401,7 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75): matchRadius: Radius to match points and fiducials (mm) nMatchMin: Minimum number of permissible matches (if <= 1, interpreted as the fraction of the number of fiducials) + fig matplotlib figure for displays; or None Returns: fids: array of `int` @@ -417,6 +422,9 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75): xd, yd = self.mcsToPfi(mcs_x_pix, mcs_y_pix) fid, dmin = matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) nMatch = sum(fid > 0) + + self._plotMatches(fig, x_fid_mm, y_fid_mm, xd, yd, fid, matchRadius, nMatch) + if nMatch < nMatchMin: raise RuntimeError(f"I only matched {nMatch} out of {len(fiducialId)} fiducial fibres") From b73782257a55a1e75fbb30b57428057babee4b29 Mon Sep 17 00:00:00 2001 From: Robert Lupton the Good Date: Tue, 16 Nov 2021 15:00:59 -1000 Subject: [PATCH 4/5] Remove trailing whitespace --- python/pfs/utils/coordinates/transform.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/python/pfs/utils/coordinates/transform.py b/python/pfs/utils/coordinates/transform.py index 1a39154..48b449b 100644 --- a/python/pfs/utils/coordinates/transform.py +++ b/python/pfs/utils/coordinates/transform.py @@ -52,7 +52,7 @@ def matchIds(u, v, x, y, fid, matchRadius=2): continue fid_out[fid_out == matched_fid] = -1 # invalidate old match - + matched_distances[matched_fid] = distance[i] fid_out[i] = matched_fid @@ -225,7 +225,7 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=2, nMatchMin=0.75, fi nMatch = sum(fid > 0) self._plotMatches(fig, x_fid_mm, y_fid_mm, xd, yd, fid, matchRadius, nMatch) - + if nMatch < nMatchMin: raise RuntimeError(f"I only matched {nMatch} out of {len(fiducialId)} fiducial fibres") @@ -424,7 +424,7 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75, fi nMatch = sum(fid > 0) self._plotMatches(fig, x_fid_mm, y_fid_mm, xd, yd, fid, matchRadius, nMatch) - + if nMatch < nMatchMin: raise RuntimeError(f"I only matched {nMatch} out of {len(fiducialId)} fiducial fibres") @@ -436,7 +436,7 @@ def updateTransform(self, mcs_data, fiducials, matchRadius=1, nMatchMin=0.75, fi xd, yd = self.mcsToPfi(mcs_x_pix, mcs_y_pix) return matchIds(xd, yd, x_fid_mm, y_fid_mm, fiducialId, matchRadius=matchRadius) - + def mcsToPfi(self, x, y): """transform ASRD 71M camera pixels to pfi mm x, y: position in mcs pixels From 8ee8b2773f308b394f318c881fc17a65bcadd9fc Mon Sep 17 00:00:00 2001 From: Robert Lupton the Good Date: Tue, 16 Nov 2021 15:06:26 -1000 Subject: [PATCH 5/5] Removed code commented out by @yukimoritani post-INSTRM-1398 --- python/pfs/utils/coordinates/transform.py | 9 --------- 1 file changed, 9 deletions(-) diff --git a/python/pfs/utils/coordinates/transform.py b/python/pfs/utils/coordinates/transform.py index 48b449b..2ebaad4 100644 --- a/python/pfs/utils/coordinates/transform.py +++ b/python/pfs/utils/coordinates/transform.py @@ -299,15 +299,6 @@ def pfiToMcs(self, x, y, niter=5, lam=1.0): xyout = CoordinateTransform(xyin, 90.0 - self.altitude, "pfi_mcs", inr=self.insrot) # first guess at MCS coordinates # - # Apparently there's a missing rotation by pi/2 in the pfs_utils code - # Implemented under INSRTM-1398 - # - """ - tmp = xyout.copy() - xyout[0], xyout[1] = tmp[1], -tmp[0] - del tmp - """ - # # Deal with coordinate transformations on the MCS # xyout = -xyout # rotate 180 about centre