InsightSoftwareConsortium / InsightSoftwareConsortium/ITK

NrrdImageIO does not stream: a small region costs a whole volume decode

Open
#6,805 4 comments 0 reactions 0 assignees View on GitHub
Dominant language
C++
Stars
1.7k
Forks
748
Avg merge
1d 1h
Merged PRs (30d)
64

Description

`NrrdImageIO` does not implement streamed reads, so extracting a small region of a NRRD file costs a full
decode of the volume. For uncompressed NRRD (`encoding: raw`) the pixel data is a flat block, exactly like
MetaImage, so the region could be served by seeking.

256^3 uint16 volume (32 MiB), region extracted with `SetExtractIndex` / `SetExtractSize`, minimum of 15
runs, SimpleITK 2.5.5, Linux, local NVMe, page cache warm:

| format | 16^3 | 64^3 | 128^3 | whole 256^3 |
| --- | --- | --- | --- | --- |
| .nii | 0.75 ms | 1.81 ms | 6.10 ms | 41.34 ms |
| .mha | 0.52 ms | 2.28 ms | 8.46 ms | 21.87 ms |
| .nrrd | 24.71 ms | 25.17 ms | 32.91 ms | 26.69 ms |

NIfTI and MetaImage scale with the size of the region. NRRD does not: it stays around 25 to 33 ms whatever
is asked of it. A 16^3 region is 0.02% of the volume and costs 93% of what reading the whole volume costs,
and 33x what NIfTI charges for the same region.

The file is genuinely uncompressed. The header says `encoding: raw` and it is 33554740 bytes for 33554432
bytes of pixel data.

The mechanism is in `Read()`. `NiftiImageIO::Read()` takes `this->GetIORegion()` and reads that region.
`MetaImageIO::Read()` compares `m_IORegion` against the largest region and reads the sub-block.
`NrrdImageIO::Read()` refers to `GetIORegion` or `m_IORegion` nowhere at all, 0 occurrences against 1 in
`itkNiftiImageIO.cxx` and 13 in `itkMetaImageIO.cxx`, and calls `nrrdLoad(nrrd, this->GetFileName(),
nullptr)` on the whole file. Consistently with that, `itkNrrdImageIO.h` overrides neither
`CanStreamRead()` nor `GenerateStreamableReadRegionFromRequestedRegion()`, where NIfTI overrides the
second and MetaImage both.

```python
import numpy as np, SimpleITK as sitk, time, tempfile, pathlib

tmp = pathlib.Path(tempfile.mkdtemp()); N = 256
arr = np.random.randint(0, 4000, (N, N, N), dtype=np.uint16)
for ext in ("nii", "mha", "nrrd"):
p = str(tmp / f"v.{ext}")
w = sitk.ImageFileWriter(); w.SetFileName(p); w.SetUseCompression(False)
w.Execute(sitk.GetImageFromArray(arr))
for R in (16, 64, 128, 256):
def read():
r = sitk.ImageFileReader(); r.SetFileName(p)
r.SetExtractIndex([0, 0, 0]); r.SetExtractSize([R, R, R])
return sitk.GetArrayViewFromImage(r.Execute())
read()
ts = []
for _ in range(15):
t = time.perf_counter()
read()
ts.append(time.perf_counter() - t)
print(f"{ext:5} {R:>4}^3 {min(ts) * 1e3:8.2f} ms")
```

Where this matters: patch based deep learning reads a volume as many small overlapping windows, one per
training sample or per inference patch. NRRD is 3D Slicer's native format, so a dataset exported from
Slicer takes the whole volume path on every patch. In KonfAI we detect the format and read the volume
whole rather than pay this per patch, which is a workaround for something that looks fixable.

Happy to work on streamed reads for the raw encoding if that would be welcome, and glad to hear if there
is a reason it was left out that I am not seeing.

Contributor guide

Open the contributing guide

Research direction

Start by reading itkNrrdImageIO.cxx and itkNrrdImageIO.h, then compare Read(), CanStreamRead(), and GenerateStreamableReadRegionFromRequestedRegion() with the corresponding NIfTI and MetaImage implementations. Use the supplied Python benchmark to verify that raw NRRD extraction reads the requested region rather than decoding the whole volume, while full-volume reads remain correct.

Written by the indexing model from the issue text.

Assessment

Tech stack
cpp
Domain
computer-vision
Issue type
Bug
Difficulty
4/5
Estimated time
3-5 days
Activity status
Active
Clarity
Clearly specified
Newbie friendliness
68/100

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.