Merge external CSM extension
This commit is contained in:
Joel Collins
parent
86e326a1e4
commit
55e423c159
10 changed files with 25 additions and 1311 deletions
|
|
@ -23,6 +23,6 @@ with handle_extension_error("zip builder"):
|
||||||
with handle_extension_error("autostorage"):
|
with handle_extension_error("autostorage"):
|
||||||
from .autostorage import autostorage_extension_v2
|
from .autostorage import autostorage_extension_v2
|
||||||
with handle_extension_error("camera stage mapping"):
|
with handle_extension_error("camera stage mapping"):
|
||||||
from .camera_stage_mapping import csm_extension
|
from camera_stage_mapping.ofm_extension import csm_extension
|
||||||
with handle_extension_error("lens shading calibration"):
|
with handle_extension_error("lens shading calibration"):
|
||||||
from .picamera_autocalibrate import lst_extension_v2
|
from .picamera_autocalibrate import lst_extension_v2
|
||||||
|
|
|
||||||
|
|
@ -1 +0,0 @@
|
||||||
from .extension import csm_extension
|
|
||||||
|
|
@ -1,102 +0,0 @@
|
||||||
# -*- coding: utf-8 -*-
|
|
||||||
"""
|
|
||||||
Created on Tue May 26 08:08:14 2015
|
|
||||||
|
|
||||||
@author: rwb27
|
|
||||||
"""
|
|
||||||
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
|
|
||||||
class AttributeDict(dict):
|
|
||||||
"""This class extends a dictionary to have a "create" method for
|
|
||||||
compatibility with h5py attrs objects."""
|
|
||||||
|
|
||||||
def create(self, name, data):
|
|
||||||
self[name] = data
|
|
||||||
|
|
||||||
def modify(self, name, data):
|
|
||||||
self[name] = data
|
|
||||||
|
|
||||||
def copy_arrays(self):
|
|
||||||
"""Replace any numpy.ndarray in the dict with a copy, to break any unintentional links."""
|
|
||||||
for k in list(self.keys()):
|
|
||||||
if isinstance(self[k], np.ndarray):
|
|
||||||
self[k] = np.copy(self[k])
|
|
||||||
|
|
||||||
|
|
||||||
def ensure_attribute_dict(obj, copy=False):
|
|
||||||
"""Given a mapping that may or not be an AttributeDict, return an
|
|
||||||
AttributeDict object that either is, or copies the data of, the input."""
|
|
||||||
if isinstance(obj, AttributeDict) and not copy:
|
|
||||||
return obj
|
|
||||||
else:
|
|
||||||
out = AttributeDict(obj)
|
|
||||||
if copy:
|
|
||||||
out.copy_arrays()
|
|
||||||
return out
|
|
||||||
|
|
||||||
|
|
||||||
def ensure_attrs(obj):
|
|
||||||
"""Return an ArrayWithAttrs version of an array-like object, may be the
|
|
||||||
original object if it already has attrs."""
|
|
||||||
if hasattr(obj, "attrs"):
|
|
||||||
return obj # if it has attrs, do nothing
|
|
||||||
else:
|
|
||||||
return ArrayWithAttrs(obj) # otherwise, wrap it
|
|
||||||
|
|
||||||
|
|
||||||
class ArrayWithAttrs(np.ndarray):
|
|
||||||
"""A numpy ndarray, with an AttributeDict accessible as array.attrs.
|
|
||||||
|
|
||||||
This class is intended as a temporary version of an h5py dataset to allow
|
|
||||||
the easy passing of metadata/attributes around nplab functions. It owes
|
|
||||||
a lot to the ``InfoArray`` example in `numpy` documentation on subclassing
|
|
||||||
`numpy.ndarray`.
|
|
||||||
"""
|
|
||||||
|
|
||||||
def __new__(cls, input_array, attrs={}):
|
|
||||||
"""Make a new ndarray, based on an existing one, with an attrs dict.
|
|
||||||
|
|
||||||
This function adds an attributes dictionary to a numpy array, to make
|
|
||||||
it work like an h5py dataset. It doesn't copy data if it can be
|
|
||||||
avoided."""
|
|
||||||
# the input array should be a numpy array, then we cast it to this type
|
|
||||||
obj = np.asarray(input_array).view(cls)
|
|
||||||
# next, add the dict
|
|
||||||
# ensure_attribute_dict always returns an AttributeDict
|
|
||||||
obj.attrs = ensure_attribute_dict(attrs)
|
|
||||||
# return the new object
|
|
||||||
return obj
|
|
||||||
|
|
||||||
def __array_finalize__(self, obj):
|
|
||||||
# this is called by numpy when the object is created (__new__ may or
|
|
||||||
# may not get called)
|
|
||||||
if obj is None:
|
|
||||||
return # if obj is None, __new__ was called - do nothing
|
|
||||||
# if we didn't create the object with __new__, we must add the attrs
|
|
||||||
# dictionary. We copy this from the source object if possible (while
|
|
||||||
# ensuring it's the right type) or create a new, empty one if not.
|
|
||||||
# NB we don't use ensure_attribute_dict because we want to make sure the
|
|
||||||
# dict object is *copied* not merely referenced.
|
|
||||||
self.attrs = ensure_attribute_dict(getattr(obj, "attrs", {}), copy=True)
|
|
||||||
|
|
||||||
|
|
||||||
def attribute_bundler(attrs):
|
|
||||||
"""Return a function that bundles the supplied attributes with an array."""
|
|
||||||
|
|
||||||
def bundle_attrs(array):
|
|
||||||
return ArrayWithAttrs(array, attrs=attrs)
|
|
||||||
|
|
||||||
|
|
||||||
class DummyHDF5Group(dict):
|
|
||||||
def __init__(self, dictionary, attrs={}, name="DummyHDF5Group"):
|
|
||||||
super(DummyHDF5Group, self).__init__()
|
|
||||||
self.attrs = attrs
|
|
||||||
for key in dictionary:
|
|
||||||
self[key] = dictionary[key]
|
|
||||||
self.name = name
|
|
||||||
self.basename = name
|
|
||||||
|
|
||||||
file = None
|
|
||||||
parent = None
|
|
||||||
|
|
@ -1,307 +0,0 @@
|
||||||
"""
|
|
||||||
1D calibration of the relationship between a stage and a camera
|
|
||||||
|
|
||||||
The `Tracker` class in this file is used to simplify code for tasks that involve moving the
|
|
||||||
stage, and tracking the corresponding motion with the camera.
|
|
||||||
|
|
||||||
(c) Richard Bowman 2019, released under GNU GPL v3
|
|
||||||
"""
|
|
||||||
import numpy as np
|
|
||||||
import time
|
|
||||||
from numpy.linalg import norm
|
|
||||||
from .camera_stage_tracker import Tracker, move_until_motion_detected
|
|
||||||
import logging
|
|
||||||
|
|
||||||
|
|
||||||
def displacements(positions):
|
|
||||||
"""Calculate the absolute distance of each point from the first point."""
|
|
||||||
return norm(positions - positions[0, :][np.newaxis, :], axis=1)
|
|
||||||
|
|
||||||
|
|
||||||
def direction_from_points(points):
|
|
||||||
"""Given an Nx2 array of points, figure out the principal component.
|
|
||||||
|
|
||||||
The return value is a normalised vector that points along the
|
|
||||||
direction with the most motion.
|
|
||||||
"""
|
|
||||||
points = points.astype(np.float)
|
|
||||||
points -= np.mean(points, axis=0)[np.newaxis, :]
|
|
||||||
eigenvalues, eigenvectors = np.linalg.eig(np.cov(points.T))
|
|
||||||
return eigenvectors[:, np.argmax(eigenvalues)]
|
|
||||||
|
|
||||||
|
|
||||||
def apply_backlash(x, backlash=0, start_unwound=True):
|
|
||||||
"""Apply a basic model of backlash to a set of coordinates.
|
|
||||||
|
|
||||||
The output (y) will lag behind the input by up to `backlash`
|
|
||||||
|
|
||||||
`start_unwound` (default: True) assumes we change direction
|
|
||||||
at the start of the time series, so you will get no motion
|
|
||||||
until `x[i]` has moved by at least `2*backlash`.
|
|
||||||
"""
|
|
||||||
y = np.zeros_like(x)
|
|
||||||
if start_unwound:
|
|
||||||
initial_direction = np.sign(x[1] - x[0])
|
|
||||||
y[0] = x[0] + initial_direction * backlash
|
|
||||||
else:
|
|
||||||
y[0] = x[0]
|
|
||||||
for i in range(1, len(x)):
|
|
||||||
d = x[i] - y[i - 1]
|
|
||||||
if np.abs(d) >= backlash:
|
|
||||||
y[i] = x[i] - np.sign(d) * backlash
|
|
||||||
else:
|
|
||||||
y[i] = y[i - 1]
|
|
||||||
return y
|
|
||||||
|
|
||||||
|
|
||||||
def fit_backlash(moves):
|
|
||||||
"""Given a set of linear moves forwards and back, estimate backlash.
|
|
||||||
|
|
||||||
The result is an estimate of the amount of backlash, and the ratio
|
|
||||||
of steps to pixels. The moves should be in the same
|
|
||||||
format as `Tracker.history`.
|
|
||||||
|
|
||||||
We use a very basic fitting method: we do a brute-force search for
|
|
||||||
the backlash value, and for each value of backlash we fit a line to
|
|
||||||
the relationship between stage position (after modelling backlash)
|
|
||||||
and image position. We then pick the value of backlash that gets
|
|
||||||
the lowest residuals. Currently the backlash values tried will
|
|
||||||
start at 0 and increase by 1 or by a factor of 1.33 each time.
|
|
||||||
|
|
||||||
The return value is a dictionary with the following keys:
|
|
||||||
backlash: float
|
|
||||||
the estimated backlash, in motor steps
|
|
||||||
pixels_per_step: float
|
|
||||||
the gradient of pixels to steps
|
|
||||||
fractional_error: float
|
|
||||||
an estimate of the goodness of fit
|
|
||||||
stage_direction: numpy.ndarray
|
|
||||||
unit vector in the direction of stage motion
|
|
||||||
image_direction: numpy.ndarray
|
|
||||||
unit vector in the direction of the motion measured
|
|
||||||
on the camera
|
|
||||||
pixels_per_step_vector: numpy.ndarray
|
|
||||||
The displacement in 2D on the camera resulting from
|
|
||||||
one step in `stage_direction`. This is equal to the
|
|
||||||
product of `pixels_per_step` and `image_direction`.
|
|
||||||
"""
|
|
||||||
all_stage_points, all_image_points = moves
|
|
||||||
|
|
||||||
# Figure out the direction of motion, and reduce everything to 1D
|
|
||||||
image_direction = direction_from_points(all_image_points)
|
|
||||||
stage_direction = direction_from_points(all_stage_points)
|
|
||||||
xfit = np.sum(all_stage_points * stage_direction[np.newaxis, :], axis=1)
|
|
||||||
yfit = np.sum(all_image_points * image_direction[np.newaxis, :], axis=1)
|
|
||||||
|
|
||||||
# We should probably use a fancy optimiser to fit the backlash, but
|
|
||||||
# brute-forcing it is reliable and doesn't take long.
|
|
||||||
def fit_motion(xfit, yfit, backlash=0):
|
|
||||||
"""Using the model of backlash, fit the observed camera motion"""
|
|
||||||
xfit_blsh = apply_backlash(xfit, backlash)
|
|
||||||
xfit_blsh -= np.mean(xfit_blsh)
|
|
||||||
m, c = np.polyfit(xfit_blsh, yfit, 1)
|
|
||||||
residuals = yfit - (xfit_blsh * m + c)
|
|
||||||
return m, c, np.std(residuals, ddof=3)
|
|
||||||
|
|
||||||
max_backlash = (np.max(xfit) - np.min(xfit)) / 3
|
|
||||||
backlash_values = []
|
|
||||||
residual_values = []
|
|
||||||
backlash = 0
|
|
||||||
while backlash < max_backlash:
|
|
||||||
m, c, residual = fit_motion(xfit, yfit, backlash)
|
|
||||||
residual_values.append(residual)
|
|
||||||
backlash_values.append(backlash)
|
|
||||||
backlash += max(1, backlash / 3)
|
|
||||||
|
|
||||||
backlash = backlash_values[np.argmin(residual_values)]
|
|
||||||
m, c, residual = fit_motion(xfit, yfit, backlash)
|
|
||||||
|
|
||||||
fractional_error = residual / norm(np.diff(yfit))
|
|
||||||
if fractional_error > 0.1:
|
|
||||||
raise ValueError("The fit didn't look successful")
|
|
||||||
|
|
||||||
return {
|
|
||||||
"backlash": backlash,
|
|
||||||
"pixels_per_step": m,
|
|
||||||
"fractional_error": fractional_error,
|
|
||||||
"stage_direction": stage_direction,
|
|
||||||
"image_direction": image_direction,
|
|
||||||
"pixels_per_step_vector": m * image_direction,
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
def calibrate_backlash_1d(tracker, move, direction=np.array([1, 0, 0])):
|
|
||||||
"""Figure out reasonable step sizes for calibration, and estimate the backlash."""
|
|
||||||
try: # Ensure that the tracker has a template set
|
|
||||||
_ = tracker.template
|
|
||||||
except:
|
|
||||||
tracker.acquire_template()
|
|
||||||
assert tracker.stage_positions.shape[0] == 1
|
|
||||||
original_stage_pos = tracker.stage_positions[-1, :]
|
|
||||||
|
|
||||||
direction = (
|
|
||||||
direction / np.sum(direction ** 2) ** 0.5
|
|
||||||
) # ensure "direction" is normalised
|
|
||||||
|
|
||||||
logging.info("Moving the stage until we see motion...")
|
|
||||||
# Move the stage until we can see a significant amount of motion
|
|
||||||
i, m = move_until_motion_detected(
|
|
||||||
tracker, move, direction, threshold=tracker.max_safe_displacement * 0.2
|
|
||||||
)
|
|
||||||
|
|
||||||
logging.info("Moving the stage to the edge of the field of view...")
|
|
||||||
i, m = move_until_motion_detected(
|
|
||||||
tracker,
|
|
||||||
move,
|
|
||||||
direction,
|
|
||||||
threshold=tracker.max_safe_displacement * 0.7,
|
|
||||||
multipliers=m / 2.0 * np.arange(20),
|
|
||||||
detect_cumulative_motion=True,
|
|
||||||
)
|
|
||||||
exponential_moves = tracker.history
|
|
||||||
|
|
||||||
# Include this final step, and make a rough estimate of the scaling from stage to image
|
|
||||||
stage_pos, image_pos = tracker.history
|
|
||||||
stage_step = stage_pos[-1, :] - stage_pos[-1 - i, :]
|
|
||||||
image_step = image_pos[-1, :] - image_pos[-1 - i, :]
|
|
||||||
steps_per_pixel = norm(stage_step) / norm(image_step)
|
|
||||||
|
|
||||||
# Calculate a step that moves roughly 0.2 times the max. displacement (i.e. 0.1 times the FoV)
|
|
||||||
sensible_step = direction * tracker.max_safe_displacement * 0.2 * steps_per_pixel
|
|
||||||
tracker.reset_history()
|
|
||||||
|
|
||||||
logging.info("Moving the stage backwards to measure backlash (1/2)")
|
|
||||||
# Now move backwards, in 10 steps that should roughly cross the field of view.
|
|
||||||
# If the stage has no backlash, this will move too far, hence the break statement to
|
|
||||||
# prevent it moving outside of the field of view.
|
|
||||||
starting_stage_pos, starting_camera_pos = tracker.append_point()
|
|
||||||
for i in range(15):
|
|
||||||
move(starting_stage_pos - sensible_step * (i + 1))
|
|
||||||
# print(".", end="")
|
|
||||||
stage_pos, image_pos = tracker.append_point()
|
|
||||||
if (
|
|
||||||
i > 3
|
|
||||||
and tracker.moving_away_from_centre
|
|
||||||
and norm(image_pos) > 0.65 * tracker.max_safe_displacement
|
|
||||||
):
|
|
||||||
break # Stop once we have moved far enough
|
|
||||||
|
|
||||||
logging.info("Moving the stage forwards to measure backlash (2/2)")
|
|
||||||
# Move forwards again, in 10 steps
|
|
||||||
starting_stage_pos, starting_camera_pos = tracker.append_point()
|
|
||||||
for i in range(15):
|
|
||||||
move(starting_stage_pos + sensible_step * (i + 1))
|
|
||||||
# print(".", end="")
|
|
||||||
stage_pos, image_pos = tracker.append_point()
|
|
||||||
if (
|
|
||||||
i > 3
|
|
||||||
and tracker.moving_away_from_centre
|
|
||||||
and norm(image_pos) > 0.65 * tracker.max_safe_displacement
|
|
||||||
):
|
|
||||||
break # Stop once we have moved far enough
|
|
||||||
linear_moves = tracker.history
|
|
||||||
|
|
||||||
try:
|
|
||||||
res = fit_backlash(linear_moves)
|
|
||||||
backlash_correction = (
|
|
||||||
sensible_step / norm(sensible_step) * res["backlash"] * 1.5
|
|
||||||
)
|
|
||||||
|
|
||||||
# Finally, move back to the starting position, doing backlash-corrected moves.
|
|
||||||
logging.info("Moving back to the start, correcting for backlash...")
|
|
||||||
tracker.reset_history()
|
|
||||||
stage_pos, camera_pos = tracker.append_point()
|
|
||||||
while np.dot(stage_pos - sensible_step - original_stage_pos, sensible_step) > 0:
|
|
||||||
move(stage_pos - sensible_step - backlash_correction)
|
|
||||||
move(stage_pos - sensible_step)
|
|
||||||
stage_pos, camera_pos = tracker.append_point()
|
|
||||||
backlash_corrected_moves = tracker.history
|
|
||||||
move(original_stage_pos - backlash_correction)
|
|
||||||
except ValueError:
|
|
||||||
return {"exponential_moves": exponential_moves, "linear_moves": linear_moves}
|
|
||||||
finally:
|
|
||||||
# Reset position
|
|
||||||
move(original_stage_pos)
|
|
||||||
|
|
||||||
logging.info(f"Estimated backlash {res['backlash']:.0f} steps")
|
|
||||||
logging.info(
|
|
||||||
f"Stage-to-image ratio {np.abs(res['pixels_per_step']):.3f} pixels/step"
|
|
||||||
)
|
|
||||||
logging.info(
|
|
||||||
f"Residuals were about {res['fractional_error']:.2f} times the step size"
|
|
||||||
)
|
|
||||||
|
|
||||||
res.update(
|
|
||||||
{
|
|
||||||
"exponential_moves": exponential_moves,
|
|
||||||
"linear_moves": linear_moves,
|
|
||||||
"backlash_corrected_moves": backlash_corrected_moves,
|
|
||||||
}
|
|
||||||
)
|
|
||||||
return res
|
|
||||||
|
|
||||||
|
|
||||||
def plot_1d_backlash_calibration(results):
|
|
||||||
"""Plot the results of a calibration run"""
|
|
||||||
from matplotlib import pyplot as plt
|
|
||||||
|
|
||||||
f, ax = plt.subplots(1, 2)
|
|
||||||
|
|
||||||
for k in ["exponential", "linear", "backlash_corrected"]:
|
|
||||||
moves = results[k + "_moves"]
|
|
||||||
if moves is not None:
|
|
||||||
ax[0].plot(moves[1][:, 0], moves[1][:, 1], "o-")
|
|
||||||
ax[0].set_aspect(1, adjustable="datalim")
|
|
||||||
|
|
||||||
image_direction = results["image_direction"]
|
|
||||||
stage_direction = results["stage_direction"]
|
|
||||||
|
|
||||||
def convert_moves(moves):
|
|
||||||
stage_pos, image_pos = moves
|
|
||||||
stage_1d = np.sum(stage_pos * stage_direction[np.newaxis, :], axis=1)
|
|
||||||
image_1d = np.sum(image_pos * image_direction[np.newaxis, :], axis=1)
|
|
||||||
return stage_1d, image_1d
|
|
||||||
|
|
||||||
ax[1].plot(*convert_moves(results["exponential_moves"]), "o-")
|
|
||||||
|
|
||||||
stage_pos, image_pos = convert_moves(results["linear_moves"])
|
|
||||||
model = apply_backlash(stage_pos, results["backlash"])
|
|
||||||
model *= results["pixels_per_step"]
|
|
||||||
model += np.mean(image_pos) - np.mean(model)
|
|
||||||
ax[1].plot(stage_pos, model, "-")
|
|
||||||
ax[1].plot(stage_pos, image_pos, "o")
|
|
||||||
if results["backlash_corrected_moves"] is not None:
|
|
||||||
ax[1].plot(*convert_moves(results["backlash_corrected_moves"]), "+")
|
|
||||||
|
|
||||||
return f, ax
|
|
||||||
|
|
||||||
|
|
||||||
def image_to_stage_displacement_from_1d(calibrations):
|
|
||||||
"""Combine X and Y calibrations
|
|
||||||
|
|
||||||
This uses the output from `calibrate_backlash_1d`, run at least
|
|
||||||
twice with orthogonal (or at least different) `direction` parameters.
|
|
||||||
The resulting 2x2 transformation matrix should map from image
|
|
||||||
to stage coordinates. Currently, the backlash estimate given
|
|
||||||
by this function is only really trustworthy if you've supplied
|
|
||||||
two orthogonal calibrations - that will usually be the case.
|
|
||||||
"""
|
|
||||||
stage_vectors = []
|
|
||||||
image_vectors = []
|
|
||||||
backlash = np.zeros(3)
|
|
||||||
for cal in calibrations:
|
|
||||||
stage_vectors.append(cal["stage_direction"][:2])
|
|
||||||
image_vectors.append(cal["pixels_per_step_vector"])
|
|
||||||
# our backlash estimate will be the maximum backlash
|
|
||||||
# measured in each direction
|
|
||||||
c_blash = np.abs(cal["backlash"] * cal["stage_direction"])
|
|
||||||
backlash[backlash < c_blash] = c_blash[backlash < c_blash]
|
|
||||||
|
|
||||||
A, res, rank, s = np.linalg.lstsq(
|
|
||||||
image_vectors, stage_vectors
|
|
||||||
) # we solve image*A = stage
|
|
||||||
return {
|
|
||||||
"image_to_stage_displacement": A,
|
|
||||||
"backlash_vector": backlash,
|
|
||||||
"backlash": np.max(backlash),
|
|
||||||
}
|
|
||||||
|
|
@ -1,90 +0,0 @@
|
||||||
"""
|
|
||||||
Camera-stage calibration, 2D
|
|
||||||
|
|
||||||
Uses 2D motion to try to calibrate the relationship between a camera and a stage.
|
|
||||||
|
|
||||||
|
|
||||||
(c) Richard Bowman 2019, released under GNU GPL v3
|
|
||||||
"""
|
|
||||||
import numpy as np
|
|
||||||
import time
|
|
||||||
from numpy.linalg import norm
|
|
||||||
from matplotlib import pyplot as plt
|
|
||||||
from camera_stage_tracker import Tracker, move_until_motion_detected
|
|
||||||
|
|
||||||
from functools import partial
|
|
||||||
|
|
||||||
|
|
||||||
def backlash_corrected_move(get_position, move, backlash_amount, pos):
|
|
||||||
"""Make two moves, arriving at `pos` from a consistent direction"""
|
|
||||||
displacement = pos - get_position()
|
|
||||||
backlash_vector = (displacement < 0).astype(np.int) * backlash_amount
|
|
||||||
if np.any(backlash_vector > 0):
|
|
||||||
move(pos - backlash_vector)
|
|
||||||
move(pos)
|
|
||||||
|
|
||||||
|
|
||||||
def bake_backlash_corrected_move(get_position, move, backlash_amount):
|
|
||||||
"""Return a function that performs backlash-corrected moves"""
|
|
||||||
return partial(backlash_corrected_move, get_position, move, backlash_amount)
|
|
||||||
|
|
||||||
|
|
||||||
def calibrate_xy_grid(tracker, move, step=100, n_steps=4, backlash_compensation=0):
|
|
||||||
"""Make a series of moves in X and Y to determine the XY components of the pixel-to-sample matrix.
|
|
||||||
|
|
||||||
Arguments:
|
|
||||||
tracker: Tracker
|
|
||||||
An initialised Tracker object, centred on the starting point. This provides position readout from the stage and the camera.
|
|
||||||
move: function
|
|
||||||
A function that accepts a 1D array and performs an absolute move to
|
|
||||||
that position. If backlash correction is needed, include it here.
|
|
||||||
step : float, optional (default 100)
|
|
||||||
The amount to move the stage by. This should move the sample by approximately 1/10th of the field of view.
|
|
||||||
"""
|
|
||||||
try: # Ensure that the tracker has a template set
|
|
||||||
_ = tracker.template
|
|
||||||
except:
|
|
||||||
tracker.acquire_template()
|
|
||||||
tracker.reset_history() # make sure we get rid of the initial (0,0) point
|
|
||||||
starting_position = tracker.get_position()
|
|
||||||
# Move the stage in a square, recording the displacement from both the stage and the camera
|
|
||||||
try:
|
|
||||||
for x in (np.arange(n_steps) - n_steps / 2.0) * step:
|
|
||||||
for y in (np.arange(n_steps) - n_steps / 2.0) * step:
|
|
||||||
move(starting_position + np.array([x, y, 0]))
|
|
||||||
tracker.append_point()
|
|
||||||
finally:
|
|
||||||
move(starting_position)
|
|
||||||
# We then use least-squares to fit the XY part of the matrix relating
|
|
||||||
# pixels to distance
|
|
||||||
# stage_positions should be the stage positions, with a zero mean.
|
|
||||||
# image_positions should be the same, but calculated from the images
|
|
||||||
stage_positions, image_positions = tracker.history
|
|
||||||
stage_positions = stage_positions.astype(np.float)
|
|
||||||
stage_positions -= np.mean(stage_positions, axis=0)
|
|
||||||
stage_positions = stage_positions[:, :2] # ensure it's 2d
|
|
||||||
image_positions -= np.mean(image_positions, axis=0)
|
|
||||||
# image_positions *= -1 # To get the matrix right, we want the position of each
|
|
||||||
# image relative to the template, rather than the other way around
|
|
||||||
A, res, rank, s = np.linalg.lstsq(
|
|
||||||
image_positions, stage_positions
|
|
||||||
) # we solve pixel_shifts*A = location_shifts
|
|
||||||
|
|
||||||
transformed_image_positions = np.dot(image_positions, A)
|
|
||||||
residuals = transformed_image_positions - stage_positions
|
|
||||||
fractional_error = norm(residuals) / stage_positions.shape[0]
|
|
||||||
logging.debug(f"Ratio of residuals to displacement is {fractional_error})")
|
|
||||||
if fractional_error > 0.05: # Check it was a reasonably good fit
|
|
||||||
logging.warning(
|
|
||||||
"Warning: the error fitting measured displacements was %.1f%%"
|
|
||||||
% (fractional_error * 100)
|
|
||||||
)
|
|
||||||
logging.info(
|
|
||||||
f"Calibrated the pixel-location matrix.\nResiduals were {fractional_error*100:.1f}% of the shift."
|
|
||||||
)
|
|
||||||
|
|
||||||
return {
|
|
||||||
"image_to_stage_displacement": A,
|
|
||||||
"moves": (stage_positions, image_positions),
|
|
||||||
"fractional_error": fractional_error,
|
|
||||||
}
|
|
||||||
|
|
@ -1,332 +0,0 @@
|
||||||
"""
|
|
||||||
Camera-stage tracker
|
|
||||||
|
|
||||||
The `Tracker` class in this file is used to simplify code for tasks that involve moving the
|
|
||||||
stage, and tracking the corresponding motion with the camera.
|
|
||||||
|
|
||||||
(c) Richard Bowman 2019, released under GNU GPL v3
|
|
||||||
"""
|
|
||||||
import numpy as np
|
|
||||||
import time
|
|
||||||
from numpy.linalg import norm
|
|
||||||
import cv2
|
|
||||||
from scipy import ndimage
|
|
||||||
|
|
||||||
|
|
||||||
def central_half(image):
|
|
||||||
"""Return the central 50% (in X and Y) of an image"""
|
|
||||||
w, h = image.shape[:2]
|
|
||||||
return image[int(w / 4) : int(3 * w / 4), int(h / 4) : int(3 * h / 4), ...]
|
|
||||||
|
|
||||||
|
|
||||||
def datum_pixel(image):
|
|
||||||
"""Get the datum pixel of an image - if no property is present, assume the central pixel."""
|
|
||||||
try:
|
|
||||||
return np.array(image.datum_pixel)
|
|
||||||
except:
|
|
||||||
return (np.array(image.shape[:2]) - 1) / 2.0
|
|
||||||
|
|
||||||
|
|
||||||
def locate_feature_in_image(image, feature, margin=0, restrict=False):
|
|
||||||
"""Find the given feature (small image) and return the position of its datum (or centre) in the image's pixels.
|
|
||||||
|
|
||||||
image : numpy.array
|
|
||||||
The image in which to look.
|
|
||||||
feature : numpy.array
|
|
||||||
The feature to look for. Ideally should be an `ImageWithLocation`.
|
|
||||||
margin : int (optional)
|
|
||||||
Make sure the feature image is at least this much smaller than the big image. NB this will take account of the
|
|
||||||
image datum points - if the datum points are superimposed, there must be at least margin pixels on each side of
|
|
||||||
the feature image.
|
|
||||||
restrict : bool (optional, default False)
|
|
||||||
If set to true, restrict the search area to a square of (margin * 2 + 1) pixels centred on the pixel that most
|
|
||||||
closely overlaps the datum points of the two images.
|
|
||||||
|
|
||||||
The `image` must be larger than `feature` by a margin big enough to produce a meaningful search area. We use the
|
|
||||||
OpenCV `matchTemplate` method to find the feature. The returned position is the position, relative to the corner of
|
|
||||||
the first image, of the "datum pixel" of the feature image. If no datum pixel is specified, we assume it's the
|
|
||||||
centre of the image. The output of this function can be passed into the pixel_to_location() method of the larger
|
|
||||||
image to yield the position in the sample of the feature you're looking for.
|
|
||||||
"""
|
|
||||||
# The line below is superfluous if we keep the datum-aware code below it.
|
|
||||||
assert (
|
|
||||||
image.shape[0] > feature.shape[0] and image.shape[1] > feature.shape[1]
|
|
||||||
), "Image must be larger than feature!"
|
|
||||||
# Check that there's enough space around the feature image
|
|
||||||
lower_margin = datum_pixel(image) - datum_pixel(feature)
|
|
||||||
upper_margin = (image.shape[:2] - datum_pixel(image)) - (
|
|
||||||
feature.shape[:2] - datum_pixel(feature)
|
|
||||||
)
|
|
||||||
assert np.all(
|
|
||||||
np.array([lower_margin, upper_margin]) >= margin
|
|
||||||
), "The feature image is too large."
|
|
||||||
# TODO: sensible auto-crop of the template if it's too large?
|
|
||||||
image_shift = np.array((0, 0))
|
|
||||||
if restrict:
|
|
||||||
# if requested, crop the larger image so that our search area is (2*margin + 1) square.
|
|
||||||
image_shift = np.array(lower_margin - margin, dtype=int)
|
|
||||||
image = image[
|
|
||||||
image_shift[0] : image_shift[0] + feature.shape[0] + 2 * margin + 1,
|
|
||||||
image_shift[1] : image_shift[1] + feature.shape[1] + 2 * margin + 1,
|
|
||||||
...,
|
|
||||||
]
|
|
||||||
|
|
||||||
corr = cv2.matchTemplate(
|
|
||||||
image, feature, cv2.TM_SQDIFF_NORMED
|
|
||||||
) # correlate them: NB the match position is the MINIMUM
|
|
||||||
corr = -corr # invert the image so we can find a peak
|
|
||||||
corr += (
|
|
||||||
corr.max() - corr.min()
|
|
||||||
) * 0.1 - corr.max() # background-subtract 90% of maximum
|
|
||||||
corr = cv2.threshold(corr, 0, 0, cv2.THRESH_TOZERO)[
|
|
||||||
1
|
|
||||||
] # zero out any negative pixels - but there should always be > 0 nonzero pixels
|
|
||||||
assert (
|
|
||||||
np.sum(corr) > 0
|
|
||||||
), "Error: the correlation image doesn't have any nonzero pixels."
|
|
||||||
peak = ndimage.measurements.center_of_mass(
|
|
||||||
corr
|
|
||||||
) # take the centroid (NB this is of grayscale values, not binary)
|
|
||||||
pos = (
|
|
||||||
np.array(peak) + image_shift + datum_pixel(feature)
|
|
||||||
) # return the position of the feature's datum point.
|
|
||||||
return pos
|
|
||||||
|
|
||||||
|
|
||||||
class Tracker:
|
|
||||||
def __init__(self, grab_image, get_position, settle=None):
|
|
||||||
"""A class to manage moving the stage and following motion in the image
|
|
||||||
|
|
||||||
Constructor Arguments:
|
|
||||||
grab_image: a function that returns the image as a numpy array
|
|
||||||
get_position: a function that returns position as a numpy array
|
|
||||||
settle: a function that waits and/or discards images
|
|
||||||
|
|
||||||
We accept functions because that seems like the easiest way to be
|
|
||||||
compatible with many different cameras/stages. Subclass and override
|
|
||||||
``__init__`` if you want to use a particular object instead.
|
|
||||||
|
|
||||||
NB the ``image_position`` that this class returns may be the negative of
|
|
||||||
what you might expect. This is because normally we are looking for
|
|
||||||
where a certain object (usually matched to a template image) is within
|
|
||||||
an image. Instead, the ``Tracker`` is following motion of the image
|
|
||||||
relative to a template. Our model is that we have a static image on
|
|
||||||
the slide, so we're tracking the slide's motion.
|
|
||||||
"""
|
|
||||||
self._grab_image = grab_image
|
|
||||||
self._get_position = get_position
|
|
||||||
self._settle = settle
|
|
||||||
self._template = None
|
|
||||||
self.margin = np.array([0, 0])
|
|
||||||
self._template_position = np.array([0.0, 0.0])
|
|
||||||
self.image_shape = None
|
|
||||||
|
|
||||||
def get_position(self):
|
|
||||||
"""Get the position of the stage"""
|
|
||||||
return np.array(self._get_position())
|
|
||||||
|
|
||||||
def settle(self):
|
|
||||||
"""Wait a short time and discard an image so the stage is no longer wobbling."""
|
|
||||||
if self._settle is not None:
|
|
||||||
self._settle()
|
|
||||||
else:
|
|
||||||
time.sleep(0.3)
|
|
||||||
self._grab_image()
|
|
||||||
|
|
||||||
@property
|
|
||||||
def template(self):
|
|
||||||
"""The template image (should be a numpy array)"""
|
|
||||||
if self._template is None:
|
|
||||||
raise ValueError("Attempt to use the tracker before setting the template")
|
|
||||||
else:
|
|
||||||
return self._template
|
|
||||||
|
|
||||||
@template.setter
|
|
||||||
def template(self, new_value):
|
|
||||||
self._template = new_value
|
|
||||||
|
|
||||||
def acquire_template(
|
|
||||||
self, settle=True, reset_history=True, relative_positions=True
|
|
||||||
):
|
|
||||||
"""Take a new image, and use it as the template. NB this records the initial point.
|
|
||||||
|
|
||||||
We will wait for the stage to settle, then acquire a new image to use as the template.
|
|
||||||
Immediately afterwards, we record the first point, so we will acquire a second image
|
|
||||||
and also read the stage's position.
|
|
||||||
|
|
||||||
The template image will be the central 50% of the starting image, which means the
|
|
||||||
maximum displacement will be 0.25 fields-of-view in all directions.
|
|
||||||
|
|
||||||
Arguments:
|
|
||||||
settle: bool, default True
|
|
||||||
Whether to wait for the stage to settle before taking the template image
|
|
||||||
reset_history: bool, default True
|
|
||||||
Whether to erase all the previously-stored positions
|
|
||||||
relative_positions: bool, default True
|
|
||||||
If true, we will define the first point (as read from the camera) to be [0,0]
|
|
||||||
and make all future measurements relative to this one. NB this won't affect
|
|
||||||
the stage positions, which are always absolute.
|
|
||||||
"""
|
|
||||||
if settle:
|
|
||||||
self.settle()
|
|
||||||
image = self._grab_image()
|
|
||||||
self.template = central_half(image)
|
|
||||||
self.image_shape = image.shape
|
|
||||||
self.margin = np.array(image.shape)[:2] - np.array(self.template.shape)[:2]
|
|
||||||
if reset_history:
|
|
||||||
self.reset_history()
|
|
||||||
self._template_position = np.array([0.0, 0.0])
|
|
||||||
if relative_positions:
|
|
||||||
self._template_position = self.track_image(
|
|
||||||
image
|
|
||||||
) # Position should be zero initially
|
|
||||||
self.append_point(settle=False)
|
|
||||||
|
|
||||||
@property
|
|
||||||
def max_displacement(self):
|
|
||||||
"""The highest position values that can be tracked"""
|
|
||||||
return (
|
|
||||||
self.margin // 2
|
|
||||||
) # TODO: be cleverer about non-trivial values of template_position
|
|
||||||
|
|
||||||
@property
|
|
||||||
def min_displacement(self):
|
|
||||||
"""The lowest position values that can be tracked"""
|
|
||||||
return (
|
|
||||||
-self.max_displacement
|
|
||||||
) # TODO: be cleverer about non-trivial template_position values
|
|
||||||
|
|
||||||
@property
|
|
||||||
def max_safe_displacement(self):
|
|
||||||
"""The biggest displacement we can safely attempt to track without knowing direction."""
|
|
||||||
return np.min(np.concatenate([self.max_displacement, -self.min_displacement]))
|
|
||||||
|
|
||||||
def track_image(self, image):
|
|
||||||
"""Find the position of the image relative to the template
|
|
||||||
|
|
||||||
NB this class is intended to track motion of the sample - most of
|
|
||||||
the time, we're interested in the motion of a (small) object that
|
|
||||||
is represented by the template, relative to the (larger) image. In
|
|
||||||
our case, we're doing the opposite - tracking motion of the image,
|
|
||||||
relative to a picture of part of the sample. That's why there is
|
|
||||||
a minus sign in front of `locate_feature_in_image` in the source
|
|
||||||
code.
|
|
||||||
"""
|
|
||||||
return -locate_feature_in_image(image, self.template) - self._template_position
|
|
||||||
|
|
||||||
def append_point(self, settle=True, image=None):
|
|
||||||
"""Find the current position using both stage and image, and append it"""
|
|
||||||
if settle:
|
|
||||||
self.settle()
|
|
||||||
if image is None:
|
|
||||||
image = self._grab_image()
|
|
||||||
image_pos = self.track_image(image)
|
|
||||||
stage_pos = self.get_position()
|
|
||||||
self._image_positions.append(image_pos)
|
|
||||||
self._stage_positions.append(stage_pos)
|
|
||||||
return stage_pos, image_pos
|
|
||||||
|
|
||||||
@property
|
|
||||||
def stage_positions(self):
|
|
||||||
"""An array of positions we have moved the stage to"""
|
|
||||||
return np.array(self._stage_positions)
|
|
||||||
|
|
||||||
@property
|
|
||||||
def image_positions(self):
|
|
||||||
"""An array of positions we have moved the stage to"""
|
|
||||||
return np.array(self._image_positions)
|
|
||||||
|
|
||||||
@property
|
|
||||||
def history(self):
|
|
||||||
"""Return arrays of stage, image positions"""
|
|
||||||
return self.stage_positions, self.image_positions
|
|
||||||
|
|
||||||
def reset_history(self, leave_first_point=False):
|
|
||||||
"""Reset the positions and displacements recorded"""
|
|
||||||
if leave_first_point:
|
|
||||||
self._stage_positions = [self._stage_positions[0]]
|
|
||||||
self._image_positions = [self._image_positions[0]]
|
|
||||||
else:
|
|
||||||
self._stage_positions = []
|
|
||||||
self._image_positions = []
|
|
||||||
|
|
||||||
@property
|
|
||||||
def moving_away_from_centre(self):
|
|
||||||
"""Whether we are moving away from [0,0] on the camera.
|
|
||||||
|
|
||||||
If we have recorded more than two steps, this property will be
|
|
||||||
`True` if the most recent point in the history is farther away from
|
|
||||||
`[0,0]` than the second most recent point. If we have recorded fewer
|
|
||||||
than 2 points, this property returns None.
|
|
||||||
"""
|
|
||||||
if len(self.image_positions) < 2:
|
|
||||||
return None
|
|
||||||
else:
|
|
||||||
return norm(self.image_positions[-1, :]) > norm(self.image_positions[-2])
|
|
||||||
|
|
||||||
|
|
||||||
def move_until_motion_detected(
|
|
||||||
tracker,
|
|
||||||
move,
|
|
||||||
displacement,
|
|
||||||
threshold=10,
|
|
||||||
multipliers=2 ** np.arange(16),
|
|
||||||
detect_cumulative_motion=False,
|
|
||||||
):
|
|
||||||
"""Move the stage until we can detect motion in the camera.
|
|
||||||
|
|
||||||
We move the stage in the direction given by ``displacement`` until the
|
|
||||||
image has shifted by at least ``threshold`` pixels. The steps will be
|
|
||||||
given by ``multipliers``, i.e. each time we move to
|
|
||||||
``displacement * multipliers[i]`` relative to the starting position.
|
|
||||||
|
|
||||||
NB we expect that the ``tracker`` object has already been initialised
|
|
||||||
with ``acquire_template``.
|
|
||||||
|
|
||||||
``detect_cumulative_motion`` will use the first point in the tracker as
|
|
||||||
the point to detect displacement relative to, rather than the last point.
|
|
||||||
The displacements are always made relative to the last point in the tracker
|
|
||||||
as it is passed in (i.e. the stage is always moved relative to where it
|
|
||||||
currently is, but motion detection may be done relative to where the tracker
|
|
||||||
was initialised). This only matters if the tracker has more than one point
|
|
||||||
in its history.
|
|
||||||
|
|
||||||
The return value `i, m` is the number of moves made, and the largest
|
|
||||||
multiplier value that was used, i.e. we moved by a total of
|
|
||||||
`displacement * m`.
|
|
||||||
"""
|
|
||||||
displacement = np.array(displacement)
|
|
||||||
starting_image_position = tracker.image_positions[
|
|
||||||
0 if detect_cumulative_motion else -1, :
|
|
||||||
]
|
|
||||||
starting_stage_position = tracker.stage_positions[-1, :]
|
|
||||||
for i, m in enumerate(multipliers):
|
|
||||||
move(starting_stage_position + displacement * m)
|
|
||||||
tracker.append_point()
|
|
||||||
if norm(tracker.image_positions[-1, :] - starting_image_position) >= threshold:
|
|
||||||
return i + 1, m
|
|
||||||
raise Exception(
|
|
||||||
"Moved the stage by {} but saw no motion.".format(
|
|
||||||
multipliers[-1] * displacement
|
|
||||||
)
|
|
||||||
)
|
|
||||||
|
|
||||||
|
|
||||||
def concatenate_tracker_histories(histories):
|
|
||||||
"""Combine a number of separate tracker history entries into one
|
|
||||||
|
|
||||||
A "tracker history" refers to the output of `Tracker.history`, i.e.
|
|
||||||
it is a tuple of `(stage_positions, image_positions)` with the two
|
|
||||||
components being a Nx3 and Nx2 `numpy.ndarray` objects respectively.
|
|
||||||
|
|
||||||
Given an array of such tuples, we will concatenate the components,
|
|
||||||
returning a single "tracker history" with the segments concatenated.
|
|
||||||
|
|
||||||
Returns: backlash, pixels_per_step, fractional_error
|
|
||||||
|
|
||||||
The return value is a tuple of 3 numbers; the estimated backlash (in
|
|
||||||
motor steps), the ratio of image_position changes to stage_position
|
|
||||||
(in units of pixels/steps), and an estimate of goodness of fit.
|
|
||||||
"""
|
|
||||||
components = zip(*histories)
|
|
||||||
return tuple(np.concatenate(c, axis=1) for c in components)
|
|
||||||
|
|
@ -1,207 +0,0 @@
|
||||||
"""
|
|
||||||
API extension for stage calibration
|
|
||||||
|
|
||||||
This file contains the HTTP API for camera/stage calibration.
|
|
||||||
"""
|
|
||||||
from labthings.server.view import View, ActionView, PropertyView
|
|
||||||
from labthings.server.find import find_component
|
|
||||||
from labthings.server.extensions import BaseExtension
|
|
||||||
from labthings.server.decorators import (
|
|
||||||
ThingAction,
|
|
||||||
use_args,
|
|
||||||
ThingProperty,
|
|
||||||
)
|
|
||||||
from labthings.server import fields
|
|
||||||
|
|
||||||
from labthings.core.utilities import get_by_path, set_by_path, create_from_path
|
|
||||||
|
|
||||||
|
|
||||||
from flask import abort, send_file
|
|
||||||
|
|
||||||
import logging
|
|
||||||
import time
|
|
||||||
import numpy as np
|
|
||||||
import PIL
|
|
||||||
import io
|
|
||||||
import os
|
|
||||||
import json
|
|
||||||
|
|
||||||
from .camera_stage_calibration_1d import (
|
|
||||||
calibrate_backlash_1d,
|
|
||||||
image_to_stage_displacement_from_1d,
|
|
||||||
)
|
|
||||||
from .camera_stage_tracker import Tracker
|
|
||||||
|
|
||||||
from openflexure_microscope.utilities import axes_to_array
|
|
||||||
from openflexure_microscope.paths import data_file_path
|
|
||||||
from openflexure_microscope.config import JSONEncoder
|
|
||||||
|
|
||||||
CSM_DATAFILE_NAME = "csm_calibration.json"
|
|
||||||
CSM_DATAFILE_PATH = data_file_path(CSM_DATAFILE_NAME)
|
|
||||||
|
|
||||||
|
|
||||||
class CSMExtension(BaseExtension):
|
|
||||||
"""
|
|
||||||
Use the camera as an encoder, so we can relate camera and stage coordinates
|
|
||||||
"""
|
|
||||||
|
|
||||||
def __init__(self):
|
|
||||||
BaseExtension.__init__(
|
|
||||||
self, "org.openflexure.camera_stage_mapping", version="0.0.1"
|
|
||||||
)
|
|
||||||
|
|
||||||
_microscope = None
|
|
||||||
|
|
||||||
@property
|
|
||||||
def microscope(self):
|
|
||||||
# TODO: does caching the microscope actually help?
|
|
||||||
if self._microscope is None:
|
|
||||||
self._microscope = find_component("org.openflexure.microscope")
|
|
||||||
return self._microscope
|
|
||||||
|
|
||||||
def update_settings(self, settings):
|
|
||||||
"""Update the stored extension settings dictionary"""
|
|
||||||
keys = ["extensions", self.name]
|
|
||||||
dictionary = create_from_path(keys)
|
|
||||||
set_by_path(dictionary, keys, settings)
|
|
||||||
logging.info(f"Updating settings with {dictionary}")
|
|
||||||
self.microscope.update_settings(dictionary)
|
|
||||||
self.microscope.save_settings()
|
|
||||||
|
|
||||||
def get_settings(self):
|
|
||||||
"""Retrieve the settings for this extension"""
|
|
||||||
keys = ["extensions", self.name]
|
|
||||||
return get_by_path(self.microscope.read_settings(), keys)
|
|
||||||
|
|
||||||
def camera_stage_functions(self):
|
|
||||||
"""Return functions that allow us to interface with the microscope"""
|
|
||||||
self.microscope.camera.start_worker() # ensure the worker thread is running, so there is an MJPEG stream
|
|
||||||
|
|
||||||
def grab_image():
|
|
||||||
jpeg = self.microscope.camera.get_frame()
|
|
||||||
return np.array(PIL.Image.open(io.BytesIO(jpeg)))
|
|
||||||
|
|
||||||
def get_position():
|
|
||||||
return self.microscope.stage.position
|
|
||||||
|
|
||||||
move = self.microscope.stage.move_abs
|
|
||||||
|
|
||||||
return grab_image, get_position, move
|
|
||||||
|
|
||||||
def calibrate_1d(self, direction):
|
|
||||||
"""Move a microscope's stage in 1D, and figure out the relationship with the camera"""
|
|
||||||
grab_image, get_position, move = self.camera_stage_functions()
|
|
||||||
|
|
||||||
def wait():
|
|
||||||
time.sleep(0.2)
|
|
||||||
|
|
||||||
tracker = Tracker(grab_image, get_position, settle=wait)
|
|
||||||
|
|
||||||
return calibrate_backlash_1d(tracker, move, direction)
|
|
||||||
|
|
||||||
def calibrate_xy(self):
|
|
||||||
"""Move the microscope's stage in X and Y, to calibrate its relationship to the camera"""
|
|
||||||
logging.info("Calibrating X axis:")
|
|
||||||
cal_x = self.calibrate_1d(np.array([1, 0, 0]))
|
|
||||||
logging.info("Calibrating Y axis:")
|
|
||||||
cal_y = self.calibrate_1d(np.array([0, 1, 0]))
|
|
||||||
|
|
||||||
# Combine X and Y calibrations to make a 2D calibration
|
|
||||||
cal_xy = image_to_stage_displacement_from_1d([cal_x, cal_y])
|
|
||||||
self.update_settings(cal_xy)
|
|
||||||
|
|
||||||
data = {
|
|
||||||
"camera_stage_mapping_calibration": cal_xy,
|
|
||||||
"linear_calibration_x": cal_x,
|
|
||||||
"linear_calibration_y": cal_y,
|
|
||||||
}
|
|
||||||
|
|
||||||
with open(CSM_DATAFILE_PATH, "w") as f:
|
|
||||||
json.dump(data, f, cls=JSONEncoder)
|
|
||||||
|
|
||||||
return data
|
|
||||||
|
|
||||||
@property
|
|
||||||
def image_to_stage_displacement_matrix(self):
|
|
||||||
try:
|
|
||||||
settings = self.get_settings()
|
|
||||||
return settings["image_to_stage_displacement"]
|
|
||||||
except KeyError:
|
|
||||||
raise ValueError("The microscope has not yet been calibrated.")
|
|
||||||
|
|
||||||
def move_in_image_coordinates(self, displacement_in_pixels):
|
|
||||||
"""Move by a given number of pixels on the camera"""
|
|
||||||
p = np.array(displacement_in_pixels)
|
|
||||||
relative_move = np.dot(p, self.image_to_stage_displacement_matrix)
|
|
||||||
self.microscope.stage.move_rel([relative_move[0], relative_move[1], 0])
|
|
||||||
|
|
||||||
|
|
||||||
csm_extension = CSMExtension()
|
|
||||||
|
|
||||||
|
|
||||||
class Calibrate1DView(ActionView):
|
|
||||||
@use_args(
|
|
||||||
{"direction": fields.List(fields.Float(), required=True, example=[1, 0, 0])}
|
|
||||||
)
|
|
||||||
def post(self, args):
|
|
||||||
"""Calibrate one axis of the microscope stage against the camera."""
|
|
||||||
|
|
||||||
direction = np.array(args.get("direction"))
|
|
||||||
|
|
||||||
return csm_extension.calibrate_1d(direction)
|
|
||||||
|
|
||||||
|
|
||||||
csm_extension.add_view(Calibrate1DView, "/calibrate_1d", endpoint="calibrate_1d")
|
|
||||||
|
|
||||||
|
|
||||||
class CalibrateXYView(ActionView):
|
|
||||||
def post(self):
|
|
||||||
"""Calibrate both axes of the microscope stage against the camera."""
|
|
||||||
return csm_extension.calibrate_xy()
|
|
||||||
|
|
||||||
|
|
||||||
csm_extension.add_view(CalibrateXYView, "/calibrate_xy", endpoint="calibrate_xy")
|
|
||||||
|
|
||||||
|
|
||||||
class MoveInImageCoordinatesView(ActionView):
|
|
||||||
@use_args(
|
|
||||||
{
|
|
||||||
"x": fields.Float(
|
|
||||||
description="The number of pixels to move in X",
|
|
||||||
required=True,
|
|
||||||
example=100,
|
|
||||||
),
|
|
||||||
"y": fields.Float(
|
|
||||||
description="The number of pixels to move in Y",
|
|
||||||
required=True,
|
|
||||||
example=100,
|
|
||||||
),
|
|
||||||
}
|
|
||||||
)
|
|
||||||
def post(self, args):
|
|
||||||
logging.debug("moving in pixels")
|
|
||||||
"""Move the microscope stage, such that we move by a given number of pixels on the camera"""
|
|
||||||
csm_extension.move_in_image_coordinates(
|
|
||||||
np.array([args.get("x"), args.get("y")])
|
|
||||||
)
|
|
||||||
|
|
||||||
return csm_extension.microscope.state["stage"]["position"]
|
|
||||||
|
|
||||||
|
|
||||||
csm_extension.add_view(MoveInImageCoordinatesView, "/move_in_image_coordinates", endpoint="move_in_image_coordinates")
|
|
||||||
|
|
||||||
|
|
||||||
class GetCalibrationFile(PropertyView):
|
|
||||||
def get(self):
|
|
||||||
"""Get the calibration data in JSON format."""
|
|
||||||
datafile_name = CSM_DATAFILE_NAME
|
|
||||||
datafile_path = CSM_DATAFILE_PATH
|
|
||||||
|
|
||||||
if os.path.isfile(datafile_path):
|
|
||||||
with open(datafile_path, "rb") as f:
|
|
||||||
return json.load(f)
|
|
||||||
else:
|
|
||||||
return {}
|
|
||||||
|
|
||||||
|
|
||||||
csm_extension.add_view(GetCalibrationFile, "/get_calibration", endpoint="get_calibration")
|
|
||||||
|
|
@ -1,270 +0,0 @@
|
||||||
"""
|
|
||||||
Image With Location
|
|
||||||
===================
|
|
||||||
|
|
||||||
This datatype supports the various operations that rely on linking a camera to a microscope stage. It is an image
|
|
||||||
along with the metadata required to relate positions in the image to positions in real life.
|
|
||||||
|
|
||||||
To create an `ImageWithLocation`, first put the image data into an `ArrayWithAttrs` and then specify the required
|
|
||||||
metadata in the `attrs` dictionary. The `pixel_to_sample_matrix` is the only required piece of metadata - the
|
|
||||||
`datum_pixel` is optional (and if missing, will be assumed to be the central pixel).
|
|
||||||
|
|
||||||
A note on coordinate systems
|
|
||||||
----------------------------
|
|
||||||
I've tried to stick to two coordinate systems: that used by the stage, generally called a "location", and pixels in an
|
|
||||||
image.
|
|
||||||
|
|
||||||
Images have a "datum pixel", specified in metadata or assumed to be the centre (i.e. pixel (N-1)/2 for a width of N).
|
|
||||||
This need not be an integer pixel position, but is specified in pixels relative to the [0,0] pixel. When considering
|
|
||||||
something within an image, *the coordinate system is always relative to pixel [0,0]*, not relative to the datum pixel.
|
|
||||||
Similarly, the transformation matrix that moves between pixel and stage coordinates uses [0,0] as its origin, not the
|
|
||||||
datum pixel. However, when considering the displacement between two images, this is usually with respect to the datum
|
|
||||||
pixels of the images - though we should generally specify this.
|
|
||||||
|
|
||||||
We transform between pixel and location coordinate systems with a matrix, the `pixel_to_sample_matrix`. Usually it
|
|
||||||
is called ``M`` in mathematical expressions. To convert a pixel coordinate to a location, we post-multiply the pixel
|
|
||||||
coordinate by the matrix, i.e. ``l = p.M`` and to convert the other way we use the inverse of ``M`` so ``p = l.M``
|
|
||||||
where the dot denotes matrix multiplication using `numpy.dot`.
|
|
||||||
|
|
||||||
Note that the calibration matrix is a 4x4 matrix, and the vectors should be ``(x, y, z, 1)`` so that we encode
|
|
||||||
the absolute position in that matrix, along with scaling and rotation. It's an entirely sensible thing to include
|
|
||||||
the stage coordinates in metadata as well as the matrix, but it is not needed.
|
|
||||||
|
|
||||||
"""
|
|
||||||
from __future__ import division
|
|
||||||
|
|
||||||
from builtins import range
|
|
||||||
from past.utils import old_div
|
|
||||||
import numpy as np
|
|
||||||
from array_with_attrs import ArrayWithAttrs, ensure_attrs
|
|
||||||
import cv2
|
|
||||||
|
|
||||||
# import cv2.cv
|
|
||||||
from scipy import ndimage
|
|
||||||
|
|
||||||
|
|
||||||
class ImageWithLocation(ArrayWithAttrs):
|
|
||||||
"""An image, as a numpy array, with attributes to provide location information
|
|
||||||
|
|
||||||
This is a functioning `numpy.ndarray` which can store the image in uncompressed format.
|
|
||||||
We require that the `attrs` dictionary (defined by `ArrayWithAttrs`) contains keys
|
|
||||||
that we use to store the crucial mapping from pixels in the image to position in the
|
|
||||||
sample.
|
|
||||||
"""
|
|
||||||
|
|
||||||
# def __array_finalize__(self, obj):
|
|
||||||
# """Ensure that the object is a properly set-up ImageWithLocation"""
|
|
||||||
# ArrayWithAttrs.__array_finalize__(self, obj) # Ensure we have self.attrs
|
|
||||||
def __getitem__(self, item):
|
|
||||||
"""Update the metadata when we extract a slice"""
|
|
||||||
try:
|
|
||||||
# Handle specially the case where we are extracting a 2D region of the image, i.e. the first and second
|
|
||||||
# indices are slices. We test for that here - and do it in a try: except block so that if, for example,
|
|
||||||
# item is not indexable,
|
|
||||||
assert isinstance(item[0], slice), "First index was not a slice"
|
|
||||||
assert isinstance(item[1], slice), "Second index was not a slice"
|
|
||||||
start = np.array([item[i].start for i in range(2)])
|
|
||||||
start = np.where(
|
|
||||||
start == np.array(None), 0, start
|
|
||||||
) # missing start points are equivalent to zero
|
|
||||||
step = np.array([item[i].step for i in range(2)])
|
|
||||||
step = np.where(
|
|
||||||
step == np.array(None), 1, step
|
|
||||||
) # missing step is equivalent to step==1
|
|
||||||
except:
|
|
||||||
# If the above doesn't work, assume we're not dealing with a 2D slice and give up.
|
|
||||||
return super(ImageWithLocation, self).__getitem__(item) # pass it on up
|
|
||||||
|
|
||||||
out = super(ImageWithLocation, self).__getitem__(item) # retrieve the slice
|
|
||||||
out.datum_pixel -= (
|
|
||||||
start
|
|
||||||
) # adjust the datum pixel so it refers to the same part of the image
|
|
||||||
# Next, we adjust the constant part of the pixel-sample matrix so pixels stay in the same place
|
|
||||||
location_shift = np.dot(ensure_3d(start), self.pixel_to_sample_matrix[:3, :3])
|
|
||||||
out.pixel_to_sample_matrix[3, :3] += location_shift
|
|
||||||
if not np.all(step == 1):
|
|
||||||
# if we're downsampling, remember to scale datum_pixel accordingly
|
|
||||||
out.datum_pixel = old_div(out.datum_pixel, step)
|
|
||||||
# Scale the pixel-to-sample matrix if we've got a non-unity step in the slice
|
|
||||||
# I don't understand why I can't do this with slicing, but it all goes wrong...
|
|
||||||
for i in range(2):
|
|
||||||
out.pixel_to_sample_matrix[i, :3] *= step[i]
|
|
||||||
return out
|
|
||||||
|
|
||||||
def pixel_to_location(self, pixel):
|
|
||||||
"""Return the location in the sample of the given pixel.
|
|
||||||
|
|
||||||
NB this returns a 3D location, including Z."""
|
|
||||||
p = ensure_2d(pixel)
|
|
||||||
l = np.dot(np.array([p[0], p[1], 0, 1]), self.pixel_to_sample_matrix)
|
|
||||||
return l[:3]
|
|
||||||
|
|
||||||
def location_to_pixel(self, location, check_bounds=False, z_tolerance=np.infty):
|
|
||||||
"""Return the pixel coordinates of a given location in the sample.
|
|
||||||
|
|
||||||
location : numpy.ndarray
|
|
||||||
A 2- or 3- element numpy array representing sample position, in units of distance.
|
|
||||||
check_bounds : bool, optional (default False)
|
|
||||||
If this is True, raise an exception if the pixel is not in the image.
|
|
||||||
z_tolerance : float, optional (defaults to infinity)
|
|
||||||
If we are checking the bounds, make sure the sample location is within this distance of the image's Z
|
|
||||||
position. The default is to allow any distance.
|
|
||||||
|
|
||||||
Returns : numpy.ndarray
|
|
||||||
A 2- or 3- element position, to match the size of location passed in.
|
|
||||||
"""
|
|
||||||
l = ensure_2d(location)
|
|
||||||
l = l[:2] - self.pixel_to_sample_matrix[3, :2]
|
|
||||||
p = np.dot(l, np.linalg.inv(self.pixel_to_sample_matrix[:2, :2]))
|
|
||||||
if check_bounds:
|
|
||||||
assert np.all(0 <= p[0:2]), "The location was not within the image"
|
|
||||||
assert np.all(
|
|
||||||
p[0:2] <= self.shape[0:2]
|
|
||||||
), "The location was not within the image"
|
|
||||||
assert (
|
|
||||||
np.abs(p[2]) < z_tolerance
|
|
||||||
), "The location was too far away from the plane of the image"
|
|
||||||
if len(location) == 2:
|
|
||||||
return p[:2]
|
|
||||||
else:
|
|
||||||
return p[:3]
|
|
||||||
|
|
||||||
def feature_at(self, centre_position, size=(100, 100), set_datum_to_centre=True):
|
|
||||||
"""Return a thumbnail cropped out of this image, centred on a particular pixel position.
|
|
||||||
|
|
||||||
This is simply a convenience method that saves typing over the usual slice syntax. Below are two equivalent
|
|
||||||
ways of extracting a thumbnail:
|
|
||||||
pos = (240,320)
|
|
||||||
size = (100,100)
|
|
||||||
thumbnail = image[pos[0] - size[0]/2:pos[0] + size[0]/2, pos[1] - size[1]/2:pos[1] + size[1]/2, ...]
|
|
||||||
thumbnail2 = image.feature_at(pos, size)
|
|
||||||
thumbnail3 = image[190:290 270:370]
|
|
||||||
|
|
||||||
``centre_position`` and ``size`` should be two-element tuples, but the intention is that this code will cope
|
|
||||||
gracefully with floating-point values.
|
|
||||||
|
|
||||||
NB the datum pixel of the returned image will be set to its centre, not the datum position of the original image
|
|
||||||
by default. Give the argument ``set_datum_to_centre=False`` to disable this behaviour.
|
|
||||||
"""
|
|
||||||
try:
|
|
||||||
float(centre_position[0])
|
|
||||||
float(centre_position[1])
|
|
||||||
float(size[0])
|
|
||||||
float(size[1])
|
|
||||||
except:
|
|
||||||
raise IndexError(
|
|
||||||
"Error: arguments of feature_at were invalid: {}, {}".format(
|
|
||||||
centre_position, size
|
|
||||||
)
|
|
||||||
)
|
|
||||||
pos = centre_position
|
|
||||||
|
|
||||||
# For now, rely on numpy to complain if the feature is outside the image. May do bound-checking at some point.
|
|
||||||
# If so, we might need to think carefully about the datum pixel of the resulting image.
|
|
||||||
thumb = self[
|
|
||||||
pos[0] - old_div(size[0], 2) : pos[0] + old_div(size[0], 2),
|
|
||||||
pos[1] - old_div(size[1], 2) : pos[1] + old_div(size[1], 2),
|
|
||||||
...,
|
|
||||||
]
|
|
||||||
if set_datum_to_centre:
|
|
||||||
thumb.datum_pixel = (
|
|
||||||
old_div(size[0], 2),
|
|
||||||
old_div(size[1], 2),
|
|
||||||
) # Make the datum point of the new image its centre.
|
|
||||||
return thumb
|
|
||||||
|
|
||||||
def downsample(self, n):
|
|
||||||
"""Return a view of the image, downsampled (sliced with a non-unity step).
|
|
||||||
|
|
||||||
In the future, an optional argument to this function may take means of blocks of the images to improve signal
|
|
||||||
to noise. Currently it just decimates (i.e. throws away rows and columns).
|
|
||||||
"""
|
|
||||||
assert n > 0, "The downsampling factor must be an integer greater than 0"
|
|
||||||
return self[
|
|
||||||
:: int(n), :: int(n), ...
|
|
||||||
] # The slicing code handles updating metadata
|
|
||||||
|
|
||||||
@property
|
|
||||||
def datum_pixel(self):
|
|
||||||
"""The pixel that nominally corresponds to where the image "is".
|
|
||||||
|
|
||||||
Usually the datum pixel is the central pixel, and if the metadata required is not present,
|
|
||||||
we will silently assume that this is the case.
|
|
||||||
"""
|
|
||||||
datum = self.attrs.get(
|
|
||||||
"datum_pixel", old_div((np.array(self.shape[:2]) - 1), 2)
|
|
||||||
)
|
|
||||||
assert len(datum) == 2, "The datum pixel didn't have length 2!"
|
|
||||||
return datum
|
|
||||||
|
|
||||||
@datum_pixel.setter
|
|
||||||
def datum_pixel(self, datum):
|
|
||||||
assert len(datum) == 2, "The datum pixel didn't have length 2!"
|
|
||||||
self.attrs["datum_pixel"] = datum
|
|
||||||
|
|
||||||
@property
|
|
||||||
def datum_location(self):
|
|
||||||
"""The location in the sample of the datum pixel"""
|
|
||||||
return self.pixel_to_location(self.datum_pixel)
|
|
||||||
|
|
||||||
@property
|
|
||||||
def pixel_to_sample_matrix(self):
|
|
||||||
"""The matrix that maps from pixel coordinates to sample coordinates.
|
|
||||||
|
|
||||||
np.dot(p, M) yields a location for the given pixel, where p is [x,y,0,1] and M is this matrix. The location
|
|
||||||
given will be 4 elements long, and will have 1 as the final element.
|
|
||||||
"""
|
|
||||||
M = self.attrs["pixel_to_sample_matrix"]
|
|
||||||
assert M.shape == (4, 4), "The pixel-to-sample matrix is the wrong shape!"
|
|
||||||
assert M.dtype.kind == "f", "The pixel-to-sample matrix is not floating point!"
|
|
||||||
return M
|
|
||||||
|
|
||||||
@pixel_to_sample_matrix.setter
|
|
||||||
def pixel_to_sample_matrix(self, M):
|
|
||||||
M = np.asanyarray(M) # ensure it's an ndarray subclass
|
|
||||||
assert M.shape == (4, 4), "The pixel-to-sample matrix must be 4x4!"
|
|
||||||
assert M.dtype.kind == "f", "The pixel-to-sample matrix must be floating point!"
|
|
||||||
self.attrs["pixel_to_sample_matrix"] = M
|
|
||||||
|
|
||||||
# TODO: split the data type out of this module and put it somewhere sensible
|
|
||||||
|
|
||||||
|
|
||||||
def add_location_metadata(image, pixel_to_sample_matrix, datum_pixel=None):
|
|
||||||
"""Wrap an image if needed, and set its pixel to sample matrix."""
|
|
||||||
awa = ensure_attrs(image) # if needed, convert the image to an ArrayWithAttrs
|
|
||||||
awa.attrs["pixel_to_sample_matrix"] = pixel_to_sample_matrix
|
|
||||||
if datum_pixel is not None:
|
|
||||||
awa.attrs["datum_pixel"] = datum_pixel
|
|
||||||
return awa
|
|
||||||
|
|
||||||
|
|
||||||
def datum_pixel(image):
|
|
||||||
"""Get the datum pixel of an image - if no property is present, assume the central pixel."""
|
|
||||||
try:
|
|
||||||
return np.array(image.datum_pixel)
|
|
||||||
except:
|
|
||||||
return (np.array(image.shape[:2]) - 1) / 2.0
|
|
||||||
|
|
||||||
|
|
||||||
def ensure_3d(vector):
|
|
||||||
"""Make sure a vector has 3 elements, appending a zero if needed."""
|
|
||||||
if len(vector) == 3:
|
|
||||||
return np.array(vector)
|
|
||||||
elif len(vector) == 2:
|
|
||||||
return np.array([vector[0], vector[1], 0])
|
|
||||||
else:
|
|
||||||
raise ValueError(
|
|
||||||
"Tried to ensure a vector was 3D, but it had neither 2 nor 3 elements!"
|
|
||||||
)
|
|
||||||
|
|
||||||
|
|
||||||
def ensure_2d(vector):
|
|
||||||
"""Make sure a vector has 3 elements, appending a zero if needed."""
|
|
||||||
if len(vector) == 2:
|
|
||||||
return np.array(vector)
|
|
||||||
elif len(vector) == 3:
|
|
||||||
return np.array(vector[:2])
|
|
||||||
else:
|
|
||||||
raise ValueError(
|
|
||||||
"Tried to ensure a vector was 2D, but it had neither 2 nor 3 elements!"
|
|
||||||
)
|
|
||||||
24
poetry.lock
generated
24
poetry.lock
generated
|
|
@ -89,6 +89,27 @@ toml = ">=0.9.4"
|
||||||
[package.extras]
|
[package.extras]
|
||||||
d = ["aiohttp (>=3.3.2)"]
|
d = ["aiohttp (>=3.3.2)"]
|
||||||
|
|
||||||
|
[[package]]
|
||||||
|
category = "main"
|
||||||
|
description = ""
|
||||||
|
name = "camera-stage-mapping"
|
||||||
|
optional = false
|
||||||
|
python-versions = "^3.6"
|
||||||
|
version = "0.1.0"
|
||||||
|
|
||||||
|
[package.dependencies]
|
||||||
|
numpy = "^1.17"
|
||||||
|
|
||||||
|
[package.extras]
|
||||||
|
all = ["opencv-python-headless (^4.1)", "scipy (^1.4)"]
|
||||||
|
correlation = ["opencv-python-headless (^4.1)", "scipy (^1.4)"]
|
||||||
|
ipython = []
|
||||||
|
ofmclient = []
|
||||||
|
|
||||||
|
[package.source]
|
||||||
|
reference = "297262362ae208f56847c823838c1d90675ef431"
|
||||||
|
type = "git"
|
||||||
|
url = "https://gitlab.com/openflexure/microscope-extensions/camera-stage-mapping.git"
|
||||||
[[package]]
|
[[package]]
|
||||||
category = "main"
|
category = "main"
|
||||||
description = "Pure Python CBOR (de)serializer with extensive tag support"
|
description = "Pure Python CBOR (de)serializer with extensive tag support"
|
||||||
|
|
@ -793,7 +814,7 @@ testing = ["coverage (>=5.0.3)", "zope.event", "zope.testing"]
|
||||||
rpi = ["RPi.GPIO"]
|
rpi = ["RPi.GPIO"]
|
||||||
|
|
||||||
[metadata]
|
[metadata]
|
||||||
content-hash = "1c923fc71ff343293798e1bd187df98791716e15e4d970fef1a8d77bf25de1e2"
|
content-hash = "38a57e780d968651d6bac204d860866ded059b4169814bc4fe00810d84a79c7f"
|
||||||
python-versions = "^3.6"
|
python-versions = "^3.6"
|
||||||
|
|
||||||
[metadata.files]
|
[metadata.files]
|
||||||
|
|
@ -825,6 +846,7 @@ black = [
|
||||||
{file = "black-18.9b0-py36-none-any.whl", hash = "sha256:817243426042db1d36617910df579a54f1afd659adb96fc5032fcf4b36209739"},
|
{file = "black-18.9b0-py36-none-any.whl", hash = "sha256:817243426042db1d36617910df579a54f1afd659adb96fc5032fcf4b36209739"},
|
||||||
{file = "black-18.9b0.tar.gz", hash = "sha256:e030a9a28f542debc08acceb273f228ac422798e5215ba2a791a6ddeaaca22a5"},
|
{file = "black-18.9b0.tar.gz", hash = "sha256:e030a9a28f542debc08acceb273f228ac422798e5215ba2a791a6ddeaaca22a5"},
|
||||||
]
|
]
|
||||||
|
camera-stage-mapping = []
|
||||||
cbor2 = [
|
cbor2 = [
|
||||||
{file = "cbor2-5.1.0.tar.gz", hash = "sha256:43ce11e8c2fe4971d386d1a60cf83bfa0a4a667b97668ba76acbf5e6398821aa"},
|
{file = "cbor2-5.1.0.tar.gz", hash = "sha256:43ce11e8c2fe4971d386d1a60cf83bfa0a4a667b97668ba76acbf5e6398821aa"},
|
||||||
]
|
]
|
||||||
|
|
|
||||||
|
|
@ -48,6 +48,7 @@ opencv-python-headless = [
|
||||||
]
|
]
|
||||||
labthings = "0.6.3"
|
labthings = "0.6.3"
|
||||||
pynpm = "^0.1.2"
|
pynpm = "^0.1.2"
|
||||||
|
camera-stage-mapping = {git = "https://gitlab.com/openflexure/microscope-extensions/camera-stage-mapping.git"}
|
||||||
|
|
||||||
[tool.poetry.extras]
|
[tool.poetry.extras]
|
||||||
rpi = ["RPi.GPIO"]
|
rpi = ["RPi.GPIO"]
|
||||||
|
|
|
||||||
Loading…
Add table
Add a link
Reference in a new issue