Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Tickets/instrm 1459 #42

Open
wants to merge 5 commits into
base: master
Choose a base branch
from
Open

Tickets/instrm 1459 #42

Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
4cabcff
Return the matched fiducials and minimum distances
RobertLuptonTheGood Nov 17, 2021
7cf9b83
If > 1 measured positions matches a fiducial only keep the best match
RobertLuptonTheGood Nov 17, 2021
3932eb9
Support plotting for both types of supported camera
RobertLuptonTheGood Nov 17, 2021
b737822
Remove trailing whitespace
RobertLuptonTheGood Nov 17, 2021
8ee8b27
Removed code commented out by @yukimoritani post-INSTRM-1398
RobertLuptonTheGood Nov 17, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
151 changes: 102 additions & 49 deletions python/pfs/utils/coordinates/transform.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,16 +10,53 @@


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)
distance = np.empty_like(u)

matched_distances = {}

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)]

return fid_out
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


def fromCameraName(cameraName, *args, **kwargs):
Expand Down Expand Up @@ -161,6 +198,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)
Expand All @@ -178,40 +221,11 @@ 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:
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")

Expand All @@ -227,6 +241,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)
Expand Down Expand Up @@ -283,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
Expand Down Expand Up @@ -325,6 +332,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):
Expand All @@ -341,7 +381,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
Expand All @@ -352,6 +392,13 @@ 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`
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)
Expand All @@ -364,8 +411,11 @@ 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)

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")

Expand All @@ -374,7 +424,10 @@ 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
x, y: position in mcs pixels
Expand Down