matlab-process-large-images

name	matlab-process-large-images
description	Patterns for using blockedImage to process large images, harness parallel compute for image processing, and write custom adapters. Use when writing code that creates, processes, or visualizes blockedImage objects, when implementing images.blocked.Adapter subclasses, or when a user needs help with large image data. Always use this skill when working with TIFF,GeoTIFF, .svs, .ndpi, .czi or other WSI, satellite imagery or microscopy volume image formats.
license	MathWorks BSD-3-Clause
metadata	{"author":"MathWorks","version":"1.0"}

When to Use

• Writing code that creates, processes, or visualizes blockedImage objects
• Implementing images.blocked.Adapter subclasses for custom file formats
• User needs help with large image data that doesn't fit in memory
• User wants to leverage parallel computing for image processing tasks
• Working with whole-slide images, geospatial rasters, 3D microscopy volumes, or tiled mosaics
• Any time a user is working with a TIFF file — load this skill and check the file first with imfinfo. Use blockedImage if the file is large (dimensions > 10,000 pixels in any axis), has multiple channels beyond RGB, or contains multiple IFDs (could be pyramid levels, time series, or Z-stack slices)

When NOT to Use

• Image fits comfortably in memory and no parallel processing is needed — use standard imread/imwrite workflows
• General image display without blockedImage — use the matlab-display-image skill instead
• Deep learning model architecture or training loop design — this skill covers data preparation with blockedImageDatastore, not network design

Large Image Processing with blockedImage

Reference files (read on demand when the topic comes up):

• references/visualization.md — imageshow, overlays, synced views, overview+detail, annotations
• references/training-and-labeling.md — Image Labeler, labeled blockedImages, training datastores, DL inference
• references/3d-volumes.md — 3D blocked images, 3D inference, 3D training data
• references/geospatial.md — GeoTIFF, shapefiles, coordinate systems, spatial referencing
• references/data-discovery.md — 7-step data discovery workflow, custom adapters, format conversion adapter selection, virtual composite adapters, adapter performance tips
• references/specialized-workflows.md — writing block-by-block, coordinate systems, mask-based block selection, ExtraImages for pixel-level masking, in-place resume
• references/whole-slide-images.md — importing WSI files via OpenSlide/Bio-Formats, choosing between libraries, display, downstream processing (segmentation, classification)

When to use blockedImage

blockedImage is not only for data that exceeds RAM. Even when the data fits in memory, wrapping it as a blockedImage lets you harness a parallel cluster: each worker can independently read and write different blocks, enabling parallel partial reads and parallel writes into separate file blocks.

Always ask the user what their intent with the data is. Ask about visualization, processing, DL training, archival, and whether they want to leverage parallel compute — even if the data fits in RAM.

When the user points you at a directory or multi-file dataset

ALWAYS read references/data-discovery.md and follow its 7-step workflow before writing any code. The workflow covers: inspecting file layout, parsing naming conventions, understanding user intent, choosing between single blockedImage vs. collection, assessing random access capabilities, checking for gaps/overlaps, and recommending an approach (built-in adapter, format conversion, or custom adapter).

Construction

% Read-only (auto-selects adapter; reads only metadata, fast)
bim = blockedImage("large_image.tif");

% Multiple images at once
[bim1, bim2] = blockedImage(["img1.tif", "img2.tif"]);

% Writable
bim = blockedImage(dest, imageSize, blockSize, initVal, Mode="w");

% Override adapter
bim = blockedImage(source, Adapter=myAdapter);

BlockSize

• Defaults to a factor of IOBlockSize (adapter's native read unit), targeting ~1024 or the minimum dimension size.
• The object's BlockSize is the default for all downstream processing and visualization.
• Override per-operation via the BlockSize name-value argument (e.g., when a DL network requires a fixed input size).
• For best performance, keep BlockSize as a multiple of IOBlockSize and as large as practical.

Reading Data

• getBlock(bim, blocksub) — read one block by block subscript. Returns [data, blockinfo].
• getRegion(bim, pixelStart, pixelEnd) — read an arbitrary pixel region.
• gather(bim) — load entire image into memory (only if it fits in RAM). Use Level to select resolution.

Processing with apply()

apply(bim, @fcn) calls fcn on each block one at a time. Scales to arbitrarily large data. Supports UseParallel=true for cluster execution.

Key name-value arguments: BlockSize, BorderSize, Level, ExtraImages, PadPartialBlocks, PadMethod, OutputLocation, BatchSize, Resume, BlockLocationSet, DisplayWaitbar.

Block struct (bs) fields

The user function receives a block struct bs with:

• bs.Data — pixel data for this block (includes border if BorderSize > 0)
• bs.Start — 1-by-N pixel subscript of the block's top-left corner (excluding border)
• bs.End — 1-by-N pixel subscript of the block's bottom-right corner (excluding border)
• bs.BorderSize — 1-by-N border size used
• bs.Level — resolution level
• bs.BlockSub — block subscript (index into SizeInBlocks)
• bs.ImageNumber — index of the source image (when processing an array)

Struct / non-image output

When the output per block is not image data (e.g., detection results, counts), use images.blocked.MATBlocks adapter to store structs:

bResults = apply(bim, @countNuclei, ...
    BlockLocationSet=bls, BorderSize=[r r], ...
    Adapter=images.blocked.MATBlocks, OutputLocation="results_dir");
allResults = gather(bResults);  % load all structs into memory

Limit work to real data

• Crop first. crop(bim, cstart, cend) is virtual — no data copy, executes instantly. Use to restrict to a spatial ROI or single channel before processing.
• Use masks with selectBlockLocations. Build a mask from a coarse level, then process only blocks overlapping the ROI. See references/specialized-workflows.md for the full pattern including InclusionThreshold and block-size tuning.

Key patterns

• Single pass. Perform all pre/post-processing inside fcn to go through data only once.
• World coordinates. fcn receives block metadata. Use sub2world() to convert block-local positions to world coordinates before returning detection results.
• BorderSize for overlap. Set BorderSize when the processing function needs context from neighboring blocks (e.g., DL network with wider input than output field of view).
• Unique detections via centroids. Set BorderSize >= half the max object diameter. In the user function, discard detections whose centroid falls in the border region, then convert remaining centroids to global coordinates using bs.Start:

  onBorder = centroids(:,1) < bs.BorderSize(1) ...
           | centroids(:,1) > size(bs.Data,1) - 2*bs.BorderSize(1) ...
           | centroids(:,2) < bs.BorderSize(2) ...
           | centroids(:,2) > size(bs.Data,2) - 2*bs.BorderSize(2);
  centroids(onBorder,:) = [];
  centroids = centroids + bs.Start(1:2);  % local -> global
  

  For inconsistent detectors, add a second apply() pass with BorderSize=BlockSize to deduplicate by bounding-box overlap in global coordinates.
• Multi-resolution output. Use makeMultiLevel2D(bResult, Scales=[1 0.5 0.3 0.1]) to generate a pyramid from a single-level result. Use "nearest" interpolation for label/mask data. Or return multiple levels from fcn and combine with concatenateLevels.
• Small block output. When per-block output is tiny (e.g., a count), use gather() to load all results into memory, or write() to rewrite with a larger block size.
• Resume. Set Resume=true to continue a previously interrupted apply() from where it left off. For in-place resume with overwrite, see references/specialized-workflows.md.
• Multi-output apply. When the processing function returns N outputs, apply returns N separate blockedImages. Use this to split channels in a single pass:

  [o1, o2, o3, o4] = apply(bim, @(bs) splitChannels(bs), ...
      BlockSize=[10 2048 2048 4], UseParallel=true, OutputLocation=outDir);
  

• ExtraImages for pixel-level masking. Pass a second blockedImage via ExtraImages to access it in the block function alongside the primary image. See references/specialized-workflows.md for the full pattern including resolution mismatch handling.
• Pyramidal TIFF in one pass. Generate multi-resolution levels and write as a single pyramidal TIFF using LevelImages:

  [b1,b2,b3,b4] = apply(bim, @resizeBlocks);
  write(bim, "pyramid.tif", LevelImages=[b1 b2 b3 b4], BlockSize=[2048 2048]);
  function [b1,b2,b3,b4] = resizeBlocks(bs)
      b1 = imresize(bs.Data, 0.5, 'nearest');
      b2 = imresize(b1, 0.5, 'nearest');
      b3 = imresize(b2, 0.5, 'nearest');
      b4 = imresize(b3, 0.5, 'nearest');
  end
  

• BatchSize. Process multiple blocks at once by setting BatchSize=N. The block function receives bs.Data with an extra trailing dimension of size N and bs.BatchSize indicating the count. Useful when per-block overhead is high.

Iterative development workflow

Processing full-resolution blocked images is slow. Shorten the feedback cycle by prototyping on smaller data first:

1. Prototype on a coarse level. Use gather(bim, Level=bim.NumLevels) to load the coarsest level into memory. Develop and tune parameters on this small image, then apply to the full image:

   imCoarse = gather(bim, Level=bim.NumLevels);
   thresh = graythresh(im2gray(imCoarse));
   bResult = apply(bim, @(bs) ~imbinarize(im2gray(bs.Data), thresh));
   

2. Prototype on an ROI. Extract a small region with representative features using getRegion, tune parameters, then apply:

   roi = getRegion(bim, [900 2400], [1700 3300], Level=1);
   thresh = graythresh(im2gray(roi));
   bResult = apply(bim, @(bs) ~imbinarize(im2gray(bs.Data), thresh));
   

Parallel processing

When writing results to disk, check for Parallel Computing Toolbox and use it if available:

outDir = tempname;
usepar = ~isempty(ver("parallel"));
bResult = apply(bim, @fcn, UseParallel=usepar, OutputLocation=outDir);

Always include UseParallel when specifying OutputLocation — set it based on toolbox availability.

Visualization

Choose 2D vs 3D based on the data, not convenience:

• 3D volumes (Z-stacks, CT, MRI, microscopy stacks): Default to volshow + viewer3d. Use OverlayData for multi-channel overlay, SlicePlanes for interactive browsing, and linkviewers for overview+detail. Read references/3d-volumes.md for full patterns. Do NOT flatten 3D data into 2D slice-by-slice viewers — that discards spatial context and ignores the 3D tooling.
• 2D images (whole-slide, satellite, mosaics): Use imageshow + viewer2d. Read references/visualization.md for overlays, synced views, overview+detail, and annotations.

If the data has a Z/depth/slice dimension, it is 3D — use the 3D path even if each slice is a separate file.

Format Conversion

Use write(bim, destination) to convert between formats or re-chunk data:

write(bim, "output.tif", BlockSize=[512 512 3]);

For adapter selection guidance when converting (JPEG vs PNG vs H5Blocks), see "Adapter Preference for Format Conversion" in references/data-discovery.md.

---

Copyright 2026 The MathWorks, Inc.