rahix/openflexure-microscope-server
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Merge branch 'scanning-stability' into 'v3'

Browse source 
Improve the stability of scanning

Closes #599 and #600

See merge request openflexure/openflexure-microscope-server!434
This commit is contained in:
Joe Knapper 2025-12-02 15:47:03 +00:00
commit 99afe89c10
10 changed files with 874 additions and 81 deletions
Download patch file Download diff file Expand all files Collapse all files

137
src/openflexure_microscope_server/background_detect.py
View file

@ -10,7 +10,6 @@ import cv2
import numpy as np
from pydantic import BaseModel, Field, ConfigDict
from pydantic.errors import PydanticUserError
from scipy.stats import norm
from labthings_fastapi.thing_description import type_to_dataschema
@ -18,6 +17,14 @@ class MissingBackgroundDataError(RuntimeError):
"""An error raised if checking for sample without background data set."""
class ChannelBlankError(RuntimeError):
"""An error raised if a channel has no measured standard deviation.
This is not physical and usually means the camera has not yet booted or changed
mode fully.
"""
class BackgroundDetectorStatus(BaseModel):
"""The status information about a background detector instance needed for the GUI.
@ -136,7 +143,7 @@ class BackgroundDetectAlgorithm:
"Each background detect algorithm must implement an image_is_sample method."
)
def set_background(self, image: np.ndarray) -> BaseModel:
def set_background(self, image: np.ndarray) -> None:
"""Use the input image to update the background data.
Background data must be a Pydantic BaseModel.
@ -188,6 +195,10 @@ class ColourChannelDetectLUV(BackgroundDetectAlgorithm):
background_data_model: BaseModel = ChannelDistributions
settings_data_model: BaseModel = ColourChannelDetectSettings
# These are the same as those used for ChannelDeviationLUV. More detail is
# provided there.
min_stds = [0.5, 0.3, 0.5]
def background_mask(self, image: np.ndarray) -> np.ndarray:
"""Calculate a binary image, showing whether each pixel is background.
@ -196,11 +207,6 @@ class ColourChannelDetectLUV(BackgroundDetectAlgorithm):
The image should be in LUV format, the output will be binary with the
same shape in the first two dimensions.
"""
if not self.background_data:
raise MissingBackgroundDataError(
"Background is not set: you need to calibrate background detection."
)
# The ``[1:]`` selects only the U and V channels of the image.
# Only U and V are used as brightness (L) often changes as
# the height of the sample changes.
@ -225,6 +231,10 @@ class ColourChannelDetectLUV(BackgroundDetectAlgorithm):
:returns: A value (between 0 and 100) is the percentage of the image that is
sample.
"""
if not self.background_data:
raise MissingBackgroundDataError(
"Background is not set: you need to calibrate background detection."
)
image_luv = cv2.cvtColor(image, cv2.COLOR_RGB2LUV)
mask = self.background_mask(image_luv)
@ -250,15 +260,114 @@ class ColourChannelDetectLUV(BackgroundDetectAlgorithm):
"""Use the input image to update the background distributions."""
image_luv = cv2.cvtColor(image, cv2.COLOR_RGB2LUV)
ch1 = (image_luv.T[0]).flatten()
ch2 = (image_luv.T[1]).flatten()
ch3 = (image_luv.T[2]).flatten()
mu = np.mean(image_luv, axis=(0, 1))
std = np.std(image_luv, axis=(0, 1))
points = np.array([np.asarray(ch1), np.asarray(ch2), np.asarray(ch3)]).T
# Get the mean and standard deviation of values in each channel
mu, std = np.apply_along_axis(norm.fit, 0, points)
if np.any(std == 0):
raise ChannelBlankError("Some LUV channels have no standard deviation.")
std = np.maximum(std, self.min_stds)
self.background_data = ChannelDistributions(
means=mu.tolist(), standard_deviations=std.tolist()
)
class ChannelDeviationLUV(BackgroundDetectAlgorithm):
"""Compare the standard deviations of the LUV channels in a grid to background data.
Using an LUV colour space, each image is divided into an 8x8 grid of images.
The standard deviation of each channel of each image is calculated and compared
to the median standard deviation for a grid of background images.
"""
# Note we don't use the means in this algorithm but we use the same channel
# distributions model
background_data_model: BaseModel = ChannelDistributions
settings_data_model: BaseModel = ColourChannelDetectSettings
# Empirically, 0.5 seems to be approximate the standard deviation for a good image
# in L and V. U appears to be about 60% of this value. U is about 65% of V when
# converting white-noise in RGB into LUV (note that we are using cv2's internal
# LUV colour space not converting to the CIELUV numbers)
min_stds = [0.5, 0.3, 0.5]
def get_sample_coverage(self, image: np.ndarray) -> float:
"""Return the percentage of the input image that is background.
Evaluate whether it is foreground or background by comparing the standard
deviations of an 8x8 grid of sub-images to the median standard deviation
from a background image.
:returns: A value (between 0 and 100) that is the percentage of the image that is
sample.
"""
if not self.background_data:
raise MissingBackgroundDataError(
"Background is not set: you need to calibrate background detection."
)
image_luv = cv2.cvtColor(image, cv2.COLOR_RGB2LUV)
stds = _chunked_stds(image_luv, 8, 8)
bg_stds = self.background_data.standard_deviations
l_cut = bg_stds[0] * self.settings.channel_tolerance
u_cut = bg_stds[1] * self.settings.channel_tolerance
v_cut = bg_stds[2] * self.settings.channel_tolerance
populated_regions = (
(stds[:, :, 0] > l_cut) | (stds[:, :, 1] > u_cut) | (stds[:, :, 2] > v_cut)
)
return float(100 * np.sum(populated_regions) / 64)
def image_is_sample(self, image: np.ndarray) -> tuple[bool, str]:
"""Label the current image as either background or sample.
:returns: A tuple of the result (boolean), and explanation string. The
explanation string is formatted so it can be added into a sentence such as
``An action was taken because the image is {message}.``
"""
sample_coverage = self.get_sample_coverage(image)
# Use bool otherwise get numpy variants of True and False.
is_sample = bool(sample_coverage > self.settings.min_sample_coverage)
message = f"{sample_coverage:0.1f}% sample"
if not is_sample:
message = "only " + message
return is_sample, message
def set_background(self, image: np.ndarray) -> None:
"""Use the input image to update the background distributions."""
image_luv = cv2.cvtColor(image, cv2.COLOR_RGB2LUV)
mu = np.zeros(3)
c_stds = _chunked_stds(image_luv, 8, 8)
channel_blank = np.all(c_stds == 0, axis=(0, 1))
if np.any(channel_blank):
raise ChannelBlankError("Some LUV channels have no standard devaition.")
std = np.median(c_stds, axis=(0, 1))
std = np.maximum(std, self.min_stds)
self.background_data = ChannelDistributions(
means=mu.tolist(), standard_deviations=std.tolist()
)
def _chunked_stds(img: np.ndarray, n_rows: int = 8, n_cols: int = 8) -> np.ndarray:
"""Split image into a grid and calculate std of each channel in each chunk.
:param img: The image to analyse
:param n_rows: The number of rows in the grid
:param n_cols: The number of cols in the grid
:return: A numpy array of the grid of standard deviations.
"""
h, w = img.shape[:2]
row_height = h // n_rows
col_width = w // n_cols
out = np.zeros((n_rows, n_cols, 3))
for i in range(n_rows):
for j in range(n_cols):
chunk = img[
i * row_height : (i + 1) * row_height,
j * col_width : (j + 1) * col_width,
]
out[i, j, :] = np.std(chunk, axis=(0, 1))
return out

195
src/openflexure_microscope_server/things/autofocus.py
View file

@ -28,6 +28,7 @@ from .stage import StageDependency as Stage
LOGGER = logging.getLogger(__name__)
MIN_TEST_IMAGE_COUNT = 3
MAX_TEST_IMAGE_COUNT = 9
EXTRA_STACK_CAPTURES = 15
class NotStreamingError(RuntimeError):
@ -133,7 +134,7 @@ class StackParams(BaseModel):
This is 15 images more then the minimum number that are captured.
"""
return self.min_images_to_test + 15
return self.min_images_to_test + EXTRA_STACK_CAPTURES
def slice_to_save(self, sharpest_index: int) -> slice:
"""Return the slice of images to save given the index of the sharpest image."""
@ -421,38 +422,38 @@ class AutofocusThing(lt.Thing):
is close to focus, but not quite within ``dz/2``. It will attempt to autofocus
up to 10 times.
"""
repeat = True
attempts = 0
attempt = 0
backlash = 200
with sharpness_monitor.run():
while repeat and attempts < 10:
while attempt < 10:
attempt += 1
if start == "centre":
stage.move_relative(x=0, y=0, z=-(backlash + dz / 2))
stage.move_relative(x=0, y=0, z=backlash)
# Always start centrally for future runs
start = "centre"
focus_data_index, _ = sharpness_monitor.focus_rel(
dz, block_cancellation=True
)
_times, heights, sizes = sharpness_monitor.move_data(focus_data_index)
peak_height = heights[np.argmax(sizes)]
height_min = np.min(heights)
height_max = np.max(heights)
target_min = np.min(heights) + dz / 5
target_max = np.max(heights) - dz / 5
if (
peak_height - height_min < dz / 5
or height_max - peak_height < dz / 5
):
attempts += 1
start = "centre"
stage.move_absolute(z=peak_height - backlash)
stage.move_absolute(z=peak_height)
else:
repeat = False
stage.move_relative(x=0, y=0, z=-(dz + backlash))
stage.move_absolute(z=peak_height)
return heights.tolist(), sizes.tolist()
# move to the peak
stage.move_absolute(z=peak_height - backlash)
stage.move_absolute(z=peak_height)
if target_min < peak_height < target_max:
# If it is within the target range then return
return heights.tolist(), sizes.tolist()
raise NoFocusFoundError(
"Looping autofocus couldn't converge on a focus location."
)
stack_images_to_save = lt.ThingSetting(
initial_value=1,
@ -566,7 +567,9 @@ class AutofocusThing(lt.Thing):
stage: Stage,
sharpness_monitor: SharpnessMonitorDep,
stack_parameters: StackParams,
) -> tuple[bool, int]:
save_on_failure: bool = False,
check_turning_points: bool = True,
) -> tuple[bool, Optional[int]]:
"""Run a smart stack.
A smart stack captures images offset in z, testing whether the sharpest image
@ -582,17 +585,22 @@ class AutofocusThing(lt.Thing):
supplied by LabThings dependency injection
:param stack_parameters: A StackParams object containing the required
parameters to run a stack.
:param save_on_failure: Whether to save an image even if no focus was found.
:param check_turning_points: Whether to check the number of turning points in
the sharpnesses of the images in the stack is exactly 1. (May fail with
thick samples)
:returns: A tuple containing:
* A boolean, True if stack was successfully
* The z position of the sharpest image
"""
tries = 0
# Loop until a stack is successful
while tries < stack_parameters.max_attempts:
attempt = 0
while True:
attempt += 1
success, captures, sharpest_id = self.z_stack(
stack_parameters=stack_parameters,
check_turning_points=check_turning_points,
cam=cam,
stage=stage,
)
@ -600,28 +608,33 @@ class AutofocusThing(lt.Thing):
if success:
break
if attempt >= stack_parameters.max_attempts:
break
# The z position of the first images in the previous attempt.
initial_z_pos = captures[0].position["z"]
# If a stack is not successful, move to the start and autofocus
self.reset_stack(
initial_z_pos,
stack_parameters.autofocus_dz,
stage,
sharpness_monitor,
)
try:
self.reset_stack(
initial_z_pos,
stack_parameters.autofocus_dz,
stage,
sharpness_monitor,
)
except NoFocusFoundError:
break
# Save stack_parameters.image_to_save images centred on the sharpest capture.
# If the smart_stack failed the exact number of images saved may not be
# stack_parameters.image_to_save
self.save_stack(
sharpest_id=sharpest_id,
captures=captures,
stack_parameters=stack_parameters,
cam=cam,
)
if success or save_on_failure:
self.save_stack(
sharpest_id=sharpest_id,
captures=captures,
stack_parameters=stack_parameters,
cam=cam,
)
# Return whether or not the smart stack was successful, and the z position of
# the sharpest image, for path planning and tracking
return success, _get_capture_by_id(captures, sharpest_id).position["z"]
def reset_stack(
@ -681,9 +694,10 @@ class AutofocusThing(lt.Thing):
def z_stack(
self,
stack_parameters: StackParams,
check_turning_points: bool,
cam: CameraClient,
stage: Stage,
) -> tuple[bool, list[CaptureInfo], Optional[int]]:
) -> tuple[bool, list[CaptureInfo], int]:
"""Capture a series of images checking that sharpest image central.
The images are separated in z offset by stack_parameters.stack_dz, as they
@ -692,6 +706,9 @@ class AutofocusThing(lt.Thing):
completes.
:param stack_parameters: a StackParams object holding stack parameters
:param check_turning_points: Whether to check the number of turning points in
the sharpnesses of the images in the stack is exactly 1. (May fail with
thick samples)
:param cam: Camera Dependency to be passed through from the calling action
:param stage: Stage Dependency to be passed through from the calling action
@ -699,7 +716,7 @@ class AutofocusThing(lt.Thing):
* the stack result (True for successful stack, False for failed stack),
* a list of CaptureInfo objects,
* the buffer_id of the sharpest image (or None if the stack failed).
* the buffer_id of the sharpest image
"""
# Move down by the height of the z stack, plus an overshoot
# Better to start too low and take too many images than too high and need to refocus
@ -718,7 +735,7 @@ class AutofocusThing(lt.Thing):
ims_to_check = slice(-stack_parameters.min_images_to_test, None)
# If the sharpest image isn't found within the maximum number of images
# end the loop and return "restart"
# end the loop and return False indicating the stack failed.
while len(captures) < stack_parameters.max_images_to_test:
time.sleep(stack_parameters.settling_time)
@ -733,16 +750,18 @@ class AutofocusThing(lt.Thing):
# If the number of images is enough to test, test them
if len(captures) >= stack_parameters.min_images_to_test:
result, capture_id = self.check_stack_result(captures[ims_to_check])
result, capture_id = self.check_stack_result(
captures[ims_to_check], check_turning_points=check_turning_points
)
if result == "success":
return True, captures, capture_id
if result == "restart":
return False, captures, None
return False, captures, capture_id
# If reached here the result was "continue"
stage.move_relative(z=stack_parameters.stack_dz)
return False, captures, None
return False, captures, capture_id
def capture_stack_image(
self,
@ -774,12 +793,15 @@ class AutofocusThing(lt.Thing):
# unlikely to improve readability. This function is basically a complex switch
# statement, having an explicit return after each option is clear.
def check_stack_result( # noqa: PLR0911
self, captures: list[CaptureInfo]
self, captures: list[CaptureInfo], check_turning_points: bool
) -> tuple[Literal["success", "continue", "restart"], int]:
"""Check if the sharpest image in a list of captures is central enough.
:param captures: a list of the capture objects to for testing if the
sharpness has converged in the centre
:param check_turning_points: Whether to check the number of turning points in
the sharpnesses of the images in the stack is exactly 1. (May fail with
thick samples)
:returns: A tuple with two values:
@ -793,28 +815,91 @@ class AutofocusThing(lt.Thing):
* capture_id - the buffer id of the sharpest image
"""
sharpest_index = np.argmax([capture.sharpness for capture in captures])
sharpnesses = np.array([capture.sharpness for capture in captures])
sharpest_index = np.argmax(sharpnesses)
# The buffer id of the sharpest image
capture_id = captures[sharpest_index].buffer_id
sharpness_length = len(captures)
n_imgs = len(captures)
# If only testing one image, then by definition the sharpest is central
if sharpness_length == 1:
if n_imgs == 1:
return "success", capture_id
# If testing three images, test if the centre is the sharpest
if sharpness_length == 3:
if n_imgs == 3:
if sharpest_index == 1:
return "success", capture_id
if sharpest_index == 0:
return "restart", capture_id
return "continue", capture_id
# For larger stacks, test if the best image is not within two of the edge of the stack
# ie for a stack of 7 images, best image must be between 3rd and and 5th
exclusion_range = 2
if sharpest_index < exclusion_range:
return "restart", capture_id
if sharpest_index >= sharpness_length - exclusion_range:
try:
turning = _get_peak_turning_point(sharpnesses)
except NotAPeakError:
return "continue", capture_id
# For larger stacks, test if the best image is not within two of the edge of
# the stack ie for a stack of 7 images, best image must be between 3rd and 5th
edge_size = 2
# Check both the peak from fitting and the sharpest image are not at the
# edge of the stack.
if turning < edge_size - 0.5 or sharpest_index < edge_size:
return "restart", capture_id
if turning > n_imgs - edge_size - 0.5 or sharpest_index >= n_imgs - edge_size:
return "continue", capture_id
if check_turning_points:
turning_points = _count_turning_points(sharpnesses)
if turning_points > 1:
return "continue", capture_id
return "success", capture_id
class NotAPeakError(RuntimeError):
"""The data to fit isn't a peak."""
class NoFocusFoundError(RuntimeError):
"""No focus found during looping Autofocus."""
def _get_peak_turning_point(sharpnesses: np.ndarray) -> float:
"""Get the turning point for a sharpnesses in a z-stack.
:param sharpnesses: A numpy array of sharpnesses
:return: The x value of the turning point where x-axis is 0 to N-1 for the N
sharpness values
:raise NotAPeakError: If the fit doesn't have a maximum within 95% confidence.
"""
# Fit the peak
x = range(len(sharpnesses))
coeffs, cov = np.polyfit(x, sharpnesses, deg=2, cov=True)
# a in the equation a*x^2 + bx + c
a = float(coeffs[0])
fit_func = np.poly1d(coeffs)
# find the peaks's x-position
turning = float(fit_func.deriv().roots[0])
# sigma_a the standard error of a
sigma_a = float(np.sqrt(cov[0, 0]))
# Estimate the upper 95% confidence bound for the x^2 term
ci_high = a + 2 * sigma_a
# If the high 95% confindence interval of the peak is not negative then the
# fit isn't sure if this is a peak or a U shape, so we continue.
# -1e-9 is used to guard against weird effects for flat and straight data
if ci_high > -1e-9:
raise NotAPeakError("Not a peak to within 95% confidence.")
return turning
def _count_turning_points(sharpnesses: np.ndarray) -> int:
"""Count the number of turing points, after rejecting those from noise."""
# Also take the difference of neighbouring points to check for turning points
d_sharpnesses = sharpnesses[1:] - sharpnesses[:-1]
# Filter out any points that are not prominent
prominent = abs(d_sharpnesses) / np.mean(abs(d_sharpnesses)) > 0.5
d_sharpnesses = d_sharpnesses[prominent]
# count the sign changes
return int(np.sum(d_sharpnesses[1:] * d_sharpnesses[:-1] < 0))

10
src/openflexure_microscope_server/things/camera/__init__.py
View file

@ -28,6 +28,7 @@ from labthings_fastapi.types.numpy import NDArray
from openflexure_microscope_server.ui import ActionButton, PropertyControl
from openflexure_microscope_server.background_detect import (
ColourChannelDetectLUV,
ChannelDeviationLUV,
BackgroundDetectAlgorithm,
BackgroundDetectorStatus,
)
@ -182,8 +183,11 @@ class BaseCamera(lt.Thing):
dictionary in this function. Configuration will be added at a later date.
"""
super().__init__()
self.background_detectors = {"Colour Channels (LUV)": ColourChannelDetectLUV()}
self._detector_name = "Colour Channels (LUV)"
self.background_detectors = {
"Colour Channels (LUV)": ColourChannelDetectLUV(),
"Channel Deviations (LUV)": ChannelDeviationLUV(),
}
self._detector_name = "Channel Deviations (LUV)"
def __enter__(self) -> None:
"""Open hardware connection when the Thing context manager is opened."""
@ -637,7 +641,7 @@ class BaseCamera(lt.Thing):
def detector_name(self, name: str) -> None:
"""Validate and set detector_name."""
if name not in self.background_detectors:
raise ValueError(f"{name} is not a valid background detector name")
LOGGER.warning(f"{name} is not a valid background detector name.")
self._detector_name = name
@property

14
src/openflexure_microscope_server/things/camera/picamera.py
View file

@ -42,6 +42,7 @@ from openflexure_microscope_server.ui import (
action_button_for,
property_control_for,
)
from openflexure_microscope_server.background_detect import ChannelBlankError
from . import picamera_recalibrate_utils as recalibrate_utils
from . import picamera_tuning_file_utils as tf_utils
@ -767,7 +768,6 @@ class StreamingPiCamera2(BaseCamera):
* ``auto_expose_from_minimum``
* ``set_static_green_equalisation`` to set geq offset to max
* ``calibrate_lens_shading`` (also sets colour gains for white balance)
* ``set_background``
"""
self.flat_lens_shading()
@ -775,8 +775,16 @@ class StreamingPiCamera2(BaseCamera):
self.set_static_green_equalisation()
self.set_ce_enable_to_off()
self.calibrate_lens_shading()
time.sleep(0.5)
self.set_background(portal)
for _i in range(3):
try:
time.sleep(self._sensor_info.long_pause)
self.set_background(portal)
# Return if background is set
return
except ChannelBlankError:
# If channel is blank, sleep a second and try again.
pass
raise RuntimeError("Couldn't set background")
@lt.thing_property
def primary_calibration_actions(self) -> list[ActionButton]:

8
src/openflexure_microscope_server/things/smart_scan.py
View file

@ -467,13 +467,17 @@ class SmartScanThing(lt.Thing):
continue
focused, focused_height = self._autofocus.run_smart_stack(
stack_parameters=self._stack_params
stack_parameters=self._stack_params,
save_on_failure=not self._scan_data.skip_background,
)
current_pos_xyz = (new_pos_xyz[0], new_pos_xyz[1], focused_height)
# An image was captured if we are focussed or we are not skipping background.
imaged = focused or not self._scan_data.skip_background
route_planner.mark_location_visited(
current_pos_xyz, imaged=True, focused=focused
current_pos_xyz, imaged=imaged, focused=focused
)
# increment capture counter as thread has completed