diff --git a/deapi/simulated_server/initialize_server.py b/deapi/simulated_server/initialize_server.py index ce533217..f9e3c858 100644 --- a/deapi/simulated_server/initialize_server.py +++ b/deapi/simulated_server/initialize_server.py @@ -28,13 +28,42 @@ def _recv_exact(conn, n): return buf +def serve_twin(port, twin_args=()): + """Serve frames rendered by the de-twin digital twin instead of the built-in fake data. + + Needs the optional ``twin`` extra (``pip install "deapi[twin]"``). ``twin_args`` are + passed to the twin's server, e.g. ``["--specimen", "Apoferritin in ice", "--camera", + "Celeritas", "--soap-port", "5002"]``; see ``python -m de_twin.faces.deapi_server --help``. + """ + try: + from de_twin.faces import deapi_server + except ImportError: + sys.stderr.write( + 'pydeserver --twin needs the digital twin: pip install "deapi[twin]"\n' + ) + return 2 + return deapi_server.main([str(port), *twin_args]) + + # Defining main function def main(port=13240): parser = argparse.ArgumentParser() parser.add_argument("--port", type=int, help="Port to listen on") - args, _ = parser.parse_known_args() + parser.add_argument( + "--twin", + action="store_true", + help='serve frames from the de-twin digital twin (pip install "deapi[twin]"); ' + "other options are passed to the twin, e.g. --specimen, --camera, --soap-port", + ) + args, rest = parser.parse_known_args() if args.port: port = args.port + if args.twin: + # the twin may call back into this loop (its --deapi-loop option): don't recurse + sys.argv = [a for a in sys.argv if a != "--twin"] + if rest and rest[0].isdigit(): # the positional port, already in ``port`` + rest = rest[1:] + return serve_twin(port, rest) HOST = "127.0.0.1" # Standard loopback interface address (localhost) PORT = port # Port to listen on (non-privileged ports are > 1023) diff --git a/deapi/tests/test_fake_server/test_twin_server.py b/deapi/tests/test_fake_server/test_twin_server.py new file mode 100644 index 00000000..eef9c9e3 --- /dev/null +++ b/deapi/tests/test_fake_server/test_twin_server.py @@ -0,0 +1,96 @@ +"""``pydeserver --twin``: the simulated server backed by the de-twin digital twin.""" + +import socket +import subprocess +import sys +import threading +import time + +import numpy as np +import pytest + +from deapi.client import Client +from deapi.simulated_server import initialize_server + + +def _free_port(): + with socket.socket() as s: + s.bind(("127.0.0.1", 0)) + return s.getsockname()[1] + + +def test_twin_flag_without_the_twin_installed_explains_how_to_get_it( + monkeypatch, capsys +): + monkeypatch.setitem(sys.modules, "de_twin", None) # import de_twin -> ImportError + monkeypatch.setitem(sys.modules, "de_twin.faces", None) + monkeypatch.setattr(sys, "argv", ["pydeserver", "--twin"]) + assert initialize_server.main(_free_port()) == 2 + assert 'pip install "deapi[twin]"' in capsys.readouterr().err + + +def test_twin_server_serves_twin_frames(): + pytest.importorskip("de_twin") + port = _free_port() + proc = subprocess.Popen( + [ + sys.executable, + "-u", + "-m", + "deapi.simulated_server.initialize_server", + str(port), + "--twin", + "--camera", + "DESim", + "--specimen", + "Dense Au on holey C", + ], + stdout=subprocess.PIPE, + stderr=subprocess.STDOUT, + text=True, + ) + lines = [] + started = threading.Event() + + def read(): + for line in proc.stdout: + lines.append(line) + if "started" in line: + started.set() + + threading.Thread(target=read, daemon=True).start() + try: + assert started.wait(120), "".join(lines) + client = Client() + client.usingMmf = False + client.connect(port=port) + assert client["Sensor Size X (pixels)"] == 1024 # the twin's DESim model + client["Frame Count"] = 2 + client.start_acquisition(1) + deadline = time.time() + 60 + while client.acquiring and time.time() < deadline: + time.sleep(0.05) + image, *_ = client.get_result("singleframe_integrated") + assert image.shape == (1024, 1024) + assert np.asarray(image).std() > 0 # a rendered specimen, not a constant + assert client["Instrument Project Magnification"] # the twin's column metadata + client.disconnect() + finally: + proc.terminate() + proc.wait(10) + + +def test_twin_options_are_passed_through(monkeypatch): + seen = {} + monkeypatch.setattr( + initialize_server, + "serve_twin", + lambda port, args: seen.update(port=port, args=args), + ) + monkeypatch.setattr( + sys, + "argv", + ["pydeserver", "13241", "--twin", "--soap-port", "5002", "--camera", "DE16"], + ) + initialize_server.main(13241) + assert seen == {"port": 13241, "args": ["--soap-port", "5002", "--camera", "DE16"]} diff --git a/doc/_static/digital_twin_holes.png b/doc/_static/digital_twin_holes.png new file mode 100644 index 00000000..cca569a5 Binary files /dev/null and b/doc/_static/digital_twin_holes.png differ diff --git a/doc/help/pyDEServer.rst b/doc/help/pyDEServer.rst index 3c49e8b6..57d0d0c8 100644 --- a/doc/help/pyDEServer.rst +++ b/doc/help/pyDEServer.rst @@ -62,4 +62,27 @@ The pyDEServer "Cheats" in a couple of ways: points of interest. This is implemented with the `BaseFakeData` class. Which implements a `__getitem__` method that returns the -data at the index in the navigation data. This can be used to return a single frame or a set of frames. \ No newline at end of file +data at the index in the navigation data. This can be used to return a single frame or a set of frames. +Realistic data from the digital twin +------------------------------------ + +For data that behaves like a real instrument, the pyDEServer can serve frames rendered by +`de-twin `_, a digital twin of a Direct Electron +camera on a TEM (column, specimen, in-situ holder and detector; Python 3.10+). Install the optional extra and +start the server with ``--twin``: + +.. code-block:: + + pip install "deapi[twin]" + pydeserver --port 13241 --twin --camera DE16 --specimen "Dense Au on holey C" + +Clients connect exactly as before (``client.usingMmf = False``). The twin supplies: + +* raw detector frames with dark offset, noise and gain structure, and references that + work like DE-Server's; +* images that follow the simulated column (stage, magnification, defocus); +* the ``Instrument ...`` metadata properties. + +Other options are passed to the twin. For example, ``--soap-port 5002`` also serves the +twin's microscope as a DE-TEM-Channel, and ``--seed``, ``--holder`` and ``--time-scale`` are +available too; see ``python -m de_twin.faces.deapi_server --help``. diff --git a/examples/digital_twin/README.rst b/examples/digital_twin/README.rst new file mode 100644 index 00000000..caaeb93d --- /dev/null +++ b/examples/digital_twin/README.rst @@ -0,0 +1,6 @@ +Digital twin +------------ + +Examples that run against the `de-twin `_ digital +twin (``pydeserver --twin``): a simulated camera, microscope and specimen that respond to each +other like a real instrument. diff --git a/examples/digital_twin/imaging_holes_with_the_digital_twin_sgskip.py b/examples/digital_twin/imaging_holes_with_the_digital_twin_sgskip.py new file mode 100644 index 00000000..fd9cfabd --- /dev/null +++ b/examples/digital_twin/imaging_holes_with_the_digital_twin_sgskip.py @@ -0,0 +1,147 @@ +""" +Imaging holes on a grid with the digital twin +============================================= + +A small automation loop that drives the microscope and the camera together, the way it +would run on an instrument: + +1. take a low-magnification atlas of a grid square, +2. find the holes in the holey carbon film, +3. move the stage to each hole with ``de_microscope`` (the DE-TEM-Channel client) and take + an image at higher magnification. + +Here the microscope and the camera are the `de-twin `_ +digital twin, so the example runs on any computer. Start it in a terminal first; it serves +a DE Server (deapi) on port 13250 and a DE-TEM-Channel on port 5002 for the same simulated +instrument, so stage moves change the images: + +.. code-block:: bash + + pip install "deapi[twin]" + pydeserver --port 13250 --twin --soap-port 5002 --camera DE16 --specimen "Dense Au on holey C" + +Against a real instrument only the addresses change: DE-Server's port and the +DE-TEM-Channel host. + +.. image:: /_static/digital_twin_holes.png + :alt: Atlas with the holes found, and an image of each hole +""" + +import time + +import matplotlib.pyplot as plt +import numpy as np +from scipy import ndimage + +import deapi +from de_microscope import Microscope + +client = deapi.Client() +client.usingMmf = False # the twin sends images over the socket +client.connect(port=13250) +scope = Microscope(host="127.0.0.1", port=5002) + + +def acquire(frames=10, fps=40): + """One integrated, dark/gain-corrected image (like DE-MC's 'Single').""" + client["Frames Per Second"] = fps + client["Frame Count"] = frames + client.start_acquisition(1) + while client.acquiring: + time.sleep(0.05) + image, *_ = client.get_result("singleframe_integrated") + return np.asarray(image, dtype=float) + + +def pixel_size_um(): + """Specimen pixel size reported by the server for the current magnification.""" + return client["Specimen Pixel Size X (nanometers)"] / 1000.0 + + +# %% +# A low-magnification atlas +# ------------------------- +# Spread the beam and go to low magnification so a few holes of the holey carbon fit in the +# field of view. Holes are where the beam passes through vacuum, so they are the brightest +# regions of the image. + +scope["Intensity"] = 0.8 +scope.set("Magnification", 2000, wait=True) +scope.set( + "StagePosition", {"x": 40.0, "y": 0.0}, wait=True +) # into a grid square, off the bar +start = scope["StagePosition"] +atlas = acquire() +atlas_px_um = pixel_size_um() + +# %% +# Find the holes +# -------------- +# Smooth, threshold the bright regions, and keep blobs of a plausible hole size +# (about 2 um across on this film) that lie wholly inside the atlas; a hole cut by the +# image edge would give a biased centre. Each hole's offset from the image centre, +# converted to micrometres, is how far the stage has to move to centre it. + +smooth = ndimage.gaussian_filter(atlas, 8) +bright = smooth > np.percentile(smooth, 80) +labels, n = ndimage.label(bright) +areas = ndimage.sum(np.ones_like(atlas), labels, range(1, n + 1)) * atlas_px_um**2 +centres = ndimage.center_of_mass(bright, labels, range(1, n + 1)) +ny, nx = atlas.shape +margin = 1.5 / atlas_px_um # pixels: a hole radius plus a little +holes = [ + ((c - nx / 2) * atlas_px_um, (r - ny / 2) * atlas_px_um, (r, c)) + for (r, c), a in zip(centres, areas) + if 1.0 < a < 8.0 # um^2: about one hole + and margin < r < ny - margin + and margin < c < nx - margin +] +print(f"found {len(holes)} holes") + +# %% +# Visit each hole +# --------------- +# Stage moves carry the specimen with them, so a feature at +dx in the image is centred by +# moving the stage by -dx. On a real column, calibrate the sign and rotation once +# (for example with ``de_microscope``'s stage calibration) and apply them here. + +scope.set("Magnification", 8000, wait=True) +scope["Intensity"] = 0.6 # converge the beam again for the close-ups +hole_images = [] +for dx_um, dy_um, _ in holes[:4]: + scope.set( + "StagePosition", {"x": start["x"] - dx_um, "y": start["y"] - dy_um}, wait=True + ) + hole_images.append(acquire(frames=40)) # 1 s +scope.set("StagePosition", {"x": start["x"], "y": start["y"]}, wait=True) + +# %% +# Plot the atlas and the holes +# ---------------------------- + + +def show(ax, img): + # the carbon film only dims the beam by ~25 %, while the gold particles are nearly + # black: set the grey range from the film and the holes, not the particles + lo, hi = np.percentile(img, [10, 99.8]) + ax.imshow(img[::4, ::4], cmap="gray", vmin=lo, vmax=hi) + + +fig, axes = plt.subplots( + 1, 1 + len(hole_images), figsize=(3.2 * (1 + len(hole_images)), 3.4) +) +axes = np.atleast_1d(axes) +show(axes[0], atlas) +for i, (_, _, (r, c)) in enumerate(holes[: len(hole_images)]): + axes[0].plot(c / 4, r / 4, "o", mfc="none", mec="C1", ms=14) + axes[0].text(c / 4, r / 4, str(i + 1), color="C1", ha="center", va="center") +axes[0].set_title("atlas, 2000x") +for i, (ax, img) in enumerate(zip(axes[1:], hole_images)): + show(ax, img) + ax.set_title(f"hole {i + 1}, 8000x") +for ax in axes: + ax.axis("off") +fig.tight_layout() +plt.show() + +client.disconnect() diff --git a/pyproject.toml b/pyproject.toml index 70a19de9..c85bb646 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -60,6 +60,10 @@ tests = [ service = [ "pyqt6", ] +# Serve frames from the de-twin digital twin: `pydeserver --twin` +twin = [ + "de-twin>=0.1.1; python_version >= '3.10'", +] doc = [ "sphinx", "pydata_sphinx_theme", diff --git a/upcoming_changes/59.doc.rst b/upcoming_changes/59.doc.rst new file mode 100644 index 00000000..0915ae29 --- /dev/null +++ b/upcoming_changes/59.doc.rst @@ -0,0 +1 @@ +Added an example that moves the stage with ``de_microscope`` and images holes in a grid, run against the digital twin. diff --git a/upcoming_changes/59.new_feature.rst b/upcoming_changes/59.new_feature.rst new file mode 100644 index 00000000..1f0e83e9 --- /dev/null +++ b/upcoming_changes/59.new_feature.rst @@ -0,0 +1 @@ +``pydeserver --twin`` serves frames rendered by the `de-twin `_ digital twin instead of the built-in fake data: a simulated microscope, specimen and detector, with a DE-TEM-Channel (``--soap-port``) so stage and optics changes show up in the images. Install it with ``pip install "deapi[twin]"``.