Examples

2D Refocusing of an HL60 cell

The data show a live HL60 cell imaged with quadriwave lateral shearing interferometry (SID4Bio, Phasics S.A., France). The diameter of the cell is about 20µm.

refocus_cell.py

import matplotlib.pylab as plt
from skimage.restoration import unwrap_phase
import numpy as np

import nrefocus

from example_helper import load_cell

# load initial cell
cell1 = load_cell("HL60_field.zip")

# refocus to two different positions
cell2 = nrefocus.refocus(cell1, 15, 1, 1)  # forward
cell3 = nrefocus.refocus(cell1, -15, 1, 1)  # backward

# amplitude range
vmina = np.min(np.abs(cell1))
vmaxa = np.max(np.abs(cell1))
ampkw = {"cmap": plt.get_cmap("gray"),
         "vmin": vmina,
         "vmax": vmaxa}

# phase range
cell1p = unwrap_phase(np.angle(cell1))
cell2p = unwrap_phase(np.angle(cell2))
cell3p = unwrap_phase(np.angle(cell3))
vminp = np.min(cell1p)
vmaxp = np.max(cell1p)
phakw = {"cmap": plt.get_cmap("coolwarm"),
         "vmin": vminp,
         "vmax": vmaxp}

# plots
fig, axes = plt.subplots(2, 3, figsize=(8, 4.5))
axes = axes.flatten()
for ax in axes:
    ax.xaxis.set_major_locator(plt.NullLocator())
    ax.yaxis.set_major_locator(plt.NullLocator())

# titles
axes[0].set_title("focused backward")
axes[1].set_title("original image")
axes[2].set_title("focused forward")

# data
mapamp = axes[0].imshow(np.abs(cell3), **ampkw)
axes[1].imshow(np.abs(cell1), **ampkw)
axes[2].imshow(np.abs(cell2), **ampkw)
mappha = axes[3].imshow(cell3p, **phakw)
axes[4].imshow(cell1p, **phakw)
axes[5].imshow(cell2p, **phakw)

# colobars
cbkwargs = {"fraction": 0.045}
plt.colorbar(mapamp, ax=axes[0], label="amplitude [a.u.]", **cbkwargs)
plt.colorbar(mapamp, ax=axes[1], label="amplitude [a.u.]", **cbkwargs)
plt.colorbar(mapamp, ax=axes[2], label="amplitude [a.u.]", **cbkwargs)
plt.colorbar(mappha, ax=axes[3], label="phase [rad]", **cbkwargs)
plt.colorbar(mappha, ax=axes[4], label="phase [rad]", **cbkwargs)
plt.colorbar(mappha, ax=axes[5], label="phase [rad]", **cbkwargs)

plt.suptitle("Refocused cell on CPU with Numpy")
plt.tight_layout()
plt.show()

Compare available Metrics for Refocusing

The HL60 cell data used is described in refocus_cell.py.

compare_metrics.py

import nrefocus
from nrefocus.metrics import METRICS
import matplotlib.pyplot as plt

from examples.example_helper import load_cell

nrefocus.set_ndarray_backend("cupy")

dataset = load_cell("HL60_field.zip")

# load initial cell and create rf object
rf = nrefocus.iface.RefocusCupy(field=dataset,
                                wavelength=647e-9,
                                pixel_size=0.139e-6,
                                kernel="helmholtz")

# autofocus the image for each metric
my_metrics = list(METRICS.keys())

fig, axes = plt.subplots(1, len(my_metrics), figsize=(12, 5))

for i, mt in enumerate(my_metrics):
    af_vals = rf.autofocus(metric=mt,
                           minimizer="legacy",
                           interval=(-5e-6, 5e-6),
                           ret_field=True,
                           ret_grid=True)
    d, grid, field = af_vals
    d = d.get()
    field = field.get()
    grid_0 = grid[0].get()
    grid_1 = grid[1].get()

    axes[i].plot(grid_0, grid_1)
    axes[i].axvline(d, color='k', ls='--')
    axes[i].set_title(f"Metric {mt}")

fig.suptitle("Comparison of Metrics")
fig.tight_layout()
plt.show()

Refocus a stack of 2D fields (3D array input).

This example shows how to refocus a stack shaped (n, y, x). nrefocus interprets the last two axes as spatial axes and broadcasts over any leading “batch” axes.

refocus_stack.py

import numpy as np

import nrefocus

rng = np.random.default_rng(0)

# Synthetic complex field stack: (n, y, x)
n, ny, nx = 5, 128, 96
amp = 1 + 0.1 * rng.normal(size=(n, ny, nx))
pha = 0.3 * rng.normal(size=(n, ny, nx))
field = amp * np.exp(1j * pha)

pixel_size = 1e-6
dist = 2.13 * pixel_size

# CPU (NumPy) refocus of the entire stack in one call
rf = nrefocus.RefocusNumpy(
    field=field,
    wavelength=8.25 * pixel_size,
    pixel_size=pixel_size,
    medium_index=1.533,
    distance=0,
    kernel="helmholtz",
    padding=True,
)
refocused = rf.propagate(distance=dist)
print("nrefocus now accepts stacks of 2D arrays (3D input)!\n"
      "Input shape:", field.shape, "| Output shape-> ", refocused.shape)