Montaging¶
Instamatic is a tool for automated electron diffraction data collection. It has interfaces for interfacing with the TEM (JEOL/TFS) and several cameras (ASI Timepix/Cheetah/TVIPS/Merlin).
This notebook shows how to process a grid montage using instamatic
, pick grid squares, and set up an acquisition (acquire_at_items
). The data were collected on a zeolite sample (2020-02-12), using a JEOL JEM-1400 @ 120 kV in combination with a TVIPS F-416 camera.
The data for this demo are available from zenodo: https://doi.org/10.5281/zenodo.3924089
Make sure to change the work directory below to point at the right location.
from instamatic.montage import *
import numpy as np
from pathlib import Path
np.set_printoptions(suppress=True)
# work directory
work = Path(r"C:/s/data/montage_1")
Setting up the montage¶
Load the montage.yaml
file and the associated images.
m = InstamaticMontage.from_montage_yaml(work / 'montage.yaml')
m.gridspec
First, we can check what the data actually look like. To do so, we can simply stitch
and plot
the data using a binning=4
to conserve a bit of memory. This naively plots the data at the expected positions. Although the stitching is not that great, it's enough to get a feeling for the data.
m.calculate_montage_coords()
m.stitch(binning=4)
m.plot()
To get better stitching, we need to:
- Better estimate the difference vectors between each tile using cross correlation
- Optimize the coordinates of the difference vectors using least-squares minimization
This approach is based on Globally optimal stitching of tiled 3D microscopic image acquisitions by Preibish et al., Bioinformatics 25 (2009), 1463–1465 (https://doi.org/10.1093/bioinformatics/btp184).
Some metrics, such as the obtained shifts and FFT scores are plotted to evaluate the stitching.
# Use cross correlation to get difference vectors
m.calculate_difference_vectors(threshold='auto',
method='skimage',
plot=False)
# plot the fft_scores
m.plot_fft_scores()
# plot the pixel shifts
m.plot_shifts()
# get coords optimized using cross correlation
m.optimize_montage_coords(plot=True)
# stitch image, use binning 4 for speed-up and memory conservation
m.stitch(binning=4)
# plot the stitched image
m.plot()
We can save the stitched image:
m.export(work / "stitched.tiff")
When the image has been stitched (with or without optimization), we can look for the positions of the grid squares/squircles. To do so, call the method .find_holes
. The grid squares are identified as objects roughly sized diameter
with a tolerance of 10%. The median as well as 5/95 percentiles are printed to evaluate the hole size distribution.
diameter = 45_000 # nm
stagecoords, imagecoords = m.find_holes(plot=True, tolerance=0.1)
It is possible to optimize the stage coordinates for more efficient navigation. In this example, the total stage movement can be reduced by about 75%, which will save a lot of time. The function uses the two-opt algorithm for finding the shortest path: https://en.wikipedia.org/wiki/2-opt.
from pyserialem.navigation import sort_nav_items_by_shortest_path
stagecoords = sort_nav_items_by_shortest_path(stagecoords, plot=True);
The stagecoords
can be used to set up an automated acquire at items. This is useful to take medium magnification images from the regions of interest, in our case, the grid squares. First, initialize a connection to the microscope.
from instamatic import TEMController
ctrl = TEMController.initialize()
Next, we should set up an acquisition function for each stage position. This should:
- Center the grid square by aligning it with a reference image
- Take an image at high mag
- Store the image and the corresponding stage position in a buffer
In this preparation step, we much first obtain a reference image from a grid square. The magnification should be so that the grid square fits in the view of the image. In this example, we use 300x
in lowmag
.
import mrcfile
# set microscope conditions
ctrl.mode.set('lowmag')
ctrl.magnification.value = 300
binsize = 4
# reference image of a centered grid square
ref_img = ctrl.get_raw_image()
f = mrcfile.new(work / 'template.mrc', data=ref_img.astype(np.int16), overwrite=True)
f.close()
Then, we can define our acquisition function. We will align the stage to the reference image using cross correlation (ctrl.align_to
), and then taken an image of the centered grid square. Although this step is optional, it makes sure that the grid square is centered. This helps when looking for particles, and reduces errors related to stage translation and calibration. We also acquire an image and store the new stage position to a buffer.
buffer = []
stagepos = []
def acquire_func(ctrl):
# Align to template
# ctrl.align_to(ref_img, apply=True)
# obtain image
img, h = ctrl.get_image(binsize=binsize)
buffer.append(img)
# store stage position and image somewhere
pos = ctrl.stage.get()
stagepos.append(pos)
When the function is defined, we can pass it and the list of grid square stage coordinates to the function ctrl.acquire_at_items
, which will automate the function at each stage position.
sel = stagecoords[0:10] # Acquire at the first 10 items
ctrl.acquire_at_items(sel, acquire=acquire_func)
Here is a minimal example of how the acquire functions can be changed to collect data can be saved to a .nav
file which can be read by SerialEM
.
This makes use of the ability to pass a post_acquire
function to .acquire_at_items
. The post acquisition can be used to save the images as well as the required metadata to SerialEM
format, making use of the instamatic.serialem
module.
The post_acquire
function saves the data to .mrc
format, and writes an input file for SerialEM: instamatic.nav
.
from pyserialem import MapItem, write_nav_file
from instamatic import config
# reference image of a centered grid square
ref_img = ctrl.get_raw_image()
f = mrcfile.new(work / "template.mrc", data=ref_img.astype(np.int16), overwrite=True)
f.close()
# empty buffers
buffer = []
stagepos = []
def write_mrc_stack(fn:str, data: list, overwrite: bool=True, mmap:bool = True):
"""Write a stack of images to an mrc file."""
if mmap:
shape = (len(buffer), *buffer[0].shape)
with mrcfile.new_mmap(fn, shape=shape, overwrite=True, mrc_mode=1) as f:
for i, im in enumerate(buffer):
f.data[i] = im
else:
data = np.array(data)
# mrc can only be saved as a 16-bit integer
data = data.astype(np.int16)
try:
f = mrcfile.new(fn, data=data, overwrite=overwrite)
except OSError:
f.close()
def post_acquire(ctrl):
fn_nav = work / 'instamatic.nav'
fn_mrc = work / 'mmm.mrc'
write_mrc_stack(fn_mrc, buffer)
items = []
magnification = ctrl.magnification.value
mode = ctrl.mode.get()
mapscalemat = config.calibration[mode]['stagematrix'][magnification]
mapscalemat = [item/binsize for item in mapscalemat]
for i, image in enumerate(buffer):
x, y, z, _, _ = stagepos[i]
shape = image.shape
# binsize = ctrl.cam.get_binning()
d = {}
d['StageXYZ'] = x / 1000, y / 1000, z / 1000
d['MapFile'] = fn_mrc
d['MapSection'] = i
d['MapBinning'] = binsize
d['MapMagInd'] = ctrl.magnification.absolute_index + 1 # SerialEM is 1-based
d['MapScaleMat'] = mapscalemat
d['MapWidthHeight'] = shape
map_item = MapItem.from_dict(d)
items.append(map_item)
write_nav_file(fn_nav, *items)
print(f"Data saved to `{fn_nav}` and `{fn_mrc}` ({len(buffer)} images)")
Next we call ctrl.acquire_at_items
as before with the new post_acquire
function.
sel = stagecoords[0:10] # Acquire at the first 10 items
ctrl.acquire_at_items(sel,
acquire=acquire_func,
post_acquire=post_acquire)
Now load up the instamatic.nav
file in SerialEM
to see the result!