Skip to content

rayolddog/calibrated-ct-encoding

Repository files navigation

calibrated-ct-encoding

Calibrated CT Preprocessing: Fixed Hounsfield Unit Windows for Deep Learning

A principled approach to encoding CT DICOM pixel data as float16 tensors for AI/deep learning, preserving the physical calibration of the Hounsfield scale. Includes a pure Python implementation and a hybrid Python + Zig implementation for high-throughput cache building.


The Problem with Per-Image Normalization

The standard image preprocessing idiom — normalizing each image by its own minimum and maximum pixel value — is destructive when applied to CT data:

# WRONG for CT — destroys Hounsfield calibration
tensor = (image - image.min()) / (image.max() - image.min())

CT pixel values are not arbitrary intensity values. They are Hounsfield Units (HU), a calibrated physical scale defined so that:

  • Water = 0 HU (by definition)
  • Air = −1000 HU (by definition)

Every CT scanner in clinical use is calibrated to this scale. A pixel value of 65 HU means the same thing — the density of acute blood — regardless of which scanner, which patient, or which slice produced it.

Per-image normalization maps that value to a different float every time, depending on what else is in the image. A bright slice (large hemorrhage) and a dim slice (thin subdural) will encode the same 65 HU blood pixel at completely different tensor values. The physical information is discarded.


The Solution: Fixed HU Window Encoding

Map a fixed, globally-defined HU range linearly to [0.0, 1.0]:

tensor = clip( (HU − hu_low) / (hu_high − hu_low),  0.0,  1.0 )

Values below hu_low0.0. Values above hu_high1.0. Values in between → linear. The same HU value always produces the same tensor value, across all slices and all patients.


Hounsfield Unit Reference

Tissue HU range
Air (external) −1000
Fat −100 to −50
Water / CSF 0 to 15
White matter 20 to 30
Gray matter 30 to 42
Acute ICH / blood 50 to 80
Hyperacute clot 80 to 100
Cortical bone 300 to 1000
Metal / hardware > 1000

Three Window Presets

Defined in hu_windows.py:

Wide window — WINDOW_WIDE

HU range : −1024 to 3071  (4095 HU span)

Full practical CT range. Covers air, fat, water, soft tissue, bone, and metal without saturation. The contrast between similar soft tissues is compressed to a small fraction of [0, 1].

Gray matter → blood delta: ~0.007

Use for: multi-tissue classification, bone and implant tasks.

Medium window — WINDOW_MEDIUM

HU range : −200 to 200  (400 HU span)

Soft-tissue and early-bone range. Excludes bulk air and dense bone. Captures fat, water, CSF, brain parenchyma, blood, and early cortical bone.

Gray matter → blood delta: ~0.068

Use for: brain anatomy, ICH in anatomical context.

Narrow window — WINDOW_NARROW

HU range : 48 to 90  (42 HU span)

Acute ICH detection range. Gray matter, CSF, white matter, fat, and air all map to 0.0. Only the acute blood density band occupies [0, 1].

Gray matter → blood delta: ~0.316 (5.5× wider than per-image)

Use for: maximum ICH vs. normal brain pixel contrast.


Experimental Validation

The figure below compares tensor value distributions for gray matter (HU 28–42) and acute blood (HU 50–80) pixels across 1,000 ICH-positive CT slices from the RSNA Intracranial Hemorrhage Detection dataset (n = 28.3M gray matter pixels, 5.6M blood pixels).

Normalization comparison

Scheme Gray matter (mean ± std) Acute blood (mean ± std) Δ mean GM std ratio
Per-image (min/max) 0.695 ± 0.022 0.752 ± 0.022 0.057 1.0× (baseline)
Fixed medium (−200 to 200 HU) 0.586 ± 0.009 0.653 ± 0.020 0.068 2.4× tighter
Fixed narrow (48 to 90 HU) 0.000 ± 0.000 0.316 ± 0.190 0.316 ∞ (constant)

Key observations:

  • Δ mean: the fixed narrow window achieves 5.5× greater separation between gray matter and acute blood than per-image normalization.
  • Gray matter std: per-image normalization maps the same gray matter tissue to different tensor values depending on the slice content (std 0.022). The fixed medium window reduces this variance by 2.4× (std 0.009), reflecting the physical reality that gray matter density is consistent across patients.
  • Fixed narrow: gray matter is entirely below the window floor and maps exactly to 0.0 on every slice. The blood distribution's spread (std 0.190) captures genuine biological variance in clot density and age.

Experiment code: compare_normalization.py


ICH Classification Results

To validate that the calibrated encoding provides actionable signal for deep learning, a SE-ResNeXt50 classifier was trained from random weights (no ImageNet pretraining) on the narrow-window cache (HU 48–90 → [0.0, 1.0]).

Training configuration:

  • Architecture: SE-ResNeXt50 (25.5M parameters), single-channel input
  • Loss: Focal loss (γ=2.0) with label smoothing (0.05)
  • Optimizer: AdamW, LR 3×10⁻⁴ with cosine decay, batch size 32
  • Dataset: RSNA Intracranial Hemorrhage Detection (Kaggle 2019)
    • Train+val: 544,685 slices | Test: 200,000 slices (held out)
  • Cache encoding: WINDOW_NARROW (48–90 HU), corrected fixed-window

Final test set results (200,000 held-out slices, best checkpoint epoch 25):

Class Prevalence AUC Sensitivity Specificity LR+ LR−
Epidural 0.42% 0.966 0.891 0.924 11.68 0.118
Intraparenchymal 4.85% 0.963 0.881 0.917 10.61 0.130
Intraventricular 3.52% 0.976 0.914 0.931 13.19 0.092
Subarachnoid 4.79% 0.936 0.875 0.847 5.71 0.147
Subdural 6.33% 0.957 0.878 0.901 8.90 0.135
Any ICH 14.49% 0.959 0.880 0.901 8.91 0.133
Mean 0.960

Thresholds were chosen by Youden J (maximising sensitivity + specificity) on the test set. LR+ and LR− are likelihood ratios for positive and negative results respectively — directly useful for Bayesian clinical reasoning.

Significance of random-weight initialisation: A prior run using the same architecture with ImageNet-pretrained weights but an incorrectly encoded cache (wrong HU recovery formula, 0.05–0.95 window margins) reached a best validation AUC of only 0.897 at epoch 29. The corrected encoding with random weights surpassed that result at epoch 2 (AUC 0.903) and reached 0.960 on the held-out test set. This demonstrates that the calibrated HU encoding, not the pretrained weights, is the primary driver of model performance.


Repository Structure

calibrated-ct-encoding/
│
├── hu_windows.py              # Core skill: HUWindow dataclass, 3 presets,
│                              # apply_window() — pure Python, no dependencies
│
├── dicom_reader_1ch.py        # Robust DICOM → HU reader (pydicom)
│                              # Handles uint16/int16/float, big-endian,
│                              # RescaleSlope/Intercept, HU validation rules
│
├── build_1ch_cache.py         # Wide-window DICOM → .npz cache builder
├── build_medium_cache.py      # Medium-window DICOM → .npz cache builder
│                              # (pure Python; uses hu_windows.py directly)
│
├── build_cache_zig.py         # High-throughput hybrid cache builder
│                              # Python handles DICOM I/O; Zig handles windowing
│                              # Supports --window wide|medium|narrow|all
│                              # 'all' builds three caches in one DICOM pass
│
├── hu_tensor/                 # Zig shared library
│   ├── build.zig              # Zig 0.15.x build file
│   └── src/
│       └── hu_tensor.zig      # apply_window(), apply_three_windows()
│                              # Exported as C ABI; called via Python ctypes
│
├── compare_normalization.py   # Experiment: per-image vs. fixed-window
│                              # Generates normalization_comparison.png/.json
│
├── normalization_comparison.png   # Experiment figure (1,000 ICH-positive slices)
└── normalization_comparison.json  # Experiment statistics (machine-readable)

Installation

Requirements

pip install pydicom numpy pandas scipy matplotlib

Python 3.10+ recommended. The pure Python implementation has no other dependencies. The hybrid Zig implementation requires Zig 0.15.x (see below).

Zig (for hybrid implementation, AMD64 Ubuntu)

Install Zig 0.15.x via snap:

sudo snap install zig --classic --channel=0.15.x/stable

Or download a binary from ziglang.org/download and place zig on your PATH.

Build the shared library:

cd hu_tensor
zig build -Doptimize=ReleaseFast

Output: hu_tensor/zig-out/lib/libhu_tensor.so

Verify:

from build_cache_zig import _load_lib
lib = _load_lib('hu_tensor/zig-out/lib/libhu_tensor.so')
print(lib.hu_tensor_version())   # → 10000  (version 1.0.0)

Usage

Pure Python — single window

import numpy as np
from hu_windows import apply_window, WINDOW_NARROW, WINDOW_MEDIUM, WINDOW_WIDE

# hu is a float32 ndarray of Hounsfield values (from dicom_reader_1ch)
tensor = apply_window(hu, WINDOW_NARROW)   # → float16, shape (512, 512)
np.savez_compressed('slice.npz', image_norm=tensor)

Build a cache — pure Python

# Medium window (−200 to 200 HU)
python build_medium_cache.py \
    --dcm-dir  /path/to/dicom/root \
    --train-dir /path/to/cache_medium_train \
    --test-dir  /path/to/cache_medium_test \
    --splits-file ./checkpoints_1ch/data_splits.json \
    --workers 8

# Wide window (−300 to 180 HU, brain/blood window)
python build_1ch_cache.py \
    --dcm-dir    /path/to/dicom/root \
    --output-dir /path/to/cache_wide \
    --workers 8

Build a cache — hybrid Python + Zig (faster)

# Single window
python build_cache_zig.py \
    --dcm-dir /path/to/dicom/root \
    --window medium \
    --workers 8

# All three windows in one DICOM pass (most efficient)
python build_cache_zig.py \
    --dcm-dir /path/to/dicom/root \
    --window all \
    --workers 8

The --window all mode calls apply_three_windows() in the Zig library, reading each DICOM pixel array once and writing wide, medium, and narrow tensors simultaneously — one third the memory bandwidth of three separate runs.

Run the normalization experiment

python compare_normalization.py \
    --cache-dir /path/to/cache_1ch \
    --labels    /path/to/stage_2_train.csv \
    --n 1000 --workers 8

DICOM to HU Conversion

The conversion from raw DICOM pixel values to Hounsfield Units uses the metadata stored in the DICOM header, as defined by the DICOM standard:

HU = pixel_value × RescaleSlope + RescaleIntercept

RescaleSlope is 1.0 for virtually all CT scanners. RescaleIntercept is typically −1024 for modern CT scanners.

dicom_reader_1ch.py applies four validation rules to exclude invalid files:

  1. File must be readable by pydicom
  2. pixel_array must be decodable and 2-dimensional
  3. BitsAllocated must be in {8, 16, 32, 64}
  4. Resulting HU range must be plausible (rejects uncalibrated integer files, scout/localizer images, and files where the intercept was not applied)

Files lacking RescaleIntercept fall through to a scanner-convention fallback (uint16 − 1024) and are subject to rule 4 validation.


Zig Library API

The Zig shared library exports two functions via C ABI, callable from any language with a C FFI (Python ctypes, Rust, C, etc.):

// Single window: hu[n] → output[n] (float16 stored as uint16 bit patterns)
void apply_window(
    const float* hu,
    size_t n,
    float hu_low,
    float hu_high,
    uint16_t* output
);

// Three windows in one pass — 1/3 the memory bandwidth of three calls
void apply_three_windows(
    const float* hu, size_t n,
    float lo0, float hi0,    // window 0 (e.g. wide)
    float lo1, float hi1,    // window 1 (e.g. medium)
    float lo2, float hi2,    // window 2 (e.g. narrow)
    uint16_t* out0,
    uint16_t* out1,
    uint16_t* out2
);

// Returns library version as major*10000 + minor*100 + patch
uint32_t hu_tensor_version(void);

Output arrays contain float16 values stored as uint16 bit patterns. In Python: np.frombuffer(out, dtype=np.uint16).view(np.float16).


Data

Experiments use the RSNA Intracranial Hemorrhage Detection dataset (Kaggle, 2019). The dataset is not included in this repository. Labels file: stage_2_train.csv.


Credits

John B. Bramble, MD Concept, clinical domain expertise, and research direction. The principle of encoding CT pixels using fixed Hounsfield Unit windows rather than per-image normalization was developed and validated through iterative experimentation with multiple neural network architectures and window widths. This work continues a personal research interest in computer-aided diagnosis that began in the late 1980s, when Dr. Bramble developed an early CAD system for arthritis classification in C on a Commodore Amiga, working under Samuel J. Dwyer III (a principal architect of the DICOM standard) and Lawrence L. Cook at the University of Kansas Medical Center.

Claude (claude-sonnet-4-6, Anthropic) Implementation of all source files in this repository, including the pure Python HU windowing skill (hu_windows.py), the DICOM reader and cache builders, the Zig shared library, the experiment design and analysis script, and this documentation. Development was conducted interactively via Claude Code in a single-session pair-programming workflow.


License

MIT License. See LICENSE for details.


References

  1. Hounsfield, G.N. (1973). Computerized transverse axial scanning (tomography). British Journal of Radiology, 46(552), 1016–1022.

  2. RSNA Intracranial Hemorrhage Detection Challenge (2019). Radiological Society of North America. https://www.kaggle.com/c/rsna-intracranial-hemorrhage-detection

  3. DICOM Standard, PS 3.3 — Information Object Definitions. National Electrical Manufacturers Association (NEMA). https://www.dicomstandard.org/

  4. Flanders, A.E. et al. (2020). Construction of a Machine Learning Dataset through Collaboration: The RSNA 2019 Brain CT Hemorrhage Challenge. Radiology: Artificial Intelligence, 2(3).

About

Calibrated CT Preprocessing: Fixed Hounsfield Unit Windows for Deep Learning

Resources

License

Stars

Watchers

Forks

Releases

No releases published

Packages

 
 
 

Contributors