stitch — Image stitching core¶

The +utils/+stitch package holds the controller-independent, headless-testable core of the controllers.Stitching tool: layout builders, pairwise registration (phase correlation and feature-based), the global least-squares solver, canvas planning, and the three fusers (in-memory, streaming to OME-Zarr, and straight out to standard image files).

utils.stitch.autocropCanvas(layout, canvas)¶

AUTOCROPCANVAS - Shrink a planned canvas to the region every slice fully covers.

Syntax:

canvas = utils.stitch.autocropCanvas(layout, canvas)

Solved tile positions are never a perfect rectangle: the outer tiles end up a few pixels apart, so the mosaic carries a ragged background frame along all four sides. This trims it away by replacing the canvas with the LARGEST axis-aligned rectangle that is covered by tiles on EVERY output slice, and re-expressing the placement plan in the cropped frame. Both fusers then produce the cropped mosaic directly - nothing is fused and thrown away, so the streaming path benefits identically to the in-memory one.

The crop is in-plane only. Z is left alone: an output slice is either produced or it is not, and dropping end slices would silently change the depth of a stack the user asked for. Slices no tile reaches at all are SKIPPED rather than intersected in - an empty slice would otherwise veto every crop.

How the rectangle is found. Every tile footprint is a rectangle in canvas coordinates, so coverage is piecewise-constant on the grid formed by the footprint edges: at most 2N+2 distinct rows and columns for N tiles, whatever the mosaic’s pixel size. Coverage is evaluated on that compressed grid (once per distinct Z-segment - the contributing tile set and the zShifts correction only change at tile-band boundaries), intersected across segments, and the maximum-area all-covered rectangle is read off with the largest-rectangle-in-histogram stack algorithm weighted by the cell sizes. The result is exact in pixels; no coverage mask the size of the mosaic is ever allocated.

Warped (affine) plans use a CONSERVATIVE footprint: the affine image of a tile is a parallelogram, and the axis-aligned rectangle spanned by its middle two corner x-coordinates and middle two y-coordinates is inscribed in it. One further pixel is trimmed off each side because imwarp blends the outermost resampled row against the zero fill. Under-claiming coverage can only leave a slightly smaller mosaic; over-claiming would put back the black edge this exists to remove. Tiles whose transform is an integer translation take the exact rectangle, so a translation plan crops identically with or without canvas.tforms.

Input Arguments:

• layout - [struct array] tile layout (.tileSize).

• canvas - [struct] from utils.stitch.planCanvas().

Output Arguments:

• canvas - [struct] same fields, with .size(1:2), .tilePlacement, .boundingBox and (when present) .tforms / .tileBounds moved into the cropped frame, plus:

  • .cropRect - [1x4] [y0 y1 x0 x1] of the kept region in the ORIGINAL canvas frame. Absent when no fully covered region exists (the canvas is then returned untouched and a warning is issued).

Example - plan and crop:

canvas = utils.stitch.planCanvas(layout, positions);
canvas = utils.stitch.autocropCanvas(layout, canvas);
fprintf('kept %d x %d\n', canvas.size(1), canvas.size(2));

utils.stitch.blendWeights(tileHW, marginPx)¶

BLENDWEIGHTS - Linear feather (distance-ramp) blend weights for one tile.

Syntax:

weightMap = utils.stitch.blendWeights(tileHW)
weightMap = utils.stitch.blendWeights(tileHW, marginPx)

Builds a single-precision [H W] weight map that is 1 across the tile interior and ramps linearly down toward the edges, reaching a small positive value at the outermost pixel. The weight at any pixel is the minimum of its four separable edge ramps (distance to the nearest edge, normalised by marginPx), so overlapping tiles cross-fade smoothly and no pixel receives exactly zero weight (which keeps sum(w*I)/sum(w) well defined everywhere). The construction is separable (outer product of 1-D ramps) - no bwdist.

Input Arguments:

• tileHW - [1x2 double] tile size [H W].

• marginPx (optional) - [double] feather width in pixels from each edge (default: round(min(H, W) / 8), at least 1).

Output Arguments:

• weightMap - [H x W single] blend weights in (0, 1].

Example - feather weights with a 32-pixel ramp:

w = utils.stitch.blendWeights([512 512], 32);
imagesc(w); axis image; colorbar;

utils.stitch.buildLayoutAtlas(veMifPath, options)¶

BUILDLAYOUTATLAS - Build a tile layout (and optionally the stitch) from a Fibics Atlas mosaic.

Syntax:

layout = utils.stitch.buildLayoutAtlas(veMifPath)
[layout, edges, positions, atlasInfo] = utils.stitch.buildLayoutAtlas(veMifPath, options)

Reads the XML files a Fibics Atlas acquisition writes next to its tiles. Atlas records three successive stages of the same stitch in three files sharing one base name, and this function can take any prefix of that chain:

.ve-mif (always read) - the acquisition record: per-tile row/col

and NOMINAL stage position, tile size, FOV and pixel size. Becomes the layout’s gridRC / nomOrigin.

.ve-tie (options.importTies) - Atlas’s pairwise seam measurements.

Becomes the edges array, so MIB can solve without re-registering a single pixel.

.ve-updates (options.importPositions) - Atlas’s FINAL solved tile

positions. Becomes positions, so the mosaic is ready to fuse with nothing recomputed.

Why the nominal placement is only a starting guess. Under some imaging conditions the stage positions Atlas records do not describe where the tiles actually overlap (the sample data this was written against is off by ~32 px in Y - Atlas’s own ties agree). The nominal grid is therefore treated exactly like any other layout source: a rough placement that Measure overlaps refines. Import the ties (or the ties + positions) to keep Atlas’s own answer instead.

Coordinate conversion. Atlas works in micrometres on a stage frame whose Y axis points UP and whose X axis may run either way. Rather than hard-coding a vendor convention, the axis directions are DERIVED per mosaic by correlating each tile’s row/col attribute with its stage coordinate, so a mosaic acquired with a mirrored stage maps correctly without a flag. The same signs are then applied to the tie shifts and the solved positions, which share the stage frame’s orientation.

Tile files are resolved locally. The paths inside the XML are absolute paths on the acquisition machine (E:\...), which almost never exist where the data is analysed. Every tile is therefore looked up by its FILE NAME in the folder holding the .ve-mif, falling back to the recorded path only when the local file is missing.

Input Arguments:

• veMifPath - [char] full path to the MosaicInfo_*.ve-mif file.

• options (optional) - struct with fields:

  • .importTies - [logical] read the .ve-tie seam measurements into edges (default: false)

  • .importPositions - [logical] read the .ve-updates solved tile positions into positions (default: false). With no ties imported the matching edges are SYNTHESISED from the solved positions, so the state is self-consistent (a re-solve reproduces the same placement) and the seam inspector has something to review.

  • .tiePath / .updatesPath - [char] explicit sidecar paths; by default they are found next to the .ve-mif by utils.stitch.findAtlasSidecars().

  • .layout - [struct array] a layout already built from this very .ve-mif; supplied to skip re-reading every tile’s dimensions when only the sidecars are wanted (default: [] - build it).

Output Arguments:

• layout - struct array per the layout contract (see utils.stitch.buildLayoutGrid()), one entry per tile, ordered by row then col. Single Z layer - one .ve-mif is one section.

• edges - struct array of imported/synthesised seams (empty when neither sidecar was imported). Ties Atlas marked <User>true</User> (placed by hand in Atlas) carry .source = 'user', so MIB weights them like its own manual fixes and a re-measure preserves them; the rest are 'auto'.

• positions - [N x 3 double] solved origins [y x z], or [].

• atlasInfo - struct describing the mosaic: .mosaicName, .folder, .pixelSizeUm, .signX / .signY (derived axis directions), .fovUm, .overlapXpercent / .overlapYpercent, .numTilesX / .numTilesY and the per-tile .tiles records.

Example - fuse an Atlas mosaic exactly as Atlas stitched it:

opts = struct('importTies', true, 'importPositions', true);
[layout, edges, positions] = utils.stitch.buildLayoutAtlas( ...
    'D:\S_001\MosaicInfo_S_001.ve-mif', opts);
canvas = utils.stitch.planCanvas(layout, positions);
mosaic = utils.stitch.fuseInMemory(layout, canvas);

See also utils.stitch.findAtlasSidecars, utils.stitch.buildLayoutGrid

utils.stitch.buildLayoutBioFormats(inputPath, options)¶

BUILDLAYOUTBIOFORMATS - Build a tile layout from embedded Bio-Formats stage coordinates.

Syntax:

layout = utils.stitch.buildLayoutBioFormats(inputPath)
layout = utils.stitch.buildLayoutBioFormats(inputPath, options)

Reads the physical stage position stored in the OME metadata of each tile (Plane PositionX/Y/Z) and converts it into the stitcher’s pixel/slice nomOrigin frame via utils.stitch.stageCoordsToOrigins(). This is the 'Bio-Formats metadata' layout source: unlike the grid / position-file / filename-pattern sources it needs no user-supplied arrangement - the microscope already recorded where every tile sits.

inputPath may be:

• a single Bio-Formats file whose series are the tiles (the typical mosaic case, e.g. one .czi / .nd2 / .lif with N series), or

• a newline-separated list of files, or a folder - one file per tile, each carrying its own stage position.

Pixel size is assumed uniform across tiles (MIB convention) and taken from the first tile; the resulting dataset pixSize should likewise come from the first tile. Distinct stage-Z values become separate zLayer layers, so a multi-focus acquisition is jointly solved just like a position file with a Z column.

Input Arguments:

• inputPath - [char] file, newline-list of files, or folder (see above).

• options (optional) - struct with fields:

  • .flipX / .flipY - [logical] negate the stage axis when it runs opposite to the pixel axis (vendor-dependent; default false).

  • .bioFormatsMemoizerMemoDir - [char] memo dir (default: tempdir).

Output Arguments:

• layout - struct array per the layout contract (see utils.stitch.buildLayoutGrid()) with an extra .seriesIndex field (1-based Bio-Formats series) honoured by utils.stitch.makeTileReader().

Example - stitch a multi-series confocal mosaic:

layout = utils.stitch.buildLayoutBioFormats('C:\data\mosaic.czi');

utils.stitch.buildLayoutFilenamePattern(filenames, options)¶

BUILDLAYOUTFILENAMEPATTERN - Build a tile layout from MIB2 chop filename tokens.

Syntax:

layout = utils.stitch.buildLayoutFilenamePattern(filenames)
layout = utils.stitch.buildLayoutFilenamePattern(filenames, options)

Parses the Z## / X## / Y## tokens embedded in each tile filename (the format produced by MIB2’s rechop tool), e.g. myStack_Z01-X02-Y03.tif → Z=1, X=2, Y=3.

The three tokens are located INDEPENDENTLY (last occurrence of each letter in the base name), so their ORDER and the separators between them do not matter - _Z01-X02-Y03, _X02-Y03-Z01 and Y03X02Z01 all parse identically. Constraints that DO matter:

• the letter must be upper case and followed by EXACTLY two digits: exactly two characters are read, so Z1 errors and Z001 silently parses as 00 (a rename to two-digit tokens is required above 99 tiles per axis);

• indices are 1-based (Z01-X01-Y01 is the first tile);

• because the LAST occurrence wins, Z/X/Y may appear before the tokens (XYZstack_Z01-X01-Y01) but not after (..._Y01_XY).

Each entry may name a single image file OR a FOLDER holding the tile’s Z-stack (auto-detected per entry by utils.stitch.resolveTileEntry()) - for folder tiles the tokens live in the FOLDER name and follow the same rules.

Nominal origins are computed from the grid indices and the (uniform) tile size, with an optional XY overlap (default 0% = abutting chunks, the MIB2-rechop reassembly case). With a non-zero overlap the step shrinks like the Grid source, so pattern-named tiles from an overlapping acquisition get honest nominal positions that the measurement pass can then refine. Z layers always abut (no Z-overlap control). Origins are 1-based pixels.

Input Arguments:

• filenames - [cell] cell array of full-path character vectors

• options (optional) - struct with fields:

  • .overlapX - [double] horizontal overlap in percent (default: 0)

  • .overlapY - [double] vertical overlap in percent (default: 0)

Output Arguments:

• layout - struct array per contract (see buildLayoutGrid for field list)

Example - parse a set of MIB2-chopped tiles with 12% overlap:

files = dir('C:\data\chop\*.tif');
layout = utils.stitch.buildLayoutFilenamePattern( ...
    fullfile({files.folder}, {files.name}), struct('overlapX', 12, 'overlapY', 12));

utils.stitch.buildLayoutGrid(filenames, gridOptions)¶

BUILDLAYOUTGRID - Build a tile layout struct from a rectangular grid specification.

Syntax:

layout = utils.stitch.buildLayoutGrid(filenames, gridOptions)

Tiles are natural-sorted before placement. The first tile’s pixel dimensions are read to compute nominal origins; all tiles are assumed to have the same size.

Nominal origin formula (1-based pixels):

x = (col-1) * width * (1 - overlapX/100) + 1 y = (row-1) * height * (1 - overlapY/100) + 1

Input Arguments:

• filenames - [cell] cell array of full-path character vectors for tile files

• gridOptions - struct with fields:

  • .rows - [double] number of grid rows (0 = auto)

  • .cols - [double] number of grid columns (0 = auto)

  • .tileOrder - [char] one of 'Horizontal', 'Horizontal snake', 'Vertical', 'Vertical snake'

  • .overlapX - [double] horizontal overlap in percent (0-90)

  • .overlapY - [double] vertical overlap in percent (0-90)

Output Arguments:

• layout - struct array with fields per contract:

  • .index - [double] 1-based tile index

  • .filename - [char] full path

  • .sliceFiles - [cell] {} for single-file tiles

  • .zLayer - [double] always 1 (single layer)

  • .gridRC - [double] [row col]

  • .nomOrigin - [double] [y x z] 1-based pixel origins

  • .tileSize - [double] [H W D C]

  • .dataClass - [char] MATLAB class string

Example - 2x3 grid of TIF tiles with 10% overlap:

opts.rows = 2; opts.cols = 3;
opts.tileOrder = 'Horizontal'; opts.overlapX = 10; opts.overlapY = 10;
layout = utils.stitch.buildLayoutGrid(myFiles, opts);

utils.stitch.buildLayoutMdoc(mdocPath, options)¶

BUILDLAYOUTMDOC - Build a tile layout (and optionally the stitch) from a SerialEM montage.

Syntax:

layout = utils.stitch.buildLayoutMdoc(mdocPath)
[layout, edges, positions, mdocInfo] = utils.stitch.buildLayoutMdoc(mdocPath, options)

Reads the plain-text .mdoc SerialEM writes beside a montage’s MRC stack. Unlike every other layout source, the tiles are not separate files: each is a SLICE of the one container, addressed through the .sliceIndex field this builder adds to the layout contract and honoured by utils.stitch.makeTileReader().

The .mdoc records three successive stages of the same stitch, exactly as Fibics Atlas splits them across its three XML files, and this function can take any prefix of that chain:

PieceCoordinates (always read) - the NOMINAL montage position of

each tile, in pixels. Becomes nomOrigin.

XedgeDxy / YedgeDxy (options.importEdges) - SerialEM’s own pairwise

seam measurements. Becomes edges, so MIB can solve without re-registering a pixel.

AlignedPieceCoords (options.importPositions) - SerialEM’s FINAL solved

tile positions. Becomes positions, so the mosaic is ready to fuse with nothing recomputed.

Why the nominal placement is only a starting guess. On the reference data the recorded grid is out by up to 5 px - MIB’s own phase correlation and SerialEM’s AlignedPieceCoords agree with each other to ~1 px and both disagree with PieceCoordinates by the same amount. Treat the nominal grid exactly like any other layout source: a rough placement that Measure overlaps refines. Import the edges (or the edges + positions) to keep SerialEM’s own answer instead.

Important

The Y axis is mirrored. MRC stores rows bottom-up and io.loaders.ImodLoader flips them on load, so a tile’s row runs OPPOSITE to its PieceCoordinates Y. Every Y quantity here - nominal origins, solved positions and the edge shifts - is negated together, and the result is normalised to a minimum of 1, which is why the montage’s overall height never enters the arithmetic.

Important

Edge sign convention (load-bearing). XedgeDxy/YedgeDxy are [dx dy] in the montage frame and are stated as the displacement of the LOWER piece, so the offset to add to the nominal step is their negation. In MIB’s [row col] frame the row component picks up a SECOND negation from the Y mirror, which cancels back to a plus: measured = nominal + [edgeDxy(2), -edgeDxy(1), 0]. Verified against phase correlation on real data - both directions agree to ~1.5 px, while the nominal grid is out by 5 px.

Note

A piece stores the edge to its neighbour at HIGHER X/Y, so XedgeDxy is absent on the last column and YedgeDxy on the last row.

Input Arguments:

• mdocPath - [char] full path to the .mdoc file.

• options (optional) - struct with fields:

  • .importEdges - [logical] read XedgeDxy/YedgeDxy into edges (default: false)

  • .importPositions - [logical] read AlignedPieceCoords into positions (default: false). With no edges imported the matching edges are SYNTHESISED from the solved positions by utils.stitch.synthesizeEdgesFromPositions().

  • .imagePath - [char] explicit path to the MRC container; by default it is found from the .mdoc name by utils.stitch.findMdocSidecar().

Output Arguments:

• layout - struct array per the layout contract (see utils.stitch.buildLayoutGrid()) plus .sliceIndex - the 1-based slice of .filename holding this tile. Ordered by row then column.

• edges - struct array of imported/synthesised seams, empty when neither stage was imported. SerialEM records no per-seam confidence, so every imported edge carries .quality = 1 and .valid = true; it is the caller’s utils.stitch.scoreSeams() pass that checks them against the pixels.

• positions - [N x 3 double] solved origins [y x z], or [].

• mdocInfo - struct describing the montage: .file, .imagePath, .pixelSizeUm, .tileHeight / .tileWidth, .mrcMode, .fullMontSize, .numSections and the per-tile .tiles records.

Example - fuse a SerialEM montage exactly as SerialEM stitched it:

opts = struct('importEdges', true, 'importPositions', true);
[layout, edges, positions] = utils.stitch.buildLayoutMdoc( ...
    'D:\Cell1.mrc.mdoc', opts);
canvas = utils.stitch.planCanvas(layout, positions);
mosaic = utils.stitch.fuseInMemory(layout, canvas);

See also utils.stitch.findMdocSidecar, utils.stitch.buildLayoutAtlas, utils.stitch.mrcTargetClass

utils.stitch.buildLayoutPositionFile(positionFilePath, options)¶

BUILDLAYOUTPOSITIONFILE - Build a tile layout struct from a position text file.

Syntax:

layout = utils.stitch.buildLayoutPositionFile(positionFilePath)
layout = utils.stitch.buildLayoutPositionFile(positionFilePath, options)

The position file contains one tile per line with columns filename  X  Y  [Z]. Delimiter is auto-detected among space, tab, and comma; repeated spaces are treated as a single delimiter. Filenames may be relative to the position file’s folder. X, Y and Z are 0-based pixel/slice origins in the file; they are stored 1-based in nomOrigin ([y x z], z in SLICES).

Distinct Z values are additionally ranked into zLayer 1..K in ascending order - zLayer drives layer-adjacency logic (findNeighborPairs), while nomOrigin(3) carries the actual nominal slice coordinate used by the solver and canvas.

A filename entry may name a single image file OR a FOLDER holding the tile’s Z-stack (one image per slice) - this is auto-detected per entry by utils.stitch.resolveTileEntry, so folder Z-stack tiles work here exactly as in the grid and filename-pattern sources.

Input Arguments:

• positionFilePath - [char] full path to the position file

• options (optional) - struct (reserved; no fields used currently)

Output Arguments:

• layout - struct array per contract (see buildLayoutGrid for field list)

Example - load a comma-delimited position file:

layout = utils.stitch.buildLayoutPositionFile('/data/positions.txt');

utils.stitch.canvasBackground(dataClass, canvasColor)¶

CANVASBACKGROUND - Fill value for mosaic pixels that no tile covers.

Syntax:

backgroundValue = utils.stitch.canvasBackground(dataClass)
backgroundValue = utils.stitch.canvasBackground(dataClass, canvasColor)

Tiles almost never tile the canvas rectangle exactly: the solved positions leave a ragged frame around the mosaic, and that frame is filled with this value by every fusion path (options.background of utils.stitch.fuseInMemory() / utils.stitch.fuseStreaming() / utils.stitch.fuseSliceComposite()).

Which value reads as “nothing here” depends on the imaging modality, not on the data: on transmission EM an empty field is BRIGHT, so a zero frame draws a black border around the specimen; on fluorescence it is dark, so zero is right. Hence the choice, and the 'white' default - MIB’s stitching input is predominantly EM. utils.stitch.autocropCanvas() removes the frame altogether; the colour still shows through any gap left by a missing tile.

Input Arguments:

• dataClass - [char] numeric class of the mosaic (canvas.dataClass).

• canvasColor (optional) - [char] 'black' (default) | 'white'.

Output Arguments:

• backgroundValue - [double] fill value, ready for cast(backgroundValue, dataClass).

Example - fuse an EM mosaic on a white canvas:

fuseOptions.background = utils.stitch.canvasBackground(canvas.dataClass, 'white');
imgOut = utils.stitch.fuseInMemory(layout, canvas, fuseOptions);

utils.stitch.computeOverlapRegion(layout, i, j, expandPx)¶

COMPUTEOVERLAPREGION - Nominal overlap rectangle of two tiles in local pixel coords.

Syntax:

[bboxA, bboxB] = utils.stitch.computeOverlapRegion(layout, i, j)
[bboxA, bboxB] = utils.stitch.computeOverlapRegion(layout, i, j, expandPx)

Given the nominal origins and sizes of tiles i and j (from layout), computes the rectangle where the two tiles are expected to overlap and returns it in EACH tile’s own local pixel coordinate system. The rectangle is grown by expandPx on every side to give the phase-correlation search room for the expected positioning jitter, then clamped to the respective tile bounds. Both returned boxes have identical width and height (the intersection extent) so the two crops can be correlated directly.

Coordinate convention: origins are 1-based [y x z] pixel positions of the top-left tile corner in the shared global frame. A tile of size [H W] spans global rows [oy, oy+H-1] and columns [ox, ox+W-1].

Input Arguments:

• layout - [struct array] tile layout with .nomOrigin ([y x z]) and .tileSize ([H W D C]) fields.

• i - [double] index of the first tile.

• j - [double] index of the second tile.

• expandPx (optional) - [double] pixels to expand the overlap on each side to absorb jitter (default: 64).

Output Arguments:

• bboxA - [2x2] overlap rectangle in tile i local coords, [yMin yMax; xMin xMax] (1-based, inclusive).

• bboxB - [2x2] overlap rectangle in tile j local coords, same shape; both boxes have equal height and width.

Note

When the nominal boxes do not overlap at all (even after expansion), the returned rectangles fall back to the full extent of the smaller shared region clamped to both tiles, guaranteeing a non-empty, equal-sized crop.

Example - overlap crop of a horizontal neighbour pair:

[bboxA, bboxB] = utils.stitch.computeOverlapRegion(layout, 1, 2, 32);
readerFcn = utils.stitch.makeTileReader(layout);
cropA = readerFcn(1, bboxA);
cropB = readerFcn(2, bboxB);

utils.stitch.estimateIntensityCorrection(layout, options)¶

ESTIMATETILECORRECTION - Estimate a per-mosaic intensity correction for the tiles.

Syntax:

correction = utils.stitch.estimateIntensityCorrection(layout)
correction = utils.stitch.estimateIntensityCorrection(layout, options)

Measures how the tiles’ brightness differs from one another and returns the correction that removes it, for utils.stitch.makeTileReader() to apply. Nothing is written and no pixels are modified here - the result is a small description that every reader can be built with, so measurement, seam scoring and fusion all see the same corrected pixels.

Why the default is a shading FIELD and not per-tile means. On the data this was written against (TEM, where the beam profile makes one side of every tile brighter than the other), the visible seam steps are NOT a tile-to-tile mean difference. Each tile carries an in-plane gradient, so at a vertical seam one tile’s dark right edge meets its neighbour’s bright left edge and they disagree by 4 - 5 % even when the two tiles’ overall means agree to 0.16 %. Measured on a 3x3 reference montage, as rms mismatch across the 12 seams:

correction	mismatch
None	3.50 %
Match tile means	3.13 %
Flat-field (shared)	1.04 %
flat-field AND mean matching	1.36 %

(Read through utils.stitch.makeTileReader() on the same MRC the numbers come out ~2.1x larger, because the header rescale removes a big offset and that amplifies ratios: 7.50 / 6.73 / 2.31 %. Flat-field (overlap-solved) scores 1.01 % on that scale, i.e. better than twice the shared field.)

Two things follow, and both are deliberate: matching tile means barely helps (it cannot touch a gradient that lives INSIDE each tile), and stacking mean matching on top of a flat field makes it WORSE - tiles genuinely contain different amounts of material, so forcing their means together fights real signal. The methods are therefore alternatives, not layers.

Match tile means is kept because it is the right correction for a DIFFERENT fault: a detector or stain whose response drifts over a long acquisition, where the tiles really do differ by a flat factor.

Note

Estimating reads every tile once. The caller is expected to cache the result (the Stitching controller keeps it in obj.intensityCorrection and drops it when the layout or the method changes) rather than re-deriving it per stage.

Warning

``Flat-field (shared)`` assumes the tiles’ CONTENT averages out. It takes the mean of every tile normalised by its own mean, so whatever survives that average is attributed to illumination. That holds for a real mosaic - many tiles, modest overlap, different specimen under each - but it fails when the tiles are few and heavily overlapping, because then they all show nearly the same thing and the specimen’s own low-frequency structure is absorbed into the field. A 3x3 montage at 10 % overlap estimates cleanly; a 2x2 at 33 % can come out worse than no correction at all.

There is no way to tell the two apart from the field alone. Use Flat-field (overlap-solved) instead: it fits the field to the DISAGREEMENT between tiles where they overlap, and because both tiles image the same specimen there, the specimen cancels exactly. It has no equivalent failure mode and is the better answer whenever seams still show.

Warning

Seams do not observe the tile interior, and a too-flexible fit bows there. Overlaps sample the field only in narrow border strips; the polynomial carries it into the middle with nothing checking it. At degree 4 on real data the fit dipped to 0.87 at the tile centre against 1.12 at the edges - a 37 % span where the whole-tile evidence says 14 %. Dividing that out brightened every tile centre and drew a dark grid along the seams, while the seam residual read 1.01 %, better than any other method. The seam metric is blind to this by construction: it measures only the quantity the fit minimises.

The degree is therefore CHOSEN, by two tests that both have to pass (chooseFieldDegree()):

• cross-validation over held-out seams, which catches a degree fitting noise in the border strips, with a parsimony margin so extra freedom has to earn its place;

• the interior check, which compares the fitted field’s centre against its borders and vetoes anything that bows.

Neither suffices alone, and the reference montage shows why: cross-validation PREFERRED degree 4 (held-out error 8 % better than degree 1) - it is the guard that rejects it at 16.9 % interior deviation. No seam-based number can see that failure, because no seam observes the tile interior.

Do not reach for a ridge on the field coefficients to tame a high degree: it competes with the gain gauge and, past a small weight, flattens the field entirely so the per-tile gains absorb everything - fitting the seams just as well while leaving every tile’s internal gradient in place. Degree is the honest lever.

Note

The other thing ``Flat-field (overlap-solved)`` cannot see. On a regular grid a linear field tilt is indistinguishable from a matching ramp of per-tile gains - both fit the seams identically, so the seams cannot choose between them and the fit needs a rule that comes from outside the data.

Two of the three possible rules are bad. Pinning the FIELD flat leaves every tile’s internal gradient in place and makes the mosaic ripple with the tile period. Pinning the GAINS flat (what fitSeamModel’s ridge does on its own) removes that ripple but leaves the tile means untouched, so the mosaic becomes a monotone STAIRCASE - on real montages this made the overall brightness ramp worse than doing nothing at all.

The rule actually used is the third: fit the shape from the seams, then spend the leftover freedom on making the assembled MOSAIC flat (levelMosaicPlane()). It is provably free - seam residuals are unchanged to the last digit - because it only moves along the direction the seams cannot see. Pass levelMosaic = false to keep the mosaic’s own brightness trend.

Note

``Re-exposure damage`` is a different fault altogether: the specimen is harmed where an earlier tile’s scan already passed, so the LATER tile of every overlap is darker there - a sharp-edged patch in one tile, not a smooth shading shared by all. No field-plus-gain model can represent it: on the reference 2-tile SEM pair Flat-field (overlap-solved) closed the step AT the seam line but left the re-imaged band 10 grey levels dark and a 28-level dark ridge beside it. The model and its measured behaviour are described under solveReexposureDamage() below.

Input Arguments:

• layout - [struct array] tile layout (see utils.stitch.buildLayoutGrid()).

• options (optional) - struct with fields:

  • .method - [char] 'None' (default), 'Flat-field (shared)', 'Flat-field (overlap-solved)', 'Match tile means' or 'Re-exposure damage'.

  • .positions - [N x 3] solved tile origins. Optional for the overlap-solved method: nominal origins are used when absent, which is what lets the method run before anything has been solved (samples are block-averaged, so a few pixels of placement error cannot move a low-order field). Required for 'Re-exposure damage', whose footprint edges are sharp: without positions it warns and returns no correction.

  • .edges - [struct array] seams, for the overlap-solved and re-exposure methods. When given, only .valid ones are fitted, so a seam excluded in the inspector cannot steer the fit either.

  • .damageReach - [double] how far, in pixels, re-exposure damage may extend beyond the earlier tile’s recorded edge (default: 150). Also the margin kept from any OTHER footprint while sampling. Re-exposure only.

  • .damageFarWidth - [double] width, in pixels, of the undamaged strip beyond damageReach used as the outside baseline (default: 50).

  • .edgeGuard - [double] pixels next to a tile’s own border that are not trusted as reference, because every tile carries its own border roll-off (default: 25). Re-exposure only.

  • .minDarkening - [double] a seam’s later tile must be darker by at least this fraction of the overlap mean before damage is assumed there (default: 0.02). Re-exposure only.

  • .previous - [struct] a 'Re-exposure damage' correction estimated earlier for the same layout. When no tile has moved by more than maxReuseShift pixels since, its fitted model is kept and only the footprints are re-placed - no pixel is read.

  • .maxReuseShift - [double] default 20.

  • .polynomialDegree - [double] force the order of the fitted field. [] (default) lets the data choose it, which is the recommended path. A forced degree is still subject to the interior check.

  • .maxPolynomialDegree - [double] highest order the search considers (default: 4).

  • .maxInteriorDeviation - [double] how far the fitted field may differ between the tile centre and its borders before it is rejected (default: 0.15).

  • .degreeSelectionMargin - [double] a higher degree must beat the best held-out error by more than this fraction to be preferred over a simpler one (default: 0.05).

  • .levelMosaic - [logical] after solving, remove the assembled mosaic’s overall brightness plane (default: true). This is free at the seams - it moves only along the direction seam data cannot see - but it does flatten a genuine large-scale trend in the specimen along with an instrumental one. Overlap-solved method only.

  • .samplesPerSeam - [double] roughly how many blocks to sample per seam (default: 800).

  • .robustIterations - [double] IRLS passes (default: 4).

  • .smoothSigma - [double] Gaussian sigma, in pixels, that separates the illumination field from the specimen. Default: 2 % of the shorter tile side (floor 8 px) - large enough that structure averages out, small enough to follow a real beam profile.

  • .showWaitbar / .parentFigure - progress reporting.

  • .cacheSizeBytes - [double] tile-reader cache budget (default 256 MB; tiles are read once each, so a big cache buys nothing here).

Output Arguments:

• correction - struct of PLAIN NUMERICS (no handles), so it survives both parfor serialisation and a JSON round-trip into the project sidecar:

  • .method - [char] the method actually applied

  • .field - [H W] single multiplicative illumination field normalised to mean 1, or [] when the method uses none

  • .gain - [N x 1] per-tile multiplier (all ones when unused)

  • .offset - [N x 1] per-tile additive term, applied BEFORE the gain

  • .tileSize - [H W] the field belongs to, so a reader can refuse a mismatched layout instead of silently mis-scaling

  • .degree / .degreeReport - the chosen field order and what each candidate scored (overlap-solved method)

  • .mosaicPlane - [dx dy] the log-brightness plane removed from the assembled mosaic, [0 0] when none was

  • .damage - [], or for 'Re-exposure damage' a struct:

    • .gain / .offset - the plateau damage, observed = gain * truth + offset where one earlier exposure fully covers the pixel

    • .distance - [1 K] signed distance, in pixels, from a footprint edge (negative inside the earlier tile’s footprint)

    • .profiles - [4 K] damage strength along .distance for the footprint’s left, right, top and bottom side (rows in that order), 1 on the plateau

    • .footprints - [M 5] rows [tile rowMin rowMax colMin colMax]: where an earlier tile’s footprint lies in tile’s own pixel frame

    • .pairs - [P 3] rows [i j direction] for every overlapping pair, direction = +1 when j was imaged after i, -1 the other way round, 0 when no damage is assumed between them

    • .order - [N 1] acquisition rank derived from the seams, NaN for a tile no seam could place

    • .positions - [N 2] the tile origins the footprints were placed with; a caller compares it against its own to know when to re-place them

Example - correct a montage’s shading before fusing it:

correction = utils.stitch.estimateIntensityCorrection(layout, ...
    struct('method', 'Flat-field (shared)'));
mosaic = utils.stitch.fuseInMemory(layout, canvas, ...
    struct('correction', correction));

See also utils.stitch.makeTileReader, utils.stitch.scoreSeams

utils.stitch.estimateOverlap(layout, options)¶

ESTIMATEOVERLAP - Estimate the true grid overlap from the tile images themselves.

Syntax:

result = utils.stitch.estimateOverlap(layout)
result = utils.stitch.estimateOverlap(layout, options)

The user-supplied overlap percentage is often only a guess; the pairwise measurement stage tolerates jitter around the nominal positions but not a systematically wrong overlap. This function recovers the actual overlap MIST-style, exploiting grid redundancy:

1. Neighbour pairs are taken from the grid structure (.gridRC), NOT from the (possibly wrong) nominal origins.

2. Each sampled pair is registered by UNRESTRICTED full-tile phase correlation (zero-padded to twice the tile size, so no shift aliases). The top-K correlation peaks are each verified by normalised cross-correlation of the overlap they imply; the best verified candidate wins (BigStitcher-style peak verification).

3. In a regular grid every x-pair shares the same true step (up to jitter), so the MEDIAN over pairs is a robust step estimate even when individual pairs mismeasure.

Large tiles are downsampled to options.maxDim for speed; the resulting step estimate is coarse (a few px), which is fine - it only repositions the nominal layout for the subsequent tight measurement pass.

Input Arguments:

• layout - [struct array] tile layout with valid .gridRC fields (grid layout source); see utils.stitch.buildLayoutGrid().

• options (optional) - struct with fields:

  • .maxPairsPerDirection - [double] pairs to sample per direction (default: 6)

  • .topK - [double] correlation peaks to verify per pair (default: 5)

  • .maxDim - [double] downsample tiles above this size (default: 1024)

  • .minNcc - [double] minimum verification NCC to accept a pair (default: 0.20)

  • .minOverlapPx - [double] minimum implied overlap extent (default: 16)

  • .colorChannel - [double|’max’] channel used for registration (default: 1)

  • .showWaitbar - [logical] show a Cancelable progress dialog (default: false); reading full-resolution tiles for the sampled pairs can take a noticeable time

  • .parentFigure - [handle] parent for the progress dialog (default: [])

Output Arguments:

• result - [struct] with fields:

  • .overlapX / .overlapY - [double] estimated overlap in percent, NaN when the direction could not be estimated (no pairs / no verified measurement)

  • .stepX / .stepY - [1x2 double] median measured step [dy dx] per direction

  • .madX / .madY - [double] median absolute deviation of the step estimates (px)

  • .numMeasuredX / .numMeasuredY - [double] verified pairs per direction

• cancelled - [logical] true when the user pressed Cancel on the progress dialog before every sampled pair was measured; result is then returned at its all-NaN/zero defaults (a cancelled estimate is discarded, not a partial one).

Example - recover the overlap for a grid built with a wrong guess:

layout = utils.stitch.buildLayoutGrid(files, wrongGridOpts);
est = utils.stitch.estimateOverlap(layout);
gridOpts.overlapX = est.overlapX;  % rebuild layout with the real overlap

utils.stitch.featureShift(cropA, cropB, options)¶

FEATURESHIFT - Feature-based transform estimate between two overlap crops.

Syntax:

[shiftYXZ, quality] = utils.stitch.featureShift(cropA, cropB)
[shiftYXZ, quality, debugInfo] = utils.stitch.featureShift(cropA, cropB, options)

Drop-in alternative to utils.stitch.pairwiseShift() for the RegistrationMethod = 'Feature-based' path. Detects keypoints in each crop (SURF by default), matches descriptors, and RANSAC-fits an estgeotform2d model - pure translation by default, or the model selected by options.transformType; the full fitted matrix is returned in debugInfo.tformA. Unlike phase correlation it does NOT rely on a large textured overlap: a thin shared strip with a handful of matchable blobs is enough, which is why it recovers small (~10%) or unknown overlaps where phase correlation loses the peak. Its weakness is feature-poor or strongly repetitive content (few / ambiguous matches) - there phase correlation still wins.

Sign convention (identical to pairwiseShift, load-bearing). If cropB equals cropA shifted DOWN by dy rows and RIGHT by dx columns (cropB(r,c) ≈ cropA(r-dy, c-dx)), then shiftYXZ = [dy dx 0] - so the composition in utils.stitch.measureAllPairs()/measureOne is unchanged between registration methods. Derivation: estgeotform2d(A, B) returns the map A→B (B ≈ transformPointsForward(tform, A)), whose translation [tx ty] = [T(1,3) T(2,3)] places the feature at cropA(y,x) onto cropB(y+ty, x+tx) ⇒ dy = ty = T(2,3), dx = tx = T(1,3).

Input Arguments:

• cropA - [numeric] reference crop from tile i, [H W] or [H W C].

• cropB - [numeric] moving crop from tile j, same size as cropA.

• options (optional) - struct with fields (all optional). It accepts the same automaticOptions shape produced by controllers.Alignment.defaultAutomaticOptions (per-detector sub-structs detectSURFFeatures …, plus estGeomTransform and rotationInvariance) so the shared settings dialog utils.align.detectorSettingsDlg() can drive it directly:

  • .featureDetector - [char] detector name understood by utils.align.detectFeatures() (default: SURF).

  • .detectSURFFeatures / .detectSIFTFeatures / … - [struct] per-detector parameters (defaults mirror defaultAutomaticOptions).

  • .rotationInvariance - [logical] passed as extractFeatures Upright (default: true - upright descriptors, appropriate for translation).

  • .downsampleFactor - [double ≥ 1] detect on tiles resized by 1/downsampleFactor (faster on big tiles); point locations are scaled back to full resolution before fitting (default: 1).

  • .estGeomTransform - [struct] .MaxNumTrials / .Confidence / .MaxDistance for the RANSAC estgeotform2d fit.

  • .featureMinInliers - [double] minimum RANSAC inliers to trust the fit (default: 8); below this quality = 0.

  • .transformType - [char] estgeotform2d model: 'translation' (default) | 'rigid' | 'similarity' | 'affine'. Non-translation models return their translation component in shiftYXZ and the full matrix in debugInfo.tformA.

  • .allowRotation - [logical] true (default). When false and a non-translation model is requested, the measured edge is constrained to carry no rotation: 'rigid' falls back to the (equivalent) pure translation fit, while 'similarity'/'affine' fits are projected through utils.stitch.projectLinearPart() (R = I branch) and the translation re-estimated over the RANSAC inliers.

  Fields used only by utils.stitch.pairwiseShift() (padPx, expectedShift, searchRadius, window, subpixel) are ignored.

Output Arguments:

• shiftYXZ - [1x3 double] [dy dx dz] correction to the nominal offset (dz is 0; Z handled by the caller).

• quality - [double] in [0,1]: the RANSAC inlier ratio (inliers / matched), 0 when fewer than featureMinInliers inliers survive or the fit fails. Comparable to the phase-correlation quality so the same QualityThreshold gates both methods.

• debugInfo - struct with .numMatched, .numInliers, .status, and .tformA - the full fitted A→B transform as a 3x3 double ([x'; y'; 1] = tformA * [x; y; 1], crop-local pixel coordinates, empty until a fit succeeds). For transformType = 'translation' its linear part is the identity.

Example - recover a known integer shift from a textured crop:

base  = mat2gray(imgaussfilt(randn(256), 1.5));
cropA = base(20:220, 20:220);
cropB = base(15:215, 23:223);          % A is B shifted down 5, left 3
[s, q] = utils.stitch.featureShift(cropA, cropB);
% s ≈ [5 -3 0]

utils.stitch.findAtlasSidecars(atlasFilePath)¶

FINDATLASSIDECARS - Resolve a Fibics Atlas mosaic’s three XML files from any one of them.

Syntax:

sidecars = utils.stitch.findAtlasSidecars(atlasFilePath)

A Fibics Atlas mosaic folder holds up to three XML files that share one base name and describe three successive stages of the same stitch:

• MosaicInfo_<name>.ve-mif - the acquisition record: per-tile grid index and NOMINAL stage position (the rough placement);

• MosaicInfo_<name>.ve-tie - Atlas’s pairwise seam measurements;

• MosaicInfo_<name>.ve-updates - Atlas’s FINAL solved tile positions.

The last two are written only after the mosaic has been stitched in Atlas, so either may be missing. This helper reports which are present, letting the caller offer only the import modes that can actually be honoured.

Any of the three may be passed in - they sit side by side with near-identical names and a user picking one of them means the same mosaic either way. Whatever is passed, .mifPath comes back pointing at the acquisition record, which is the file everything else is read relative to.

.mifPath doubles as the “is this an Atlas input at all?” test: a path whose extension is none of the three (a plain position .txt, an image, …) yields an all-empty struct rather than an error, so a caller can branch on it.

Input Arguments:

• atlasFilePath - [char] full path to any of the mosaic’s three XML files.

Output Arguments:

• sidecars - struct with fields:

  • .mifPath - [char] full path to the .ve-mif, '' if the input is not an Atlas file or the .ve-mif itself is missing

  • .tiePath - [char] full path to the .ve-tie file, '' if absent

  • .updatesPath - [char] full path to the .ve-updates file, '' if absent

Example - branch on whether a picked position file is really an Atlas mosaic:

sidecars = utils.stitch.findAtlasSidecars(selectedFile);
if ~isempty(sidecars.mifPath)
    layout = utils.stitch.buildLayoutAtlas(sidecars.mifPath);
else
    layout = utils.stitch.buildLayoutPositionFile(selectedFile);
end

See also utils.stitch.buildLayoutAtlas, utils.stitch.buildLayoutPositionFile

utils.stitch.findMdocSidecar(montageFilePath)¶

FINDMDOCSIDECAR - Resolve a SerialEM montage’s .mdoc and image container from either one.

Syntax:

sidecar = utils.stitch.findMdocSidecar(montageFilePath)

A SerialEM montage is TWO files: an MRC stack holding every tile as a slice, and a plain-text .mdoc next to it recording where those slices go. The .mdoc is named after the image it describes (Cell1.mrc → Cell1.mrc.mdoc), so either file identifies the pair and a user picking one of them means the same montage. Whatever is passed, .mdocPath comes back pointing at the text file - the one that has to be parsed - and .imagePath at the container the tiles are read from.

.mdocPath doubles as the “is this a SerialEM montage at all?” test: a path whose extension is neither .mdoc nor an MRC image extension yields an all-empty struct rather than an error, so a caller can branch on it exactly as it branches on utils.stitch.findAtlasSidecars().

Important

A non-empty .mdocPath does NOT mean the file describes a mosaic - SerialEM writes the same format for tilt series and single-shot acquisitions, which have no PieceCoordinates and are not stitchable. .isMontage is the flag for that, and it is reported rather than raised so the caller can still route the file to utils.stitch.buildLayoutMdoc() and get one clear error instead of a confusing parse of the wrong file kind.

The three counts describe how far SerialEM got with its own stitch, mirroring the .ve-tie / .ve-updates presence checks on the Atlas side: they let a caller offer only the import modes the file can actually honour.

Input Arguments:

• montageFilePath - [char] full path to either the .mdoc or the MRC image (.mrc, .st, .rec, .ali, .preali, .pre).

Output Arguments:

• sidecar - struct with fields:

  • .mdocPath - [char] full path to the .mdoc; '' when the input is not a SerialEM montage file or the .mdoc is missing

  • .imagePath - [char] full path to the image container; '' if absent

  • .isMrcImage - [logical] the input path is an MRC IMAGE (by extension), whatever came of the lookup. Lets a caller tell “an MRC whose .mdoc is missing” (worth an explanatory error) from “not a SerialEM file at all” (fall through to the next format) - the two are otherwise identical, since both leave .mdocPath empty.

  • .isMontage - [logical] the .mdoc describes a mosaic (it declares Montage = 1 or carries PieceCoordinates)

  • .numTiles - [double] tiles the .mdoc places

  • .numEdges - [double] measured seams (XedgeDxy + YedgeDxy)

  • .numAligned - [double] tiles with a solved AlignedPieceCoords

Example - branch on whether a picked file is a SerialEM montage:

sidecar = utils.stitch.findMdocSidecar(selectedFile);
if ~isempty(sidecar.mdocPath)
    layout = utils.stitch.buildLayoutMdoc(sidecar.mdocPath);
end

See also utils.stitch.buildLayoutMdoc, utils.stitch.findAtlasSidecars

utils.stitch.findNeighborPairs(layout, options)¶

FINDNEIGHBORPAIRS - Find overlapping tile pairs from nominal origins and sizes.

Syntax:

pairs = utils.stitch.findNeighborPairs(layout)
pairs = utils.stitch.findNeighborPairs(layout, options)

Tests all tile pairs for rectangle overlap using their nominal origins and tile sizes. Within-layer pairs are tagged 'x' (primarily side-by-side, abs(dx) >= abs(dy)) or 'y' (primarily top-bottom). Pairs in adjacent Z-layers whose XY footprints overlap are tagged 'z'.

Pairs with overlap smaller than options.minOverlapPx pixels in both dimensions are excluded.

Input Arguments:

• layout - struct array as returned by buildLayoutGrid etc.

• options (optional) - struct with fields:

  • .minOverlapPx - [double] minimum overlap in pixels (default: 16)

Output Arguments:

• pairs - struct array with fields:

  • .i - [double] index of first tile in the pair

  • .j - [double] index of second tile in the pair

  • .direction - [char] 'x', 'y', or 'z'

  • .nominal - [double] [dy dx dz] = layout(j).nomOrigin - layout(i).nomOrigin

Example - find all neighbor pairs in a simple 2x2 grid:

opts.rows = 2; opts.cols = 2;
opts.tileOrder = 'Horizontal'; opts.overlapX = 10; opts.overlapY = 10;
layout = utils.stitch.buildLayoutGrid(files, opts);
pairs  = utils.stitch.findNeighborPairs(layout);

utils.stitch.fuseInMemory(layout, canvas, options)¶

FUSEINMEMORY - Fuse all tiles into a resident mosaic array.

Syntax:

imgOut = utils.stitch.fuseInMemory(layout, canvas)
imgOut = utils.stitch.fuseInMemory(layout, canvas, options)

Allocates the full [H W Z C T] mosaic of canvas.dataClass (initialised to options.background) and fills it by compositing every tile at its planned placement, one output Z-slice at a time via utils.stitch.fuseSliceComposite(), so the transient blending accumulators stay bounded to a single slice even though the final array is resident. This is the fast path for mosaics that fit in RAM; larger jobs use utils.stitch.fuseStreaming().

Input Arguments:

• layout - [struct array] tile layout.

• canvas - [struct] from utils.stitch.planCanvas().

• options (optional) - struct with fields:

  • .blendMode - [char] 'Feather' (default) | 'Average' | 'Max' | 'Min' | 'Overwrite'

  • .background - [double] background fill value (default: 0)

  • .marginPx - [double] feather margin (default: derived from tile size)

  • .cacheSizeBytes - [double] LRU tile-cache budget (default: 2*1024^3)

  • .readerFcn - [function_handle] reuse an existing tile reader (optional)

  • .tileStack - [1 x N] drawing order for 'Overwrite', bottom first (default: [], see utils.stitch.tileDrawOrder())

  • .showWaitbar - [logical] show progress (default: false)

  • .parentFigure - [handle] progress-dialog parent (default: [])

Output Arguments:

• imgOut - [H x W x Z x C x T] fused mosaic of class canvas.dataClass.

Example - feather-blend a solved layout:

opts.blendMode = 'Feather';
imgOut = utils.stitch.fuseInMemory(layout, canvas, opts);

utils.stitch.fuseSliceComposite(layout, canvas, zGlobal, ~, readerFcn, options)¶

FUSESLICECOMPOSITE - Composite all tiles intersecting one output Z-slice.

Syntax:

outSlice = utils.stitch.fuseSliceComposite(layout, canvas, zGlobal, t, readerFcn, options)

Shared per-slice fusion kernel used by utils.stitch.fuseInMemory(), utils.stitch.fuseStreaming() and io.savers.StitchSliceProvider. It builds one output slice [H W C] of canvas.dataClass at global depth zGlobal (1-based) and time t by placing every tile whose Z-band covers zGlobal at its canvas.tilePlacement origin using the requested blend mode. Feather/Average use single-precision accumulators and divide once at the end; Max/Min/Overwrite write in place. When the plan carries per-slice mosaic corrections (canvas.zShifts from the inspector’s Fix Z), every tile on this slice is additionally shifted by canvas.zShifts(zGlobal, :).

When canvas.tforms is present (affine plan from utils.stitch.planCanvas() + utils.stitch.solveGlobalAffine()), each tile whose transform is not an integer translation is RESAMPLED into its warped output footprint with imwarp (bilinear, zero fill); its blend weights are warped identically so the feather follows the warped footprint, and Max/Min/Overwrite only touch pixels the warped tile actually covers. Tiles whose transform reduces to an integer translation - and every tile when there are no canvas.tforms - take the resampling-free integer-placement fast path, which is bit-identical to the translation-only kernel.

Input Arguments:

• layout - [struct array] tile layout.

• canvas - [struct] from utils.stitch.planCanvas().

• zGlobal - [double] 1-based global output slice index.

• t - [double] time frame (1-based).

• readerFcn - [function_handle] tile reader from utils.stitch.makeTileReader().

• options - [struct] with fields:

  • .blendMode - [char] 'Feather' | 'Average' | 'Max' | 'Min' | 'Overwrite'

  • .background - [double] background fill value

  • .marginPx - [double] feather margin (Feather mode; default derived from tile size)

  • .correction - [struct] (optional) the intensity correction the tiles are read with (see utils.stitch.estimateIntensityCorrection()). Only its .damage.order is used here, and only by 'Overwrite': see below.

  • .tileStack - [1 x N] (optional) explicit drawing order, bottom first, set in the seam inspector. 'Overwrite' only.

Note

Which tile wins an ``’Overwrite’`` overlap is decided by :func:`utils.stitch.tileDrawOrder`: an explicit tileStack if set, otherwise - with a re-exposure damage model - the tile imaged FIRST (the later one shows specimen the beam had already hit; the correction evens its brightness but cannot give back destroyed structure), otherwise the highest index. The seam inspector colours the winning tile magenta from the same function. Max/Min/Feather/Average do not depend on the drawing order and are left alone.

Output Arguments:

• outSlice - [H x W x C] fused slice of class canvas.dataClass.

utils.stitch.fuseStreaming(layout, canvas, outputZarrPath, options)¶

FUSESTREAMING - Fuse tiles straight to an OME-Zarr v3 store (out-of-core).

Syntax:

utils.stitch.fuseStreaming(layout, canvas, outputZarrPath)
utils.stitch.fuseStreaming(layout, canvas, outputZarrPath, options)

Primary streaming fusion path. When a single output XY slice fits in options.maxSliceBytes (the common case), a io.savers.StitchSliceProvider is fed to io.savers.Zarr3Saver.saveStream(), which writes level 0 and the full downsampled pyramid with sharding - the whole mosaic is never resident. When even one slice is too large, a chunk-wise fallback creates the zarr array directly (mirroring Zarr3Saver chunk/shard defaults), iterates output chunks, loads only the intersecting tiles (bounded cache), blends, and writes each block with io.zarr.Array.write(); a manual block-downsample pass then adds pyramid levels. Either way io.savers.Zarr3Saver.patchMetadata() stamps the physical bounding box and pixel size so the store reopens as a MIB BigData dataset.

Input Arguments:

• layout - [struct array] tile layout.

• canvas - [struct] from utils.stitch.planCanvas().

• outputZarrPath - [char] destination .zarr3 folder (overwritten).

• options (optional) - struct with fields:

  • .blendMode - [char] 'Feather' (default) | 'Average' | 'Max' | 'Min' | 'Overwrite'

  • .background - [double] background fill value (default: 0)

  • .marginPx - [double] feather margin (default: derived from tile size)

  • .maxSliceBytes - [double] slice-fits threshold in bytes (default: 4*1024^3)

  • .ChunkSize - [1x3 double] level-0 chunk shape (default: [256 256 16])

  • .Compressors - [char] codec name (default: 'zstd')

  • .Levels - [double] explicit pyramid level count (default: auto)

  • .ShardSize - [1x3 double] per-axis chunk multipliers for zarr v3 sharding ([] = no sharding)

  • .DownsampleMethod - [char] pyramid downsampling method (default: 'bilinear')

  • .DownsampleStrategy - [char] 'XY only' (default) | 'Anisotropy-preserving'

  • .cacheSizeBytes - [double] LRU tile-cache budget (default: 2*1024^3)

  • .tileStack - [1 x N] drawing order for 'Overwrite', bottom first (default: [], see utils.stitch.tileDrawOrder())

  • .showWaitbar - [logical] show progress (default: false)

  • .parentFigure - [handle] progress-dialog parent (default: [])

  • .pixSize - [struct] override canvas.pixSize for metadata (optional)

Output Arguments:

(none) - writes outputZarrPath to disk.

Example - stream a mosaic to zarr and reopen as BigData:

utils.stitch.fuseStreaming(layout, canvas, 'C:\out\mosaic.zarr3', ...
    struct('blendMode', 'Feather'));
loader = io.loaders.Zarr3VirtualSetupLoader(struct('datasetMode', 'BigData'));

utils.stitch.fuseToFiles(layout, canvas, outputPath, options)¶

FUSETOFILES - Fuse tiles straight into standard image files (TIF/PNG/AM).

Syntax:

fnOut = utils.stitch.fuseToFiles(layout, canvas, outputPath)
fnOut = utils.stitch.fuseToFiles(layout, canvas, outputPath, options)

Third fusion path beside utils.stitch.fuseInMemory() (a resident MIB dataset) and utils.stitch.fuseStreaming() (an OME-Zarr BigData store): writes the mosaic as ordinary image files, so it can be handed to software that reads nothing else.

A io.savers.StitchSliceProvider is fed to the saver’s saveStream, which is the same slice-by-slice contract the zarr path uses - so for the 2-D sequence formats (TIF, PNG, Amira file sequence) only ONE output slice is ever resident, whatever the mosaic’s depth. Amira Mesh binary writes a single 3-D file and has no streaming writer, so io.savers.BaseSaver.saveStream gathers the volume first: that format needs the whole mosaic in RAM, exactly as In memory does.

Per-slice filenames are the saver’s own (utils.generateSequentialFilename via io.savers.BaseSaver.buildSliceNames): <stem>_001.tif … <stem>_NNN.tif, zero-padded to the digit count the slice total needs, so the files sort in acquisition order. No naming or format dialog is raised - every choice arrives through options.

Input Arguments:

• layout - [struct array] tile layout.

• canvas - [struct] from utils.stitch.planCanvas().

• outputPath - [char] full destination path INCLUDING the extension; for a 2-D sequence it is the stem the numbered files are derived from.

• options (optional) - struct with fields:

  • .Format - [char] io.SaverFactory format string (default: 'TIF format uncompressed (*.tif)')

  • .Saving3DPolicy - [char] '2D sequence' (default) | '3D stack'

  • .blendMode - [char] 'Feather' (default) | 'Average' | 'Max' | 'Min' | 'Overwrite'

  • .background - [double] background fill value (default: 0)

  • .marginPx - [double] feather margin (default: derived from tile size)

  • .correction - [struct] intensity correction from utils.stitch.estimateIntensityCorrection() (default: [])

  • .tileStack - [1 x N] drawing order for 'Overwrite', bottom first (default: [], see utils.stitch.tileDrawOrder())

  • .cacheSizeBytes - [double] LRU tile-cache budget (default: 2*1024^3)

  • .readerFcn - [function_handle] reuse an existing tile reader (optional)

  • .pixSize - [struct] override canvas.pixSize for the file metadata

  • .filename - [char] source name recorded in the metadata (default: outputPath)

  • .showWaitbar - [logical] show progress (default: false)

  • .parentFigure - [handle] progress-dialog parent (default: [])

  • .mibPath - [char] MIB installation directory, for saver dialogs

Output Arguments:

• fnOut - [char] the single written file for a 3-D stack, or a FLAT [cell of char] of the per-slice paths for a 2-D sequence; [] when the saver declined the job (e.g. more colour channels than the format holds).

Example - write a mosaic as a numbered TIF sequence:

opts.Format = 'TIF format uncompressed (*.tif)';
opts.Saving3DPolicy = '2D sequence';
opts.blendMode = 'Overwrite';
fnOut = utils.stitch.fuseToFiles(layout, canvas, 'C:\out\mosaic.tif', opts);
% -> C:\out\mosaic_01.tif ... C:\out\mosaic_NN.tif

See also: utils.stitch.fuseInMemory, utils.stitch.fuseStreaming, io.savers.StitchSliceProvider, io.SaverFactory

utils.stitch.layoutPixSize(layout)¶

LAYOUTPIXSIZE - The acquisition pixel size a layout carries, if it carries one.

Syntax:

pixSize = utils.stitch.layoutPixSize(layout)

Layout sources that read a real acquisition record - SerialEM .mdoc, Fibics Atlas .ve-mif, Bio-Formats OME metadata - know the physical pixel size and attach it to every tile as layout(k).pixSize. Sources that only see a folder of images (Grid, filename pattern, MIB position file) do not, and say so by not having the field at all.

This resolves that into one answer for the whole mosaic, so callers do not each have to re-implement the “is it there, and is it usable?” test. It returns [] rather than a default when nothing is known, which lets the caller leave planCanvas to apply its own 1 um fallback instead of inventing a number here that would then be indistinguishable from a measured one.

Input Arguments:

• layout - [struct array] tile layout.

Output Arguments:

• pixSize - [struct] .x .y .z in µm plus .units, or [] when the layout carries no usable pixel size.

Example

pixSize = utils.stitch.layoutPixSize(layout);
if ~isempty(pixSize); canvasOptions.pixSize = pixSize; end

See also utils.stitch.planCanvas, utils.stitch.buildLayoutMdoc

utils.stitch.loadProject(filePath)¶

LOADPROJECT - Load a stitching project from a JSON sidecar file.

Syntax:

[layout, edges, positions, solverInfo, outputInfo] = utils.stitch.loadProject(filePath)
[layout, edges, positions, solverInfo, outputInfo, tforms] = utils.stitch.loadProject(filePath)
[..., tforms, zSliceFixes, settings] = utils.stitch.loadProject(filePath)

Reads the *.mibstitch.json file saved by utils.stitch.saveProject and reconstructs all stitching state structs. positions is [] when solvedOrigin was not recorded in the file.

Input Arguments:

• filePath - [char] full path to the .mibstitch.json file

Output Arguments:

• layout - struct array with tile fields (index, filename, etc.)

• edges - struct array with pair edge fields (i, j, direction, nominal, …; tform carries the 3x3 pairwise transform when one was measured)

• positions - [double] N-by-3 solved origins [y x z], or []

• solverInfo - struct with solver settings / RMSE

• outputInfo - struct with blend / output settings

• tforms - [N x 1 cell] solved per-tile 3x3 transforms ({} for translation-only projects; tiles saved without one fall back to a pure translation synthesised from solvedOrigin)

• zSliceFixes - [K x 3] per-slice mosaic corrections [z dy dx] from the seam inspector’s Fix Z ([] when none were saved)

• settings - struct of flattened tool settings (schema v3 and newer); an empty struct for older files that carry no settings block

• tileStack - [1 x N] 'Overwrite' drawing order, bottom first, as set in the seam inspector; [] when none was saved (default order)

Example - round-trip save / load:

utils.stitch.saveProject('C:\data\exp.mibstitch.json', layout, edges, [], [], []);
[layout2, edges2, pos, si, oi] = utils.stitch.loadProject('C:\data\exp.mibstitch.json');

utils.stitch.localCorrelate(tileA, tileB, clickXY, currentOffsetYX, options)¶

LOCALCORRELATE - Click-seeded local registration of a tile pair (ROI NCC).

Syntax:

[newOffsetYX, score, confident] = utils.stitch.localCorrelate(tileA, tileB, clickXY, currentOffsetYX)
[newOffsetYX, score, confident, debugInfo] = utils.stitch.localCorrelate(..., options)

The seam inspector’s “human picks WHERE, machine finds EXACTLY” tool (see development/stitching/plan_inspector.md): the user clicks a distinctive spot in the overlap; a small ROI around the click is cut from tile A and matched by normxcorr2 against tile B’s neighbourhood of the corresponding location (within a search radius of the current offset). The NCC peak gives a corrected pair offset with subpixel (parabolic) refinement; the peak height and its prominence over the second-best peak gate a confident flag so a weak/ambiguous match never silently moves a tile.

Offset convention. currentOffsetYX/newOffsetYX are the pair displacement positions(j,1:2) - positions(i,1:2) - tile-A pixel (r, c) corresponds to tile-B pixel (r - dy, c - dx) (the solver/edge.measured convention).

Input Arguments:

• tileA - [numeric] full tile i, 2D or [H W D C] (depth is mean-projected, channel selected per options.colorChannel).

• tileB - [numeric] full tile j, same conventions.

• clickXY - [1x2 double] [x y] click location in tile-A local pixels.

• currentOffsetYX - [1x2 double] current solved [dy dx].

• options (optional) - struct with fields:

  • .roiSize - [double] ROI edge length around the click (default: 128)

  • .searchRadius - [double] search extent around the corresponding location in B, per side (default: 64)

  • .colorChannel - [double|char] channel / 'max' (default: 1)

  • .subpixel - [logical] parabolic peak refinement (default: true)

  • .minPeak - [double] minimum peak NCC for confidence (default: 0.5)

  • .minProminence - [double] minimum peak − second-peak separation (second peak sampled outside a 5-px exclusion zone; default: 0.05)

  • .smoothSigma - [double] sigma of the Gaussian blur of the fallback match; 0 disables the fallback (default: 1). The raw template and search crop are matched first; only when that match is not confident is it repeated on both crops blurred with this sigma (the tiles themselves are never blurred). Pixel noise caps the NCC of correctly aligned raw tiles: on a noisy uint8 EM pair the true offset peaked at 0.27-0.41, under minPeak, and never snapped; with sigma 1 the same spots gave 0.84-0.89 and the same offset. The 'flat-template' test runs on the raw template.

Output Arguments:

• newOffsetYX - [1x2 double] corrected [dy dx]; equals currentOffsetYX when not confident (never a silent bad move).

• score - [double] peak NCC in [-1, 1] (0 when no match ran).

• confident - [logical] peak and prominence above the thresholds.

• debugInfo - [struct] .secondPeak, .templateBBox / .searchBBox ([rowStart rowEnd; colStart colEnd]), .reason (why not confident: '' | 'roi-too-small' | 'flat-template' | 'search-too-small' | 'weak-peak'), .smoothSigma (sigma of the reported match: 0 = raw crops, options.smoothSigma = the blurred fallback, NaN = no match ran). When neither match is confident, the reported score/.secondPeak are those of the last attempt.

Example - snap a seam from a click at a landmark:

[offset, score, ok] = utils.stitch.localCorrelate(tileI, tileJ, [412 88], [2 -158]);
if ok; edge.measured(1:2) = offset; edge.source = 'user'; end

utils.stitch.makeTileReader(layout, options)¶

MAKETILEREADER - Build a cached reader closure that returns tile images on demand.

Syntax:

readerFcn = utils.stitch.makeTileReader(layout)
[readerFcn, isCachedFcn] = utils.stitch.makeTileReader(layout, options)
img = readerFcn(tileIndex)
img = readerFcn(tileIndex, pixelRegion)
tf  = isCachedFcn(tileIndex)

Returns a function handle that loads and caches individual tiles described by the layout struct array. Only the first time point is returned. Single-file tiles are read via io.loadImagesWrapper(); subfolder tiles (with a non-empty .sliceFiles list) are read slice-by-slice and stacked along the depth dimension; MRC-container tiles (with a .sliceIndex, from utils.stitch.buildLayoutMdoc()) are read one slice at a time out of the shared stack. A bounded least-recently-used (LRU) cache holds decoded full tiles so repeated reads (e.g. a tile appearing in several pairwise registrations) do not hit disk again. The cache is a plain cell/struct ring - no containers.Map - so it is safe to serialise into parfor workers.

The optional pixelRegion second argument requests a sub-rectangle of a tile. For single-file TIFF/PNG tiles this uses imread(..., 'PixelRegion', ...) and for MRC-container tiles a ranged getVolume, both avoiding a full decode; for every other case the full tile is loaded (and cached) and then cropped.

Important

An MRC container is opened afresh on every read rather than kept open in the closure: an open file handle cannot cross into a parfor worker, and the header read it costs is negligible beside the pixels. The intensity scaling is likewise taken from the FILE HEADER, so every tile of a montage shares one scale - per-slice statistics would give each tile its own, injecting exactly the intensity mismatch a stitch must not have.

Input Arguments:

• layout - [struct array] tile layout; each element has fields .filename (char, full path - a folder for subfolder tiles, the shared container for MRC tiles), .sliceFiles (cellstr, {} for single-file tiles), .sliceIndex (optional) (1-based slice inside an MRC container) and .tileSize.

• options (optional) - struct with fields:

  • .cacheSizeBytes - [double] LRU cache budget in bytes (default: utils.stitch.tileCacheBudget(), which sizes it to the layout and to the memory this machine has - a fixed 2 GB could not hold even ONE PAIR of large tiles, so the pair view re-decoded both on every revisit)

  • .mibBioformatsCheck - [logical] force the BioFormats reader (default: false)

  • .correction - [struct] intensity correction from utils.stitch.estimateIntensityCorrection(), applied to every tile as it is read (default: none). Building the reader with it is what makes the correction reach measurement, seam scoring and fusion identically - there is no second place pixels enter the pipeline. A .damage model ('Re-exposure damage') is inverted first, per tile, from that tile’s own footprint list; the strength map is built for exactly the region read, so a cropped read and a crop of a full read agree.

Output Arguments:

• readerFcn - [function_handle] img = readerFcn(tileIndex) returns the tile as [H, W, D, C] (first time point); img = readerFcn(tileIndex, pixelRegion) returns a sub-region where pixelRegion = [yMin yMax; xMin xMax].

• isCachedFcn - [function_handle] tf = isCachedFcn(tileIndex): is that tile resident, i.e. would a whole-tile read return immediately? Exists so a caller can tell a free read from one that will stall on disk and put a progress dialog around only the latter - a decode of a large tile is several seconds, and a dialog flashed on every cached read would be worse than none. Never treat it as a promise: a later read can evict the tile.

Example - read two tiles with a shared cache:

readerFcn = utils.stitch.makeTileReader(layout);
tileA = readerFcn(1);
tileB = readerFcn(2);
crop  = readerFcn(1, [10 200; 10 200]);   % sub-region of tile 1

utils.stitch.measureAllPairs(layout, pairs, options)¶

MEASUREALLPAIRS - Measure the residual shift of every neighbour pair.

Syntax:

edges = utils.stitch.measureAllPairs(layout, pairs)
edges = utils.stitch.measureAllPairs(layout, pairs, options)

For each candidate neighbour pair, extracts the nominal overlap crops from both tiles (utils.stitch.computeOverlapRegion() + a cached tile reader), runs utils.stitch.pairwiseShift(), and returns an edges struct array that carries the original pair fields plus the measured offset, its quality and a valid flag. The measured global offset obeys the pairwiseShift sign convention:

edge.measured = edge.nominal + shiftYXZ

Multi-channel tiles are registered on a single channel (options.colorChannel) or their max-projection (options.colorChannel = 'max').

Input Arguments:

• layout - [struct array] tile layout (see utils.stitch.makeTileReader()).

• pairs - [struct array] neighbour pairs; each has .i, .j, .direction ('x'/'y'/'z') and .nominal ([dy dx dz]).

• options (optional) - struct with fields:

  • .expandPx - [double] jitter expansion for the overlap (default: 64)

  • .qualityThreshold - [double] valid = quality >= threshold (default: 0.30)

  • .colorChannel - [double|char] channel index to register on, or 'max' for a max-projection over channels (default: 1)

  • .subpixel - [logical] subpixel refinement in pairwiseShift (default: true)

  • .registrationMethod - [char] 'Phase correlation' (default) or 'Feature-based'; selects utils.stitch.pairwiseShift() or utils.stitch.featureShift() as the per-pair estimator (both share the same sign convention, so the displacement composition is identical).

  • .transformType - [char] 'Translation' (default) | 'Rigid' | 'Similarity' | 'Affine'. Phase correlation can only measure translation, so any non-translation model implies the feature-based estimator regardless of registrationMethod. Non-translation edges additionally carry the full fitted transform in .tform.

  • .allowRotation - [logical] true (default). false constrains every pairwise fit to carry no rotation (forwarded to utils.stitch.featureShift()); pair it with the same option on utils.stitch.solveGlobalAffine() so measurement and solve agree.

  • .preserveEdges - [struct array] previously measured edges whose USER-made fixes (.source = 'user', from the seam inspector) must survive this re-measure: after measuring, any output edge whose (i, j) pair matches a preserved user edge is replaced by it. Without this, one re-measure silently discards a QC session. Preserved user edges whose pair no longer exists in pairs are dropped.

  • .featureOptions - [struct] detector settings forwarded to utils.stitch.featureShift() when registrationMethod is 'Feature-based' (detector type, per-detector params, downsampling, RANSAC; the automaticOptions shape). Ignored for phase correlation.

  • .cacheSizeBytes - [double] LRU tile-cache budget (default: sized to the layout by utils.stitch.tileCacheBudget(), divided by the pool size on the parfor path where each worker caches separately)

  • .useParallel - [logical] measure pairs with parfor (default: false)

  • .showWaitbar - [logical] show a progress dialog (default: false); only the sequential path (useParallel = false) makes it Cancelable - a parfor batch cannot poll the dialog mid-iteration

  • .parentFigure - [handle] parent for the progress dialog (default: [])

Output Arguments:

• edges - [struct array] one per pair with fields .i .j .direction .nominal (copied) plus .measured ([dy dx dz]), .quality ([0,1]), .valid (logical), .tform - the tile-local A→B transform as a 3x3 double in xy pixel coordinates ([x_j; y_j; 1] = tform * [x_i; y_i; 1]), filled only when a non-translation transformType was fitted ([] otherwise; the translation solver reads .measured alone) - and the seam-inspector bookkeeping fields .source ('auto' here; 'user'/ 'confirmed' are set by the inspector) and .seamScore ([] here; filled by utils.stitch.scoreSeams()). When cancelled partway, only the pairs measured before the cancel are included - the caller must check cancelled rather than assume edges covers every pair in pairs.

• cancelled - [logical] true when the user pressed Cancel on the progress dialog before all pairs were measured; false otherwise (always false when showWaitbar is off or no dialog was shown).

Example - measure all pairs sequentially:

opts.qualityThreshold = 0.3;
edges = utils.stitch.measureAllPairs(layout, pairs, opts);
solvable = edges([edges.valid]);

utils.stitch.moveInTileStack(stack, tileIdx, action, positions, tileSizes)¶

MOVEINTILESTACK - Move one tile in the Overwrite drawing order.

Syntax:

stack = utils.stitch.moveInTileStack(stack, tileIdx, action, positions, tileSizes)

stack lists the tiles bottom first (see utils.stitch.tileDrawOrder()): the last one wins every overlap it takes part in. The four actions the seam inspector offers:

• 'top' / 'bottom' - to the very end / start of the stack.

• 'up' / 'down' - one step past the nearest tile above / below it that it actually overlaps. Stepping past a tile it does not touch changes no pixel of the mosaic, so it would be an action with no visible effect; skipping those makes every press count.

An action that changes nothing ('top' on the top tile, 'up' with no overlapping tile above) returns the stack unchanged, which is how the caller tells which actions are worth offering.

Input Arguments:

• stack - [1 x N double] tile indices, bottom first.

• tileIdx - [double] the tile to move.

• action - [char] 'top' | 'up' | 'down' | 'bottom'.

• positions - [N x 3 double] tile origins [y x z] (solved or nominal).

• tileSizes - [N x 3+ double] [H W D ...] per tile, as in layout.tileSize.

Output Arguments:

• stack - [1 x N double] the new order, bottom first.

Example - put tile 4 below tile 2:

stack = utils.stitch.moveInTileStack(stack, 4, 'down', positions, ...
    reshape([layout.tileSize], 4, []).');

See also utils.stitch.tileDrawOrder

utils.stitch.mrcTargetClass(mrcMode)¶

MRCTARGETCLASS - Pixel class a stitched MRC montage is read into, from the MRC mode.

Syntax:

[targetClass, needsScaling] = utils.stitch.mrcTargetClass(mrcMode)

MRC stores pixels as unsigned bytes, SIGNED 16-bit integers or 32-bit floats, while MIB works in unsigned integer classes only. This function is the single place that decides which class a montage’s tiles become, so the two callers - utils.stitch.buildLayoutMdoc() (which records layout.dataClass, and therefore what utils.stitch.planCanvas() allocates) and utils.stitch.makeTileReader() (which does the actual conversion) - cannot drift apart.

Why every scaled mode lands on uint16 rather than the widest class that fits. io.loaders.ImodLoader picks the target class from the header’s density RANGE, which sends a typical float montage to uint32 and leaves MIB to stretch it to uint16 on load. A mosaic is built once and kept, so the intermediate uint32 would only double the canvas and every fused slice for a range that is going to be squeezed anyway.

Important

Scaling uses the range from the FILE HEADER, never per-slice statistics. Every tile of a montage is a slice of one container, so a per-slice range would give each tile its own intensity scale - injecting exactly the tile-to-tile mismatch a stitch is supposed to be free of.

Input Arguments:

• mrcMode - [double] MRC mode field: 0 = unsigned byte, 1 = int16, 2 = float32, 6 = uint16 (3/4 complex and 16 RGB are treated as scaled modes).

Output Arguments:

• targetClass - [char] MATLAB class the tiles are returned as.

• needsScaling - [logical] true when pixels must be mapped from the header’s [minDensity maxDensity] onto the full range of targetClass; false when the stored values are already unsigned and pass through.

Example - allocate a canvas for a float32 montage:

[tileClass, mustScale] = utils.stitch.mrcTargetClass(2);
% tileClass = 'uint16', mustScale = true

See also utils.stitch.buildLayoutMdoc, utils.stitch.makeTileReader

utils.stitch.naturalSortFiles(cellstrFiles)¶

NATURALSORTFILES - Sort a cell array of filenames in natural (alphanumeric) order.

Syntax:

sortedFiles = utils.stitch.naturalSortFiles(cellstrFiles)

Natural sort treats embedded digit sequences as numbers, so that tile2.tif sorts before tile10.tif. Comparison is case-insensitive and operates on the full path string.

Input Arguments:

• cellstrFiles - [cell] cell array of character vectors (filenames or full paths)

Output Arguments:

• sortedFiles - [cell] input cell sorted in natural/alphanumeric order

Example - sort a mixed list of tile names:

files = {'tile10.tif', 'tile2.tif', 'tile1.tif'};
sorted = utils.stitch.naturalSortFiles(files);
% sorted = {'tile1.tif', 'tile2.tif', 'tile10.tif'}

utils.stitch.pairwiseShift(cropA, cropB, options)¶

PAIRWISESHIFT - Windowed FFT phase correlation between two overlap crops.

Syntax:

[shiftYXZ, quality] = utils.stitch.pairwiseShift(cropA, cropB)
[shiftYXZ, quality, debugInfo] = utils.stitch.pairwiseShift(cropA, cropB, options)

Estimates the residual translation between two same-size overlap crops taken from a pair of neighbouring tiles at their NOMINAL relative position. A Hann window is applied to suppress FFT edge effects, then normalised cross-power spectrum (phase correlation) locates the peak. Optional parabolic subpixel refinement fits a quadratic to the three samples straddling the integer peak on each axis. Registration is 2D; shiftYXZ(3) (dz) is returned as 0 in Phase 1 (Z handled by the caller via facing-slice / mean-projection crops).

Sign convention (load-bearing). If cropB equals cropA shifted DOWN by dy rows and RIGHT by dx columns (i.e. cropB(r,c) ≈ cropA(r-dy, c-dx)), then shiftYXZ = [dy dx 0]. Verified in tests/utils/StitchCoreTest.m.

When the crops start at local positions bboxA(:,1) / bboxB(:,1) (from utils.stitch.computeOverlapRegion()), the tile displacement follows as

% P_j - P_i in the shared global frame:
measured = (bboxA(:,1) - bboxB(:,1))' - [dy dx]

This composition (including the sign NEGATION of the raw shift) lives in utils.stitch.measureAllPairs/measureOne.

Input Arguments:

• cropA - [numeric] reference crop from tile i, [H W] or [H W C].

• cropB - [numeric] moving crop from tile j, same size as cropA.

• options (optional) - struct with fields:

  • .subpixel - [logical] enable parabolic subpixel refinement (default: true)

  • .window - [logical] apply a separable Hann window (default: true)

Output Arguments:

• shiftYXZ - [1x3 double] [dy dx dz] correction to the nominal offset (dz is 0 in Phase 1).

• quality - [double] peak prominence in [0, 1]: the phase-correlation peak height relative to the surrounding field. Textured overlaps score high (> 0.5), flat/noise overlaps score low (< 0.1).

Example - recover a known integer shift:

cropA = single(imfilter(randn(128), fspecial('gaussian', 9, 2)));
cropB = circshift(cropA, [5 -3]);   % down 5, left 3
[s, q] = utils.stitch.pairwiseShift(cropA, cropB);
% s ≈ [5 -3 0], q > 0.5

utils.stitch.planCanvas(layout, positions, options)¶

PLANCANVAS - Compute the fused mosaic canvas from solved tile positions.

Syntax:

canvas = utils.stitch.planCanvas(layout, positions)
canvas = utils.stitch.planCanvas(layout, positions, options)

Turns the fractional solved origins into an integer placement plan for the output mosaic: it shifts all origins so the minimum origin lands at pixel 1, rounds to integer tile placements (keeping the fractional remainder as a per-tile subpixel residual for later resampled placement), and derives the total canvas size and physical bounding box. All three axes use the SOLVED positions - positions(:,3) is a slice coordinate, so overlapping or jittered Z-stacks land where the global solve put them (tiles within a 2D layer share one z by construction of the within-layer dz constraints).

When per-tile transforms from utils.stitch.solveGlobalAffine() are passed in options.tforms, the in-plane extent comes from the WARPED tile footprints instead: the corner points of every tile are pushed through its transform, the union bounding box defines the canvas, and the transforms are re-expressed in canvas coordinates (tile-local xy → canvas xy, pixel centres at integers) in canvas.tforms together with each tile’s integer output bounding box in canvas.tileBounds. The z axis keeps the translation placement logic. Without options.tforms the behaviour is bit-identical to the translation-only planner.

Input Arguments:

• layout - [struct array] tile layout with .tileSize ([H W D C]), .zLayer and .dataClass.

• positions - [N x 3 double] solved [y x z] origins (fractional).

• options (optional) - struct with fields:

  • .pixSize - [struct] .x .y .z physical voxel size (default: 1 µm iso)

  • .numChannels - [double] output channel count (default: from tileSize(4))

  • .numFrames - [double] output time frames (default: 1)

  • .tforms - [N x 1 cell] per-tile 3x3 tile-local→global xy transforms from utils.stitch.solveGlobalAffine(); enables the warped-footprint plan (default: none - translation placement)

  • .zSliceFixes - [K x 3] per-slice mosaic corrections [z dy dx] from the seam inspector’s Fix Z: every output slice >= z shifts in-plane by [dy dx] (cumulative over rows). Produces canvas.zShifts and grows the canvas so nothing is clipped (default: none)

  • .autocrop - [logical] trim the ragged background frame the solved positions leave around the mosaic, by handing the finished plan to utils.stitch.autocropCanvas() (default: false)

Output Arguments:

• canvas - [struct] with fields:

  • .size - [1x5] [H W Z C T] output dimensions

  • .tilePlacement - [N x 3] integer 1-based [y x z] origins

  • .subpixelResidual - [N x 3] fractional part discarded by rounding

  • .dataClass - [char] numeric class of the mosaic

  • .boundingBox - [1x6] [xmin xmax ymin ymax zmin zmax] physical extent

  • .pixSize - [struct] the pixel size used

  • .zLayers - [vector] sorted distinct z-layer ids

  • .tforms - [N x 1 cell] tile-local→CANVAS xy transforms (only when options.tforms was given; canvas world frame == intrinsic pixels)

  • .tileBounds - [N x 4] integer output bbox [y0 y1 x0 x1] of each warped tile footprint, clipped to the canvas (only with options.tforms)

  • .zShifts - [Z x 2] integer extra [dy dx] applied to every tile on that output slice by the fusers (only with options.zSliceFixes; baseline-shifted so all entries are >= 0 and fit the grown canvas)

  • .cropRect - [1x4] [y0 y1 x0 x1] kept region in the uncropped frame (only with options.autocrop, and only when a fully covered region exists)

Example - plan a canvas at 20 nm isotropic:

opts.pixSize = struct('x', 0.02, 'y', 0.02, 'z', 0.02);
canvas = utils.stitch.planCanvas(layout, positions, opts);
fprintf('Canvas: %d x %d x %d\n', canvas.size(1), canvas.size(2), canvas.size(3));

utils.stitch.projectLinearPart(L, transformType, allowRotation)¶

PROJECTLINEARPART - Closest member of a 2D transform group to a linear part.

Syntax:

Lp = utils.stitch.projectLinearPart(L, transformType, allowRotation)

Projects a 2x2 linear part onto the requested transform group in the Frobenius norm, via the polar decomposition L = R * P (rotation x symmetric stretch, SVD-based). This is the shared projection used by utils.stitch.solveGlobalAffine() (per-tile, after the linear global solve) and utils.stitch.featureShift() (per-edge, when the rotation lock is on) - the R = I branch implements the AllowRotation = false constraint (see development/stitching/plan_transforms.md).

Input Arguments:

• L - [2x2 double] linear part of an affine transform.

• transformType - [char] 'Rigid' | 'Similarity' | 'Affine' (case-insensitive).

• allowRotation - [logical] false locks the rotation factor to identity.

Output Arguments:

• Lp - [2x2 double] the projected linear part:

  • Rigid → R (or I when rotation is disallowed)

  • Similarity → s*R with s = trace(R'*L)/2, the Frobenius-optimal scalar (or s*I with s = trace(L)/2)

  • Affine → L unchanged when rotation is allowed; the symmetric stretch part P (scale/shear, no rotation) when disallowed.

Example - strip the rotation out of a fitted linear part:

Lp = utils.stitch.projectLinearPart(tform.A(1:2,1:2), 'Affine', false);

utils.stitch.rankSeams(edges)¶

RANKSEAMS - Worst-first review order for an already-scored edge set.

Syntax:

ranking = utils.stitch.rankSeams(edges)

Orders edges the way the seam inspector reviews them: pruned (valid = false) edges first, then by seamScore ascending with unscored / NaN scores leading (no overlap at the solved placement is the worst thing an edge can be).

Split out of utils.stitch.scoreSeams() because it reads NO PIXELS - it is a pure function of the edge fields. That is what lets a consumer re-derive the order after an edge is excluded, or adopt an edge set that was already scored, without paying for a second pass over the overlaps (which on a large mosaic is measured in minutes). utils.stitch.scoreSeams() calls this at the end, so the two orders cannot drift.

Input Arguments:

• edges - [struct array] with .valid and .seamScore. Either field may be missing or empty; missing valid counts as valid, an empty or NaN score sorts worst.

Output Arguments:

• ranking - [1 x M double] edge indices, worst first.

See also utils.stitch.scoreSeams

utils.stitch.resolveTileEntry(entryPath)¶

RESOLVETILEENTRY - Resolve a tile source path to its slices, size and class.

Syntax:

[sliceFiles, tileSize, dataClass] = utils.stitch.resolveTileEntry(entryPath)

A tile source is either a single image file or a FOLDER holding the tile’s Z-stack (one image per slice). This helper is the single place that decides which, so every layout builder (utils.stitch.buildLayoutGrid(), utils.stitch.buildLayoutPositionFile(), utils.stitch.buildLayoutFilenamePattern()) supports folder Z-stack tiles identically: whether a tile is a folder is a property of the SOURCE, orthogonal to how the tiles are ARRANGED (grid, position file, filename pattern).

Input Arguments:

• entryPath - [char] path to a single image file or a tile folder.

Output Arguments:

• sliceFiles - [cell] {} for a single-file tile, otherwise the natural-sorted list of per-slice image paths in the folder.

• tileSize - [1x4 double] [H W D C]; D is the slice count for a folder tile.

• dataClass - [char] numeric class of the tile pixels.

Example - resolve a folder Z-stack tile:

[slices, sz, cls] = utils.stitch.resolveTileEntry('C:\data\tile_01');
% numel(slices) == sz(3)

utils.stitch.saveProject(filePath, layout, edges, positions, solverInfo, outputInfo, tforms, zSliceFixes, settings, tileStack)¶

SAVEPROJECT - Save a stitching project to a JSON sidecar file.

Syntax:

utils.stitch.saveProject(filePath, layout, edges, positions, solverInfo, outputInfo)
utils.stitch.saveProject(filePath, layout, edges, positions, solverInfo, outputInfo, tforms)
utils.stitch.saveProject(filePath, layout, edges, positions, solverInfo, outputInfo, tforms, zSliceFixes)
utils.stitch.saveProject(..., tforms, zSliceFixes, settings)

Writes all stitching state to <name>.mibstitch.json for reproducibility and later use by the QC / seam checker. The file is human-readable JSON encoded with jsonencode.

Input Arguments:

• filePath - [char] full path for the output JSON file (the .mibstitch.json extension is appended if not already present)

• layout - struct array as returned by the layout builders

• edges - struct array of measured/validated pair edges (may be [])

• positions - [double] N-by-3 array of solved origins [y x z] (may be [] if not yet solved)

• solverInfo - struct with solver settings and RMSE (may be [])

• outputInfo - struct with blend mode, output path, canvas size (may be [])

• tforms (optional) - [N x 1 cell] solved per-tile 3x3 affine transforms from utils.stitch.solveGlobalAffine() (stored per tile as solvedTform); pass {}/omit for translation-only projects

• zSliceFixes (optional) - [K x 3] per-slice mosaic corrections [z dy dx] from the seam inspector’s Fix Z; pass []/omit for none

• settings (optional) - struct of flattened tool settings (one scalar / char / logical per BatchOpt field, plus the nested FeatureOptions) as produced by controllers.Stitching.collectProjectSettings(). Stored under project.settings so Load project can restore the whole dialog, or reuse the parameters alone on a different set of tiles. Omit for none.

• tileStack (optional) - [1 x N] the 'Overwrite' drawing order set in the seam inspector, bottom first (see utils.stitch.tileDrawOrder()); stored as project.tileStack. []/omit when the user never set one - the default order is then re-derived on load rather than frozen into the file.

Example - save after solving:

utils.stitch.saveProject('C:\data\experiment.mibstitch.json', ...
    layout, edges, positions, solverInfo, outputInfo);

utils.stitch.scoreSeams(layout, edges, positions, options)¶

SCORESEAMS - Pixel-based seam quality of every edge at the SOLVED positions.

Syntax:

[edges, ranking] = utils.stitch.scoreSeams(layout, edges, positions)
[edges, ranking] = utils.stitch.scoreSeams(layout, edges, positions, options)
[edges, ranking, cancelled] = utils.stitch.scoreSeams(layout, edges, positions, options)

For every edge, reads the two tiles’ overlap strips AT THE SOLVED POSITIONS and stores their zero-mean normalised cross-correlation in edges(k).seamScore (range [-1, 1]; higher = better seam). This is the inspector’s primary ranking metric (see development/stitching/plan_inspector.md): a confidently WRONG pairwise measurement - RANSAC or phase correlation locked one period off on repetitive content - satisfies the solver perfectly on a chain-like graph (zero residual), but its pixels do not agree at the solved placement, so only re-checking actual pixels catches it. Edges whose tiles do not overlap at all at the solved positions get seamScore = NaN (worst possible).

3D (z-stack tiles): the strips are aligned in depth by the solved dz and scored as the MEAN of per-slice 2D NCCs over the overlapping slab (the same metric measureAllPairs uses to choose dz - a volumetric NCC would inherit the thick-slab bias and stay high at wrong offsets). No slab overlap at the solved dz scores NaN. Cross-layer (direction == 'z') edges are additionally re-scored at dz ± dzScanRadius: when a neighbouring offset beats the solved one by dzHintMargin the difference is stored in edges(k).dzHint (e.g. +2 = the pixels prefer dz + 2); dzHint = 0 means the solved dz is the local optimum.

The returned ranking orders edges for worst-first review: pruned (valid = false) edges first, then by seam score ascending with NaN scores leading.

Non-translation solves: the strips are cut by the TRANSLATION between the solved origins; the residual rotation/scale of an otherwise good seam lowers its absolute score slightly, but the ranking (relative scores) remains meaningful. Transform-aware strip warping is a later refinement.

Input Arguments:

• layout - [struct array] tile layout (.tileSize, reader fields).

• edges - [struct array] from utils.stitch.measureAllPairs().

• positions - [N x 3 double] solved [y x z] origins.

• options (optional) - struct with fields:

  • .colorChannel - [double|char] channel to score on, or 'max' (default: 1; same semantics as measureAllPairs)

  • .maxStripPx - [double] strips longer than this are downsampled before correlation (default: 1024)

  • .maxScoreSlices - [double] at most this many z-slices of the overlap slab are correlated, evenly sampled (default: 16)

  • .dzScanRadius - [double] cross-layer edges are re-scored at dz ± radius to fill dzHint; 0 disables (default: 2)

  • .dzHintMargin - [double] a neighbouring dz must beat the solved one by this much before it is hinted (default: 0.05)

  • .cacheSizeBytes - [double] LRU tile-cache budget (default: 2*1024^3)

  • .readerFcn - [function_handle] reuse an existing tile reader (optional)

  • .showWaitbar / .parentFigure - progress dialog (default: off). The dialog is Cancelable: scoring re-reads every overlap from disk, which on a large mosaic is the slowest advisory step in the tool.

Output Arguments:

• edges - the input edges with .seamScore and .dzHint filled in. When cancelled, EVERY score is cleared (seamScore = [], dzHint = 0), including the ones already computed and any the input carried. A partial set would be worse than none: the quality chip and the inspector’s ranking both take a MINIMUM over the scored seams, so scoring half of them and stopping would rate the mosaic on its better half and silently ignore the seams nobody looked at. Whatever the input scores described, it was a different placement - this function is only ever run after the positions changed.

• ranking - [1 x M double] edge indices in worst-first review order; plain input order 1:M when cancelled (there is nothing to rank by).

• cancelled - [logical] true when the user pressed Cancel on the progress dialog before all edges were scored; false otherwise (always false when showWaitbar is off or no dialog was shown).

Example - score and list the three worst seams:

[edges, ranking] = utils.stitch.scoreSeams(layout, edges, positions);
worst = edges(ranking(1:3));

utils.stitch.solveGlobalAffine(layout, edges, options)¶

SOLVEGLOBALAFFINE - Global per-tile affine transforms from pairwise edges (linear LS).

Syntax:

[tforms, positions, stats] = utils.stitch.solveGlobalAffine(layout, edges)
[tforms, positions, stats] = utils.stitch.solveGlobalAffine(layout, edges, options)

Affine generalisation of utils.stitch.solveGlobalLeastSquares(). Every tile t gets a 2D affine map from 1-based tile-local pixel coordinates to global coordinates, x_global = L_t * x_local + p_t (L_t 2x2, p_t 2x1, xy order). A pairwise edge measured as the tile-local map x_j = M * x_i + c (edge.tform, from the feature-based estimator) pins the two tiles’ maps together: a physical point seen at x_i in tile i and x_j in tile j must land on the same global point, which for all overlap points splits into

• L_i = L_j * M (4 scalar rows)

• p_i - p_j - L_j * c = 0 (2 scalar rows)

Both are LINEAR in the stacked unknowns {L_t, p_t}, so the global solve stays one sparse weighted least-squares system (the BigStitcher affine model) - no nonlinear optimiser. Edges without a stored transform (phase-correlation or translation-model measurements) participate with M = I and c synthesised from .measured, so mixed edge sets are fine.

The linear-part rows are scaled by a characteristic tile extent so that a dimensionless error in L is weighted like the pixel error it causes at the tile edge; without this the position rows (pixels) would numerically dominate the linear rows (~1).

Gauge, springs and connectivity mirror the translation solver: tile 1 is anchored at L = I / its nominal origin; pruned edges touching a component disconnected from the anchor re-enter as weak springs toward nominal; disconnected tiles get a bridge spring; every tile gets a tiny self-spring (position toward nominal, linear part toward identity) for rank.

Rigid / Similarity - solved by projection, not by a nonlinear optimiser. options.transformType restricts the model: after the (always-linear) affine solve, each tile’s linear part is polar-decomposed L = R * S (rotation x symmetric stretch, via SVD) and replaced by the closest member of the requested group - R for Rigid, s*R (s = mean singular value, the Frobenius-optimal scalar) for Similarity. options.allowRotation = false additionally locks the rotation factor to identity (R = I branch of the same decomposition): Rigid degenerates to pure translation, Similarity to scale + translation, Affine keeps scale/shear but no rotation. After the projection the tile translations are RE-SOLVED with the projected linear parts held fixed - each edge then contributes the known offset pos_j - pos_i = (L_j - L_i)*[1;1] - L_j*c, which is exactly the translation solver’s row shape, so the refinement reuses utils.stitch.solveGlobalLeastSquares() (same springs/anchoring).

The z axis composes additively regardless of the in-plane model (2D affine acts within a slice), so z origins are solved with the existing scalar path (utils.stitch.solveGlobalLeastSquares()) and merged into positions.

Input Arguments:

• layout - [struct array] tile layout with .nomOrigin ([y x z]).

• edges - [struct array] from utils.stitch.measureAllPairs(); uses .i .j .measured .nominal .quality .valid and (when present) .tform (3x3 double, tile-local xy map i → j).

• options (optional) - struct with fields:

  • .springWeight - [double] weight of re-added pruned/bridge springs (default: 0.10)

  • .nominalSpringWeight - [double] weight of the per-tile position self-spring (default: 0.001; rank guard only - real weight biases the solution)

  • .identitySpringWeight - [double] weight of the per-tile linear-part self-spring toward identity (default: 0.001). A mild prior against scale/shear drift accumulating across the mosaic.

  • .transformType - [char] 'Affine' (default) | 'Rigid' | 'Similarity' - the group each tile’s linear part is projected onto after the linear solve (see above).

  • .allowRotation - [logical] true (default) permits per-tile rotation; false locks the rotation factor to identity in the projection (robustness prior for stage-tiled data that cannot rotate).

  • .userEdgeWeight - [double] weight of USER-fixed edges (edge.source = 'user'; default: 5.0) - never pruned, dominate conflicting automatic edges (see the translation solver’s doc).

Output Arguments:

• tforms - [N x 1 cell] per-tile 3x3 doubles mapping 1-based tile-local xy to global xy: [x'; y'; 1] = tforms{t} * [x; y; 1]. The anchor tile has an identity linear part and its nominal origin.

• positions - [N x 3 double] solved [y x z] origins - the global coordinate of each tile’s pixel (1,1) (L_t*[1;1] + p_t in yx order, z from the scalar solve). Collapses to the translation-solver result when every edge is a pure translation.

• stats - [struct] same shape as the translation solver’s: .residuals (M x 3, position mismatch per valid edge evaluated at the source tile’s centre), .rmse (1x3), .rmseTotal, .nPruned, .nComponents, .anchorComponent, .disconnectedTiles.

Example - solve an affine-jittered grid:

[tforms, positions, stats] = utils.stitch.solveGlobalAffine(layout, edges);
fprintf('RMSE = %.3f px\n', stats.rmseTotal);