diff --git a/data/DirLab-4DCT/README.md b/data/DirLab-4DCT/README.md index 686b65d..a3f9744 100644 --- a/data/DirLab-4DCT/README.md +++ b/data/DirLab-4DCT/README.md @@ -8,10 +8,16 @@ DIR-Lab distributes each case individually and may require registration. 1. Visit the DIR-Lab 4D-CT page and request/download the case archives: https://med.emory.edu/departments/radiation-oncology/research-laboratories/deformable-image-registration/downloads-and-reference-data/4dct.html -2. Extract each case's raw images into `data/DirLab-4DCT/` (see layout - below). The `.mhd` headers already committed in this directory point at - those raw files, so no conversion step is required — see - "About the Committed `.mhd` Files" below. +2. Extract each case's raw images into `data/DirLab-4DCT/downloaded_data/` + (see layout below). The `.mhd` headers already committed in that + directory point at those raw files, so no header-writing step is + required — see "About the Committed `.mhd` Files" below. +3. Run `data/DirLab-4DCT/fix_downloaded_data.py`. DIR-Lab's raw volumes are + not in Hounsfield units (see "Fixing Raw Intensities" below); this + script reads every `.mhd` header in `downloaded_data/` with its backing + pixel data present, corrects the intensities, and writes one compressed + `.mha` volume per phase directly into `data/DirLab-4DCT/`. Tutorials and + experiments read these top-level `.mha` files, not the raw `.mhd` ones. Once populated, check the layout with: @@ -21,26 +27,31 @@ from physiotwin4d import DataDownloadTools assert DataDownloadTools.VerifyDirLab4DCTData("data/DirLab-4DCT") ``` -**Directory structure after download:** +**Directory structure after download and fixing:** ``` data/DirLab-4DCT/ -├── Case1Pack/ -│ ├── Images/ # T00-T90 phase images -│ ├── ExtremePhases/ # T00 and T50 (max inhale/exhale) -│ └── Sampled4D/ # Sampled time points -├── Case1Pack_T00.mhd # Already-committed headers, per phase -├── Case1Pack_T10.mhd +├── downloaded_data/ +│ ├── Case1Pack/ +│ │ ├── Images/ # T00-T90 phase images +│ │ ├── ExtremePhases/ # T00 and T50 (max inhale/exhale) +│ │ └── Sampled4D/ # Sampled time points +│ ├── Case1Pack_T00.mhd # Already-committed headers, per phase +│ ├── Case1Pack_T10.mhd +│ ... +│ ├── Case10Pack/ +│ └── Convert4DCTToMHD.py # Documents how the .mhd headers were generated +├── Case1Pack_T00.mha # Written by fix_downloaded_data.py +├── Case1Pack_T10.mha ... -├── Case10Pack/ -├── Convert4DCTToMHD.py # Documents how the .mhd headers were generated +├── fix_downloaded_data.py └── README.md (this file) ``` ### About the Committed `.mhd` Files -The `Case*Pack_T*.mhd` files in this directory are already committed to the -repository, but they are only MetaImage **headers** (a few hundred bytes -each) — for example: +The `Case*Pack_T*.mhd` files in `downloaded_data/` are already committed to +the repository, but they are only MetaImage **headers** (a few hundred +bytes each) — for example: ``` ObjectType = Image @@ -60,6 +71,16 @@ until you complete the manual download above and the referenced raw DIR-Lab archives. It is included for provenance/documentation only — you do not need to run it; the `.mhd` files are already committed. +### Fixing Raw Intensities + +DIR-Lab's raw `.mhd`/`.img` volumes store intensities offset by +1024 from +Hounsfield units, not real HU. `fix_downloaded_data.py` calls +`DataDownloadTools.FixDirLab4DCTData`, which subtracts 1024 and clips the +result to `[-1024, 1024]` for every downloaded phase, writing one +compressed `.mha` file per phase into `data/DirLab-4DCT/` (`.gitignore`d, +like the raw data). Run it once after downloading; tutorials and +experiments assume the `.mha` files it produces already exist. + ## Overview Benchmark dataset for 4D CT respiratory motion analysis. Contains 10 cases @@ -104,9 +125,16 @@ If you use this dataset, please cite: ### Files in This Directory -- `Convert4DCTToMHD.py` — documents how the committed `.mhd` headers were - generated from the raw DIR-Lab archives; not needed to run -- `Case*Pack_T*.mhd` — MetaImage headers for each case/phase (see above) +- `fix_downloaded_data.py` — converts `downloaded_data/*.mhd` to + HU-corrected `Case*Pack_T*.mha` files in this directory (see "Fixing Raw + Intensities" above) +- `downloaded_data/Convert4DCTToMHD.py` — documents how the committed + `.mhd` headers were generated from the raw DIR-Lab archives; not needed + to run +- `downloaded_data/Case*Pack_T*.mhd` — MetaImage headers for each + case/phase (see above) +- `Case*Pack_T*.mha` — HU-corrected volumes written by + `fix_downloaded_data.py`; what tutorials and experiments actually read ### Additional Resources diff --git a/data/DirLab-4DCT/Case10Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case10Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case10Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case10Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case10Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case1Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case1Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case1Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case1Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case2Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case2Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case2Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case2Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case3Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case3Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case3Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case3Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case4Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case4Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case4Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case4Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case5Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case5Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case5Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case5Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case6Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case6Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case6Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case6Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case7Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case7Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case7Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case7Pack_T90.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T00.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T10.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T20.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T30.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T40.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T50.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T60.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T70.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T80.mhd diff --git a/data/DirLab-4DCT/Case8Deploy_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case8Deploy_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case8Deploy_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case8Deploy_T90.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T00.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T00.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T00.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T00.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T10.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T10.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T10.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T10.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T20.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T20.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T20.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T20.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T30.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T30.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T30.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T30.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T40.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T40.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T40.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T40.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T50.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T50.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T50.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T50.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T60.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T60.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T60.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T60.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T70.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T70.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T70.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T70.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T80.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T80.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T80.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T80.mhd diff --git a/data/DirLab-4DCT/Case9Pack_T90.mhd b/data/DirLab-4DCT/downloaded_data/Case9Pack_T90.mhd similarity index 100% rename from data/DirLab-4DCT/Case9Pack_T90.mhd rename to data/DirLab-4DCT/downloaded_data/Case9Pack_T90.mhd diff --git a/data/DirLab-4DCT/Convert4DCTToMHD.py b/data/DirLab-4DCT/downloaded_data/Convert4DCTToMHD.py similarity index 100% rename from data/DirLab-4DCT/Convert4DCTToMHD.py rename to data/DirLab-4DCT/downloaded_data/Convert4DCTToMHD.py diff --git a/data/DirLab-4DCT/fix_downloaded_data.py b/data/DirLab-4DCT/fix_downloaded_data.py new file mode 100644 index 0000000..e5f1d29 --- /dev/null +++ b/data/DirLab-4DCT/fix_downloaded_data.py @@ -0,0 +1,14 @@ +# %% +from pathlib import Path + +from physiotwin4d.data_download_tools import DataDownloadTools + +_HERE = Path(__file__).resolve().parent + +# %% +input_dir = _HERE / "downloaded_data" +output_dir = _HERE + +# Converts raw DirLab-4DCT intensities (.mhd) to clipped Hounsfield units and +# writes each result as a compressed .mha file in output_dir. +DataDownloadTools.FixDirLab4DCTData(input_dir, output_dir) diff --git a/docs/api/model_registration/index.rst b/docs/api/model_registration/index.rst index 456aed3..d612dea 100644 --- a/docs/api/model_registration/index.rst +++ b/docs/api/model_registration/index.rst @@ -42,7 +42,7 @@ Most users should access model registration through the workflow: patient_image=itk.imread("patient_ct.nii.gz"), ) - result = workflow.run_workflow() + result = workflow.process() See Also ======== diff --git a/docs/api/workflows.rst b/docs/api/workflows.rst index 4798d43..41502c5 100644 --- a/docs/api/workflows.rst +++ b/docs/api/workflows.rst @@ -147,7 +147,7 @@ Statistical Shape Modeling patient_image=itk.imread("patient_ct.nii.gz"), ) - result = workflow.run_workflow() + result = workflow.process() High-Resolution 4D CT Reconstruction ==================================== @@ -172,7 +172,7 @@ High-Resolution 4D CT Reconstruction ) workflow.set_upsample_to_fixed_resolution(True) - result = workflow.run_workflow() + result = workflow.process() See Also ======== diff --git a/docs/developer/registration_models.rst b/docs/developer/registration_models.rst index f7cada6..f18b535 100644 --- a/docs/developer/registration_models.rst +++ b/docs/developer/registration_models.rst @@ -22,7 +22,7 @@ Recommended Entry Point patient_image=itk.imread("patient_ct.nii.gz"), ) - result = workflow.run_workflow() + result = workflow.process() Lower-Level Classes =================== diff --git a/docs/developer/segmentation.rst b/docs/developer/segmentation.rst index 9f05580..3fb758e 100644 --- a/docs/developer/segmentation.rst +++ b/docs/developer/segmentation.rst @@ -72,8 +72,9 @@ produces by populating ``self.taxonomy``. The base class owns the ... The base class then provides the shared pre/post-processing, contrast-agent -fusion, and ``create_anatomy_group_masks()`` that walks the taxonomy and -emits one mask per registered group plus ``"other"``. +fusion, and ``create_anatomy_group_labelmaps()`` that walks the taxonomy and +emits one labelmap (preserving the original label ids) per registered group +plus ``"other"``. Custom Anatomy Looks ==================== diff --git a/docs/tutorials.rst b/docs/tutorials.rst index dae135e..a617f2e 100644 --- a/docs/tutorials.rst +++ b/docs/tutorials.rst @@ -315,7 +315,7 @@ Inner API usage ) if pca_model is not None: workflow.set_use_pca_registration(True, pca_model=pca_model) - result = workflow.run_workflow() + result = workflow.process() registered_surface = result["registered_template_model_surface"] Run @@ -425,7 +425,7 @@ Inner API usage registration_method=registration_method, ) workflow.set_modality("ct") - result = workflow.run_workflow() + result = workflow.process() reconstructed_images = result["reconstructed_images"] Run @@ -490,9 +490,9 @@ Inner API usage labelmap_interior_object_ids=LABELMAP_INTERIOR_OBJECT_IDS, ) ssm_fit_workflow.set_use_pca_registration( - use_pca_registration=True, pca_model=ssm_model, pca_uses_surface=False, + use_pca_registration=True, pca_model=ssm_model, use_surface=False, ) - ssm_fit_workflow_result = ssm_fit_workflow.run_workflow() + ssm_fit_workflow_result = ssm_fit_workflow.process() ssm_pca_coefficients = ssm_fit_workflow.pca_coefficients icon_registration_method = RegisterImagesICON() @@ -502,7 +502,7 @@ Inner API usage fixed_image=ref_image, registration_method=icon_registration_method, ) - reg_result = reg_workflow.run_workflow() + reg_result = reg_workflow.process() reconstructed_images = reg_result["reconstructed_images"] Run diff --git a/experiments/Heart-Create_Statistical_Model/README.md b/experiments/Heart-Create_Statistical_Model/README.md index 1ea3a03..92727e9 100644 --- a/experiments/Heart-Create_Statistical_Model/README.md +++ b/experiments/Heart-Create_Statistical_Model/README.md @@ -161,7 +161,7 @@ workflow = WorkflowFitStatisticalModelToPatient( pca_number_of_modes=10 ) -registered_mesh = workflow.run_workflow() +registered_mesh = workflow.process()["registered_template_model_surface"] ``` ## Requirements diff --git a/experiments/Heart-GatedCT_To_USD/2-generate_segmentation.py b/experiments/Heart-GatedCT_To_USD/2-generate_segmentation.py index f8859c1..7482466 100644 --- a/experiments/Heart-GatedCT_To_USD/2-generate_segmentation.py +++ b/experiments/Heart-GatedCT_To_USD/2-generate_segmentation.py @@ -81,7 +81,7 @@ ) else: labelmap_image = itk.imread(os.path.join(output_dir, f"{outname}.all_mask.mha")) - result = seg.create_anatomy_group_masks(labelmap_image) + result = seg.create_anatomy_group_labelmaps(labelmap_image) lung_mask = result["lung"] heart_mask = result["heart"] diff --git a/experiments/Heart-Statistical_Model_To_Patient/heart_model_to_patient-CHOPValve.py b/experiments/Heart-Statistical_Model_To_Patient/heart_model_to_patient-CHOPValve.py index 59ff4f8..cacf30e 100644 --- a/experiments/Heart-Statistical_Model_To_Patient/heart_model_to_patient-CHOPValve.py +++ b/experiments/Heart-Statistical_Model_To_Patient/heart_model_to_patient-CHOPValve.py @@ -67,7 +67,7 @@ ) # %% -results = registrar.run_workflow() +results = registrar.process() # %% registered_model = results["registered_template_model"] diff --git a/experiments/Lung-GatedCT_To_USD/0-register_dirlab_4dct.py b/experiments/Lung-GatedCT_To_USD/0-register_dirlab_4dct.py index b509c09..7621194 100644 --- a/experiments/Lung-GatedCT_To_USD/0-register_dirlab_4dct.py +++ b/experiments/Lung-GatedCT_To_USD/0-register_dirlab_4dct.py @@ -98,11 +98,12 @@ def register_image( os.makedirs(output_dir, exist_ok=True) for case_name in case_names: + # .mha files are DirLab-4DCT data already converted to HU by + # data/DirLab-4DCT/fix_downloaded_data.py. fixed_image_filename = ( - f"{input_dir}/{case_name}_T{fixed_image_num * 10:02d}.mhd" + f"{input_dir}/{case_name}_T{fixed_image_num * 10:02d}.mha" ) fixed_image = itk.imread(fixed_image_filename) - fixed_image = DataDirLab4DCT().fix_image(fixed_image) print("Segmenting fixed image...") fixed_result = seg_image.segment(fixed_image) @@ -151,9 +152,8 @@ def register_image( for image_num in images: if image_num != fixed_image_num: moving_image = itk.imread( - os.path.join(input_dir, f"{case_name}_T{image_num * 10:02d}.mhd") + os.path.join(input_dir, f"{case_name}_T{image_num * 10:02d}.mha") ) - moving_image = DataDirLab4DCT().fix_image(moving_image) print("***") print( diff --git a/experiments/Lung-GatedCT_To_USD/Experiment_SegReg.py b/experiments/Lung-GatedCT_To_USD/Experiment_SegReg.py index ddb88fa..d01cffa 100644 --- a/experiments/Lung-GatedCT_To_USD/Experiment_SegReg.py +++ b/experiments/Lung-GatedCT_To_USD/Experiment_SegReg.py @@ -3,8 +3,6 @@ import itk -from data_dirlab_4d_ct import DataDirLab4DCT - from physiotwin4d import RegisterImagesICON from physiotwin4d import SegmentChestTotalSegmentator @@ -18,12 +16,10 @@ _RESULTS_DIR = os.path.join(_HERE, "results_SegReg") # %% - fixed_image = DataDirLab4DCT().fix_image( - itk.imread(os.path.join(_DATA_DIR, "Case1Pack_T30.mhd")) - ) - moving_image = DataDirLab4DCT().fix_image( - itk.imread(os.path.join(_DATA_DIR, "Case1Pack_T00.mhd")) - ) + # .mha files are DirLab-4DCT data already converted to HU by + # data/DirLab-4DCT/fix_downloaded_data.py. + fixed_image = itk.imread(os.path.join(_DATA_DIR, "Case1Pack_T30.mha")) + moving_image = itk.imread(os.path.join(_DATA_DIR, "Case1Pack_T00.mha")) # %% # Register images diff --git a/experiments/Lung-GatedCT_To_USD/data_dirlab_4d_ct.py b/experiments/Lung-GatedCT_To_USD/data_dirlab_4d_ct.py index 4e4f5d2..731d36d 100644 --- a/experiments/Lung-GatedCT_To_USD/data_dirlab_4d_ct.py +++ b/experiments/Lung-GatedCT_To_USD/data_dirlab_4d_ct.py @@ -2,9 +2,13 @@ class pmDataDirLab4dCT: This module contains the pmDataDirLab4DCT class, which is used to store the data for the DirLab 4DCT dataset. -""" -import itk +DirLab-4DCT's raw ``.mhd``/``.img`` volumes are not in Hounsfield units; run +``data/DirLab-4DCT/fix_downloaded_data.py`` (backed by +``DataDownloadTools.FixDirLab4DCTData``) once to write corrected ``.mha`` +volumes before using this dataset. That script's output is what this +directory's ``case_names`` are meant to be read from. +""" class DataDirLab4DCT: @@ -30,13 +34,3 @@ def __init__(self): def get_case_names(self) -> list[str]: """Get the case names""" return self.case_names - - def fix_image(self, input_image: itk.image) -> itk.image: - """Fix DirLab_4DCT intensities to conform to HU""" - input_image_arr = itk.GetArrayViewFromImage(input_image) - input_image_arr -= 1024 - input_image_arr = input_image_arr.clip(-1024, 3071) - new_image = itk.GetImageFromArray(input_image_arr) - new_image.CopyInformation(input_image) - - return new_image diff --git a/src/physiotwin4d/cli/fit_statistical_model_to_patient.py b/src/physiotwin4d/cli/fit_statistical_model_to_patient.py index 85976a9..adc1621 100644 --- a/src/physiotwin4d/cli/fit_statistical_model_to_patient.py +++ b/src/physiotwin4d/cli/fit_statistical_model_to_patient.py @@ -262,7 +262,7 @@ def main() -> int: # Execute registration workflow print("\nStarting registration pipeline...") print("=" * 70) - result = workflow.run_workflow( + result = workflow.process( use_ICON_registration_refinement=args.use_ICON_refinement, ) diff --git a/src/physiotwin4d/cli/reconstruct_highres_4d_ct.py b/src/physiotwin4d/cli/reconstruct_highres_4d_ct.py index 3e0fc7d..0d09819 100644 --- a/src/physiotwin4d/cli/reconstruct_highres_4d_ct.py +++ b/src/physiotwin4d/cli/reconstruct_highres_4d_ct.py @@ -312,7 +312,7 @@ def main() -> int: print("\nStarting reconstruction pipeline...") print("=" * 70) workflow.set_upsample_to_fixed_resolution(True) - result = workflow.run_workflow() + result = workflow.process() # Save results print("\n" + "=" * 70) diff --git a/src/physiotwin4d/contour_tools.py b/src/physiotwin4d/contour_tools.py index 2c46267..5d0c2cd 100644 --- a/src/physiotwin4d/contour_tools.py +++ b/src/physiotwin4d/contour_tools.py @@ -31,8 +31,8 @@ def __init__(self, log_level: int | str = logging.INFO): """ super().__init__(class_name=self.__class__.__name__, log_level=log_level) + @staticmethod def extract_contours( - self, labelmap_image: itk.image, ) -> pv.PolyData: """ @@ -71,8 +71,8 @@ def extract_contours( return contours + @staticmethod def transform_contours( - self, contours: pv.PolyData, tfm: itk.Transform, with_deformation_magnitude: bool = False, @@ -135,8 +135,8 @@ def merge_meshes( return merged_mesh, pv_meshes + @staticmethod def create_reference_image( - self, mesh: pv.DataSet, spatial_resolution: float = 0.5, buffer_factor: float = 0.25, @@ -164,8 +164,8 @@ def create_reference_image( reference_image.Allocate() return reference_image + @staticmethod def create_mask_from_mesh( - self, mesh: pv.DataSet | pv.UnstructuredGrid, reference_image: itk.Image, ) -> itk.Image: @@ -360,8 +360,8 @@ def create_distance_map( return distance_image + @staticmethod def create_deformation_field( - self, points: np.ndarray, point_displacements: np.ndarray, reference_image: itk.Image, diff --git a/src/physiotwin4d/data_download_tools.py b/src/physiotwin4d/data_download_tools.py index 879352e..da83d37 100644 --- a/src/physiotwin4d/data_download_tools.py +++ b/src/physiotwin4d/data_download_tools.py @@ -16,7 +16,10 @@ import urllib.request import zipfile from pathlib import Path -from typing import Union +from typing import Optional, Union + +import itk +import numpy as np from .convert_image_4d_to_3d import ConvertImage4DTo3D @@ -355,6 +358,54 @@ def VerifyDirLab4DCTData(dirname: Union[str, Path]) -> bool: # noqa: N802 ) return has_case_dir_layout or has_pack_layout + DIRLAB_4DCT_HU_OFFSET = 1024 + DIRLAB_4DCT_HU_CLIP_RANGE = (-1024, 1024) + + @staticmethod + def FixDirLab4DCTData( # noqa: N802 + dirname: Union[str, Path], + output_dirname: Optional[Union[str, Path]] = None, + ) -> list[Path]: + """Convert DirLab-4DCT raw intensities to clipped Hounsfield units. + + DirLab-4DCT's raw ``.img``/``.mhd`` volumes store intensities offset + by +1024 from Hounsfield units. This finds every ``.mhd`` header + under ``dirname`` (recursively) whose raw data has been downloaded, + subtracts 1024, clips the result to ``DIRLAB_4DCT_HU_CLIP_RANGE``, + and writes the corrected volume with the same base filename but a + ``.mha`` extension, compressed. + + Args: + dirname: Directory to search recursively for ``.mhd`` headers. + output_dirname: Directory to write the ``.mha`` files into. + Defaults to writing each result next to its source ``.mhd`` + header. + + Returns: + Paths to the written ``.mha`` files. + """ + output_dir = Path(output_dirname) if output_dirname is not None else None + if output_dir is not None: + output_dir.mkdir(parents=True, exist_ok=True) + + output_files = [] + for mhd_file in sorted(Path(dirname).rglob("*.mhd")): + if not DataDownloadTools._MetaImageHeaderHasBackingData(mhd_file): + continue + image = itk.imread(str(mhd_file)) + image_arr = ( + itk.array_from_image(image) - DataDownloadTools.DIRLAB_4DCT_HU_OFFSET + ) + image_arr = np.clip(image_arr, *DataDownloadTools.DIRLAB_4DCT_HU_CLIP_RANGE) + fixed_image = itk.image_from_array(image_arr) + fixed_image.CopyInformation(image) + target_dir = output_dir if output_dir is not None else mhd_file.parent + output_file = target_dir / mhd_file.with_suffix(".mha").name + itk.imwrite(fixed_image, str(output_file), compression=True) + output_files.append(output_file) + _logger.info("Fixed %s -> %s", mhd_file.name, output_file.name) + return output_files + @staticmethod def VerifyKCLHeartModelData(dirname: Union[str, Path]) -> bool: # noqa: N802 """Return True when KCL-Heart-Model has its expected mesh inputs.""" diff --git a/src/physiotwin4d/segment_anatomy_base.py b/src/physiotwin4d/segment_anatomy_base.py index c0f6865..0afb8ef 100644 --- a/src/physiotwin4d/segment_anatomy_base.py +++ b/src/physiotwin4d/segment_anatomy_base.py @@ -348,51 +348,55 @@ def postprocess_after_labelmap( """ return labelmap_image - def create_anatomy_group_masks( + def create_anatomy_group_labelmaps( self, labelmap_image: itk.image ) -> dict[str, itk.image]: """ - Create binary masks for different anatomical groups from the labelmap. + Create labelmaps for different anatomical groups from the labelmap. - Generates separate binary masks for major anatomical systems by + Generates separate labelmaps for major anatomical systems by grouping related anatomical structures from the detailed labelmap. - This is useful for motion analysis and visualization. + Each group's labelmap retains the original label ids for voxels + belonging to that group and is zero elsewhere. This is useful for + motion analysis and visualization. Args: labelmap_image (itk.image): The detailed segmentation labelmap Returns: - dict[str, itk.image]: Dictionary of binary masks keyed by group + dict[str, itk.image]: Dictionary of labelmaps keyed by group name. Exactly one entry per group registered in :attr:`taxonomy` (plus ``"other"``). The returned key set is segmenter-specific — callers that need a particular - group should check membership (``"lung" in masks``) rather - than assume a fixed schema. + group should check membership (``"lung" in labelmaps``) + rather than assume a fixed schema. Example: - >>> masks = segmenter.create_anatomy_group_masks(labelmap) - >>> if "lung" in masks: - ... lung_mask = masks["lung"] + >>> labelmaps = segmenter.create_anatomy_group_labelmaps(labelmap) + >>> if "lung" in labelmaps: + ... lung_labelmap = labelmaps["lung"] """ labelmap_arr = itk.GetArrayFromImage(labelmap_image) - other_mask_arr = np.where(labelmap_arr > 0, 1, 0) + other_labelmap_arr = np.where(labelmap_arr > 0, labelmap_arr, 0) - masks: dict[str, itk.image] = {} + labelmaps: dict[str, itk.image] = {} for group_name in self.taxonomy.group_names(): if group_name == AnatomyTaxonomy.OTHER_GROUP: continue group_ids = list(self.taxonomy.labels_in_group(group_name).keys()) - group_mask_arr = np.isin(labelmap_arr, group_ids).astype(np.uint8) - other_mask_arr = np.where(group_mask_arr > 0, 0, other_mask_arr) - group_mask = itk.GetImageFromArray(group_mask_arr) - group_mask.CopyInformation(labelmap_image) - masks[group_name] = group_mask + group_labelmap_arr = np.where( + np.isin(labelmap_arr, group_ids), labelmap_arr, 0 + ) + other_labelmap_arr = np.where(group_labelmap_arr > 0, 0, other_labelmap_arr) + group_labelmap = itk.GetImageFromArray(group_labelmap_arr) + group_labelmap.CopyInformation(labelmap_image) + labelmaps[group_name] = group_labelmap - other_mask = itk.GetImageFromArray(other_mask_arr.astype(np.uint8)) - other_mask.CopyInformation(labelmap_image) - masks[AnatomyTaxonomy.OTHER_GROUP] = other_mask + other_labelmap = itk.GetImageFromArray(other_labelmap_arr) + other_labelmap.CopyInformation(labelmap_image) + labelmaps[AnatomyTaxonomy.OTHER_GROUP] = other_labelmap - return masks + return labelmaps def segmentation_method(self, preprocessed_image: itk.image) -> itk.image: """ @@ -425,7 +429,7 @@ def segment( This is the main segmentation method that coordinates preprocessing, segmentation, subclass-specific labelmap refinement, and anatomical - group mask creation. + group labelmap creation. Args: input_image (itk.image): The input 3D image to segment @@ -433,12 +437,13 @@ def segment( Returns: dict[str, itk.image]: Dictionary containing: - "labelmap": Detailed segmentation labelmap - - one binary mask image per anatomy group, keyed by group name + - one labelmap image per anatomy group, keyed by group name, + preserving the original label ids for that group Example: >>> result = segmenter.segment(image) >>> labelmap = result['labelmap'] - >>> heart_mask = result['heart'] + >>> heart_labelmap = result['heart'] """ preprocessed_image = self.preprocess_input(input_image) @@ -448,11 +453,11 @@ def segment( labelmap_image = self.postprocess_after_labelmap(input_image, labelmap_image) - masks = self.create_anatomy_group_masks(labelmap_image) + labelmaps = self.create_anatomy_group_labelmaps(labelmap_image) labelmap_image = itk.GetImageFromArray( itk.GetArrayFromImage(labelmap_image).astype(np.uint8) ) labelmap_image.CopyInformation(input_image) - return {"labelmap": labelmap_image, **masks} + return {"labelmap": labelmap_image, **labelmaps} diff --git a/src/physiotwin4d/segment_chest_total_segmentator.py b/src/physiotwin4d/segment_chest_total_segmentator.py index cd0bc75..a8ba499 100644 --- a/src/physiotwin4d/segment_chest_total_segmentator.py +++ b/src/physiotwin4d/segment_chest_total_segmentator.py @@ -14,6 +14,7 @@ import nibabel as nib import numpy as np +from .image_tools import ImageTools from .segment_anatomy_base import SegmentAnatomyBase @@ -49,7 +50,7 @@ class SegmentChestTotalSegmentator(SegmentAnatomyBase): >>> segmenter = SegmentChestTotalSegmentator() >>> result = segmenter.segment(ct_image) >>> labelmap = result['labelmap'] - >>> heart_mask = result['heart'] + >>> heart_labelmap = result['heart'] """ def __init__(self, log_level: int | str = logging.INFO): @@ -65,7 +66,7 @@ def __init__(self, log_level: int | str = logging.INFO): """ super().__init__(log_level=log_level) - self.target_spacing = 1.5 + self.target_spacing = 0.0 # TotalSegmentator class indices, grouped by anatomy. for group_name, organs in ( @@ -74,13 +75,19 @@ def __init__(self, log_level: int | str = logging.INFO): { 51: "heart", 61: "atrial_appendage_left", - 140: "heart_envelop", + 140: "highres_myocardium", + 141: "highres_atrium_left", + 142: "highres_ventricle_left", + 143: "highres_atrium_right", + 144: "highres_ventricle_right", + 146: "highres_pulmonary_artery", }, ), ( "major_vessels", { 52: "aorta", + 145: "highres_aorta", 53: "pulmonary_vein", 54: "brachiocephalic_trunk", 55: "right_subclavian_artery", @@ -202,7 +209,10 @@ def __init__(self, log_level: int | str = logging.INFO): 90: "brain", 15: "esophagus", 16: "trachea", - 133: "soft_tissue", + 133: "body", + 134: "body_trunc", + 135: "body_extremities", + 136: "body_skin", }, ), ): @@ -212,6 +222,16 @@ def __init__(self, log_level: int | str = logging.INFO): self._add_extra_taxonomy_groups() self._finalize_other_group() + self.has_highres_heart_license = False + + def set_has_highres_heart_license(self, has_highres_heart_license: bool) -> None: + """Set whether the highres heart license is available. + + Args: + has_highres_heart_license (bool): Whether the highres heart license is available + """ + self.has_highres_heart_license = has_highres_heart_license + def _add_extra_taxonomy_groups(self) -> None: """Hook for subclasses to add taxonomy groups before finalization. @@ -246,7 +266,7 @@ def segmentation_method(self, preprocessed_image: itk.image) -> itk.image: soft tissue labels from the 'body' task Note: - Requires GPU acceleration (device="gpu:0") for reasonable performance. + Requires GPU acceleration (device="gpu") for reasonable performance. The method automatically handles coordinate system conversions between ITK and nibabel formats. @@ -270,60 +290,106 @@ def segmentation_method(self, preprocessed_image: itk.image) -> itk.image: # nr_thr_resamp defaults to 1; TotalSegmentator's post-prediction # resampling back to native resolution is CPU-bound and benefits # from parallelizing across the available cores. - resamp_threads = min(8, os.cpu_count() or 1) if self.fast_mode else 1 - output_nib_image1 = totalsegmentator( + resamp_threads = min(8, os.cpu_count() or 1) + output_nib_image_total = totalsegmentator( nib_image, task="total", device="gpu", fast=self.fast_mode, nr_thr_resamp=resamp_threads, ) - labelmap_arr1 = output_nib_image1.get_fdata().astype(np.uint8) - - if self.fast_mode: - final_arr = labelmap_arr1 - else: - output_nib_image2 = totalsegmentator( - nib_image, task="body", device="gpu" + labelmap_arr_total = output_nib_image_total.get_fdata().astype(np.uint8) + + final_arr = labelmap_arr_total + + if not self.fast_mode: + if self.has_highres_heart_license: + self.log_info("Running heart chambers task") + output_nib_image_heart = totalsegmentator( + nib_image, + task="heartchambers_highres", + device="gpu", + nr_thr_resamp=resamp_threads, + ) + labelmap_arr_heart = output_nib_image_heart.get_fdata().astype( + np.uint8 + ) + # labelmap_arr_heart contains: 1=myocardium, 2=atrium_left, 3=ventricle_left, + # 4=atrium_right, 5=ventricle_right, 6=aorta, 7=pulmonary_artery + final_arr = np.where(labelmap_arr_heart == 1, 140, final_arr) + final_arr = np.where(labelmap_arr_heart == 2, 141, final_arr) + final_arr = np.where(labelmap_arr_heart == 3, 142, final_arr) + final_arr = np.where(labelmap_arr_heart == 4, 143, final_arr) + final_arr = np.where(labelmap_arr_heart == 5, 144, final_arr) + final_arr = np.where(labelmap_arr_heart == 7, 146, final_arr) + # final_arr = np.where(labelmap_arr_heart == 6, 145, final_arr) + # Aorta is not included in heart model. + # Should include only a portion of the aorta in the heart model. + + self.log_info("Running lung vessels task") + output_nib_image_lung = totalsegmentator( + nib_image, + task="lung_vessels", + device="gpu", + nr_thr_resamp=resamp_threads, ) - labelmap_arr2 = output_nib_image2.get_fdata().astype(np.uint8) - - output_nib_image3 = totalsegmentator( - nib_image, task="lung_vessels", device="gpu" + labelmap_arr_lung = output_nib_image_lung.get_fdata().astype(np.uint8) + # labelmap_arr_lung contains: 1=arteries, 2=veins, 3=airways, + # 4=airways_wall + final_arr = np.where(labelmap_arr_lung == 1, 120, final_arr) + final_arr = np.where(labelmap_arr_lung == 2, 121, final_arr) + final_arr = np.where(labelmap_arr_lung == 3, 122, final_arr) + # final_arr = np.where(labelmap_arr_lung == 4, 123, final_arr) + # Airway wall segmentation is too zealous. Fills right atrium + + self.log_info("Running body task") + output_nib_image_body = totalsegmentator( + nib_image, task="body", device="gpu", nr_thr_resamp=resamp_threads ) - labelmap_arr3 = output_nib_image3.get_fdata().astype(np.uint8) - - mask1 = labelmap_arr1 == 0 - mask2 = labelmap_arr2 > 0 - mask = mask1 & mask2 - soft_tissue_id = next( - label_id - for label_id, name in self.taxonomy.labels_in_group( - "soft_tissue" - ).items() - if name == "soft_tissue" - ) - final_arr = np.where(mask, soft_tissue_id, labelmap_arr1) - - # labelmap_arr3 contains: 1=arteries, 2=veins, 3=airways, - # 4=airways_wall - final_arr = np.where( - labelmap_arr3 == 1, 120, final_arr - ) # lung arteries - final_arr = np.where(labelmap_arr3 == 2, 121, final_arr) # lung veins - final_arr = np.where(labelmap_arr3 == 3, 122, final_arr) # lung airways - final_arr = np.where( - labelmap_arr3 == 4, 123, final_arr - ) # lung airways wall + labelmap_arr_body = output_nib_image_body.get_fdata().astype(np.uint8) + # labelmap_arr_body contains: 1=body, 2=body_trunc, 3=body_extremities, + # 4=skin + # Only overwrite the background with body labels + mask = final_arr > 0 + labelmap_arr_body[mask] = 0 + final_arr = np.where(labelmap_arr_body == 1, 133, final_arr) + final_arr = np.where(labelmap_arr_body == 2, 134, final_arr) + final_arr = np.where(labelmap_arr_body == 3, 135, final_arr) + final_arr = np.where(labelmap_arr_body == 4, 136, final_arr) + # To create an ITK image, we save the result and read it back with # ITK. This correctly handles the coordinate system and data # layout conversions. out_tmp_file = os.path.join(tmp_dir, "out.nii.gz") # Use the affine from one of the outputs to preserve spatial info - result_nib = nib.Nifti1Image(final_arr, output_nib_image1.affine) + result_nib = nib.Nifti1Image(final_arr, output_nib_image_total.affine) nib.save(result_nib, out_tmp_file) labelmap_image = itk.imread(out_tmp_file) labelmap_arr = itk.array_from_image(labelmap_image).astype(np.uint8) + + # Add heart around interior regions. + if self.has_highres_heart_license: + interior_mask = np.isin(labelmap_arr, [141, 142, 143, 144]) + # Binarize to foreground value 1 so the dilate/erode calls + # below (which use foreground=1) operate on the mask. + interior_arr = interior_mask.astype(np.uint8) + interior_image = itk.GetImageFromArray(interior_arr) + interior_image.CopyInformation(preprocessed_image) + imMath = ImageTools() + spacing = interior_image.GetSpacing() + exterior_image = imMath.binary_dilate_image( + interior_image, round(7 / spacing[0]), 1, 0 + ) + exterior_image = imMath.binary_erode_image( + exterior_image, round(4 / spacing[0]), 1, 0 + ) + exterior_arr = itk.GetArrayFromImage(exterior_image) + mask_id = 51 # Heart mask id + exterior_arr = exterior_arr * mask_id + labelmap_arr = np.where(labelmap_arr == 0, exterior_arr, labelmap_arr) + replace_arr = np.where(labelmap_arr == 133, exterior_arr, 0) + labelmap_arr = np.where(replace_arr > 0, exterior_arr, labelmap_arr) + labelmap_image = itk.image_from_array(labelmap_arr) labelmap_image.CopyInformation(preprocessed_image) diff --git a/src/physiotwin4d/segment_chest_total_segmentator_with_contrast.py b/src/physiotwin4d/segment_chest_total_segmentator_with_contrast.py index 96d5710..101e665 100644 --- a/src/physiotwin4d/segment_chest_total_segmentator_with_contrast.py +++ b/src/physiotwin4d/segment_chest_total_segmentator_with_contrast.py @@ -23,7 +23,7 @@ class SegmentChestTotalSegmentatorWithContrast(SegmentChestTotalSegmentator): Extends :class:`SegmentChestTotalSegmentator` with an additional connected-component pass that identifies contrast-enhanced blood vessels and cardiac chambers, labeling them under a ``"contrast"`` taxonomy - group (label id 135). Use this class instead of + group (label id 155). Use this class instead of :class:`SegmentChestTotalSegmentator` for contrast-enhanced studies. Attributes: @@ -34,7 +34,7 @@ class SegmentChestTotalSegmentatorWithContrast(SegmentChestTotalSegmentator): >>> segmenter = SegmentChestTotalSegmentatorWithContrast() >>> result = segmenter.segment(ct_image) >>> labelmap = result['labelmap'] - >>> contrast_mask = result['contrast'] + >>> contrast_labelmap = result['contrast'] """ def __init__(self, log_level: int | str = logging.INFO) -> None: @@ -48,8 +48,8 @@ def __init__(self, log_level: int | str = logging.INFO) -> None: self.contrast_threshold: int = 500 def _add_extra_taxonomy_groups(self) -> None: - """Register the ``"contrast"`` taxonomy group (label id 135).""" - self.taxonomy.add_organ("contrast", 135, "contrast") + """Register the ``"contrast"`` taxonomy group (label id 155).""" + self.taxonomy.add_organ("contrast", 155, "contrast") def postprocess_after_labelmap( self, input_image: itk.Image, labelmap_image: itk.Image @@ -104,7 +104,7 @@ def segment_connected_component( Example: >>> # Segment contrast-enhanced blood >>> updated_labels = segmenter.segment_connected_component( - ... preprocessed_image, labels, 700, 4000, mask_id=135 + ... preprocessed_image, labels, 700, 4000, mask_id=155 ... ) """ thresh_image = itk.binary_threshold_image_filter( diff --git a/src/physiotwin4d/segment_heart_simpleware.py b/src/physiotwin4d/segment_heart_simpleware.py index d569f35..b73f2d9 100644 --- a/src/physiotwin4d/segment_heart_simpleware.py +++ b/src/physiotwin4d/segment_heart_simpleware.py @@ -52,7 +52,7 @@ class SegmentHeartSimpleware(SegmentAnatomyBase): >>> segmenter = SegmentHeartSimpleware() >>> result = segmenter.segment(ct_image) >>> labelmap = result['labelmap'] - >>> heart_mask = result['heart'] + >>> heart_labelmap = result['heart'] See :class:`SegmentHeartSimplewareTrimmedBranches` for a variant that additionally clips pulmonary/great-vessel branches to the cardiac region. diff --git a/src/physiotwin4d/workflow_convert_image_to_usd.py b/src/physiotwin4d/workflow_convert_image_to_usd.py index d578f97..5d3dd6c 100644 --- a/src/physiotwin4d/workflow_convert_image_to_usd.py +++ b/src/physiotwin4d/workflow_convert_image_to_usd.py @@ -54,7 +54,7 @@ def __init__( registration_method: Optional[RegisterImagesBase] = None, dynamic_labelmap_ids: Optional[list[int]] = None, mask_dilation_radius: int = 10, - frames_per_second: float = 24.0, + frames_per_second: float = 1.0, log_level: int | str = logging.INFO, save_assets: bool = True, ) -> None: @@ -137,7 +137,7 @@ def __init__( # Data storage for processing pipeline self._num_time_points = len(time_series_images) - self.reference_segmentation: Optional[dict[str, itk.Image]] = None + self.reference_segmentation_results: Optional[dict[str, itk.Image]] = None self.reference_contours: dict[str, pv.PolyData] = {} self.transformed_contours: dict[str, list[pv.PolyData]] = { @@ -225,8 +225,10 @@ def _segment_and_register_frames(self) -> None: # Set up registrar with reference image self.registrar.set_fixed_image(self.reference_image) - self.reference_segmentation = self.segmenter.segment(self.reference_image) - labelmap = self.reference_segmentation["labelmap"] + self.reference_segmentation_results = self.segmenter.segment( + self.reference_image + ) + labelmap = self.reference_segmentation_results["labelmap"] if self.save_assets: itk.imwrite( labelmap, @@ -266,8 +268,8 @@ def _segment_and_register_frames(self) -> None: moving_image = self.time_series_images[i] - moving_segmentation = self.segmenter.segment(moving_image) - moving_labelmap = moving_segmentation["labelmap"] + moving_segmentation_results = self.segmenter.segment(moving_image) + moving_labelmap = moving_segmentation_results["labelmap"] if self.save_assets: itk.imwrite( moving_labelmap, @@ -345,10 +347,10 @@ def _generate_reference_contours(self) -> None: """Generate contour meshes from reference segmentation.""" self.log_info("Generating reference contours...") - assert self.reference_segmentation is not None, ( - "reference segmentation must be set" + assert self.reference_segmentation_results is not None, ( + "reference segmentation results must be set" ) - labelmap = self.reference_segmentation["labelmap"] + labelmap = self.reference_segmentation_results["labelmap"] # Generate all anatomy contours all_contours = self.contour_tools.extract_contours(labelmap) diff --git a/src/physiotwin4d/workflow_fit_statistical_model_to_patient.py b/src/physiotwin4d/workflow_fit_statistical_model_to_patient.py index 94c25ab..44041aa 100644 --- a/src/physiotwin4d/workflow_fit_statistical_model_to_patient.py +++ b/src/physiotwin4d/workflow_fit_statistical_model_to_patient.py @@ -74,7 +74,7 @@ class WorkflowFitStatisticalModelToPatient(PhysioTwin4DBase): patient_model_surface (pv.PolyData): Primary patient model surface (first in list) combined_patient_model (pv.PolyData): Merged patient models before surface - extraction; used when pca_uses_surface=False. + extraction; used when use_surface=False. patient_image (itk.Image): Reference image providing coordinate frame patient_labelmap (itk.Image): Multi-label labelmap for patient model patient_mask (itk.Image): Binary mask for patient registration region @@ -115,11 +115,11 @@ class WorkflowFitStatisticalModelToPatient(PhysioTwin4DBase): ... patient_models=[lv_model, mc_model, rv_model], ... ) >>> registrar.set_mask_dilation_mm(20) - >>> # To enable PCA registration, call before run_workflow(): + >>> # To enable PCA registration, call before process(): >>> # registrar.set_use_pca_registration(True, pca_model=pca_model_dict, pca_number_of_modes=10) >>> # To enable labelmap-to-image refinement: >>> # registrar.set_use_labelmap_to_image_registration(True, template_labelmap, organ_mesh_ids, organ_extra_ids, background_ids) - >>> result = registrar.run_workflow() + >>> result = registrar.process() """ def __init__( @@ -264,7 +264,7 @@ def __init__( self.pca_template_model: Optional[pv.DataSet] = None self.pca_template_model_surface: Optional[pv.PolyData] = None self.pca_template_labelmap: Optional[itk.Image] = None - self.pca_uses_surface: bool = False + self.use_surface: bool = False # Stage 2: Labelmap-to-labelmap registration results self.use_l2l_registration = True @@ -301,7 +301,7 @@ def set_use_pca_registration( use_pca_registration: bool, pca_model: Optional[dict[str, Any]] = None, pca_number_of_modes: int = 0, - pca_uses_surface: bool = True, + use_surface: bool = False, ) -> None: """Set whether to use PCA-based registration and provide the PCA model. @@ -314,7 +314,7 @@ def set_use_pca_registration( "eigenvalues" and "components". pca_number_of_modes: Required when use is True. Number of PCA modes to use. Default 0 means use all modes. - pca_uses_surface: Whether to use the surface of the patient model for PCA registration. + use_surface: Whether to use the surface of the patient model for PCA registration. Raises: ValueError: If use is True and pca_model is None. """ @@ -328,7 +328,7 @@ def set_use_pca_registration( else: self.pca_model = None self.pca_number_of_modes = 0 - self.pca_uses_surface = pca_uses_surface + self.use_surface = use_surface self.use_pca_registration = use_pca_registration def set_use_labelmap_to_labelmap_registration( @@ -515,7 +515,7 @@ def register_model_to_model_pca(self) -> dict: } pca_template_model: Optional[pv.DataSet] - if self.pca_uses_surface: + if self.use_surface: pca_template_model = self.icp_template_model_surface fixed_model = self.patient_model_surface fixed_distance_map = None @@ -549,7 +549,7 @@ def register_model_to_model_pca(self) -> dict: result = self.pca_registrar.register() self.pca_coefficients = result["pca_coefficients"] registered_model = cast(pv.DataSet, result["registered_model"]) - if self.pca_uses_surface: + if self.use_surface: self.pca_template_model_surface = cast(pv.PolyData, registered_model) else: self.pca_template_model_surface = registered_model.extract_surface( @@ -579,7 +579,7 @@ def register_model_to_model_pca(self) -> dict: itk.imwrite(tfm_y_img, "pca_forward_point_transform_y.nii.gz") itk.imwrite(tfm_z_img, "pca_forward_point_transform_z.nii.gz") - if self.pca_uses_surface: + if self.use_surface: assert self.icp_template_model is not None, "ICP template model must be set" self.pca_template_model = self._transform_model_dataset( self.icp_template_model, @@ -710,7 +710,7 @@ def register_labelmap_to_image( raise ValueError( "Labelmap-to-image registration requires template labelmap and label IDs. " "Call set_use_labelmap_to_image_registration(True, template_labelmap, " - "organ_mesh_ids, organ_extra_ids, background_ids) before run_workflow()." + "organ_mesh_ids, organ_extra_ids, background_ids) before process()." ) propagated_labelmap = ( self.l2l_template_labelmap @@ -886,7 +886,7 @@ def transform_model( transformed_model.points = new_points return transformed_model - def run_workflow( + def process( self, use_ICON_registration_refinement: bool = False, ) -> dict: diff --git a/src/physiotwin4d/workflow_reconstruct_highres_4d_ct.py b/src/physiotwin4d/workflow_reconstruct_highres_4d_ct.py index c074cfe..49d9ab9 100644 --- a/src/physiotwin4d/workflow_reconstruct_highres_4d_ct.py +++ b/src/physiotwin4d/workflow_reconstruct_highres_4d_ct.py @@ -346,7 +346,7 @@ def reconstruct_time_series(self) -> dict: return {"reconstructed_images": self.reconstructed_images} - def run_workflow(self) -> dict: + def process(self) -> dict: """Execute the complete high-resolution 4D CT reconstruction workflow. Runs the full pipeline: diff --git a/tests/conftest.py b/tests/conftest.py index af220d8..9f34b8b 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -515,17 +515,16 @@ def test_labelmaps( itk.imwrite(result["labelmap"], str(labelmap_file), compression=True) labelmap = itk.imread(str(labelmap_file)) - masks = segmenter_total_segmentator.create_anatomy_group_masks(labelmap) + labelmaps = segmenter_total_segmentator.create_anatomy_group_labelmaps(labelmap) results.append( { "labelmap": labelmap, - "lung": masks["lung"], - "heart": masks["heart"], - "major_vessels": masks["major_vessels"], - "bone": masks["bone"], - "soft_tissue": masks["soft_tissue"], - "other": masks["other"], - "contrast": masks["contrast"], + "lung": labelmaps["lung"], + "heart": labelmaps["heart"], + "major_vessels": labelmaps["major_vessels"], + "bone": labelmaps["bone"], + "soft_tissue": labelmaps["soft_tissue"], + "other": labelmaps["other"], } ) diff --git a/tests/test_segment_chest_total_segmentator.py b/tests/test_segment_chest_total_segmentator.py index e11c332..804bcb1 100644 --- a/tests/test_segment_chest_total_segmentator.py +++ b/tests/test_segment_chest_total_segmentator.py @@ -138,18 +138,18 @@ def test_segment_multiple_images( assert len(results) == 2, "Expected 2 segmentation results" print(f"\nSuccessfully segmented {len(results)} time points") - def test_anatomy_group_masks( + def test_anatomy_group_labelmaps( self, segmenter_total_segmentator: SegmentChestTotalSegmentator, test_images: list[Any], ) -> None: - """Test that anatomy group masks are created correctly.""" + """Test that anatomy group labelmaps are created correctly.""" input_image = test_images[0] # Run segmentation result = segmenter_total_segmentator.segment(input_image) - # Check each anatomy group mask + # Check each anatomy group labelmap anatomy_groups = [ "lung", "heart", @@ -158,26 +158,34 @@ def test_anatomy_group_masks( "soft_tissue", "other", ] + taxonomy = segmenter_total_segmentator.taxonomy for group in anatomy_groups: - mask = result[group] - assert mask is not None, f"{group} mask is None" - - # Check that mask is binary - mask_arr = itk.array_from_image(mask) - unique_values = np.unique(mask_arr) - assert len(unique_values) <= 2, f"{group} mask should be binary" - assert 0 in unique_values, f"{group} mask should contain background" - - # Check that mask has same size as input - assert itk.size(mask) == itk.size(input_image), ( - f"{group} mask size mismatch" + group_labelmap = result[group] + assert group_labelmap is not None, f"{group} labelmap is None" + + group_labelmap_arr = itk.array_from_image(group_labelmap) + unique_values = set(np.unique(group_labelmap_arr).tolist()) + assert 0 in unique_values, f"{group} labelmap should contain background" + + # "other" collects whatever ids no group claimed, so it has no + # fixed id set to check against. + if group != "other": + allowed_values = {0} | set(taxonomy.labels_in_group(group).keys()) + assert unique_values <= allowed_values, ( + f"{group} labelmap contains unexpected label ids: " + f"{unique_values - allowed_values}" + ) + + # Check that labelmap has same size as input + assert itk.size(group_labelmap) == itk.size(input_image), ( + f"{group} labelmap size mismatch" ) - print("\nAll anatomy group masks created correctly") + print("\nAll anatomy group labelmaps created correctly") for group in anatomy_groups: - mask_arr = itk.array_from_image(result[group]) - num_voxels = np.sum(mask_arr > 0) + group_labelmap_arr = itk.array_from_image(result[group]) + num_voxels = np.sum(group_labelmap_arr > 0) print(f" {group}: {num_voxels} voxels") def test_contrast_detection( diff --git a/tutorials/tutorial_01b_lung_gated_ct_to_usd.py b/tutorials/tutorial_01b_lung_gated_ct_to_usd.py new file mode 100644 index 0000000..00add8c --- /dev/null +++ b/tutorials/tutorial_01b_lung_gated_ct_to_usd.py @@ -0,0 +1,211 @@ +""" +Tutorial 1b: Lung-Gated 4D CT to Animated USD + +Purpose +------- +Convert a respiratory-gated 4D lung CT scan (multiple breathing phases) into an +animated USD model suitable for visualization in NVIDIA Omniverse. The workflow +segments the lungs and surrounding chest anatomy from a reference phase, +registers all other respiratory phases to that reference using deep-learning +registration, and assembles the resulting time-varying surface meshes into a +single USD file with anatomical materials applied. + +Inputs +------ +- A set of 3D CT volumes (``*.mha``) representing successive respiratory + phases of one DirLab-4DCT case. + Expected location: ``data/DirLab-4DCT/Case1Pack_T??.mha`` (already converted + to Hounsfield units by ``data/DirLab-4DCT/fix_downloaded_data.py``). +- The mid-inspiration phase (index ~0.7 through the series) is used as the + reference frame for segmentation and registration. + +Outputs +------- +- An animated USD file with anatomy materials, written under ``output_dir`` + and named after the workflow's ``usd_project_name``. +- Screenshots (PNG) for documentation and regression testing: + - ``slice__registered_test.png`` - axial slice of the registered + reference phase + - ``slice__labelmap_test.png`` - segmentation mask overlaid on that slice + - a rendered view of the exported USD model + +Strengths +--------- +- Single call (``WorkflowConvertImageToUSD.process()``) runs the full pipeline. +- Supports both GPU-accelerated ICON registration and CPU-capable Greedy registration. +- Output is Omniverse-ready with anatomical materials (USDAnatomyTools). + +Weaknesses / Limitations +------------------------ +- Requires a GPU for ICON registration (``registration_method=RegisterImagesICON()``); + use ``registration_method=RegisterImagesGreedy()`` for CPU-only environments + (about 10x slower). +- Segmentation quality depends on TotalSegmentator's training distribution; + unusual pathologies or pediatric anatomy may degrade results. +- Large 4D datasets (>20 phases, high resolution) can require 32 GB+ RAM. + +Classes Used +------------ +- WorkflowConvertImageToUSD (workflow_convert_image_to_usd.py): + Orchestrates the full pipeline: CT phases -> segmentation -> registration -> + contour extraction -> USD export. +- SegmentChestTotalSegmentator (segment_chest_total_segmentator.py): + Deep-learning segmentation of 117 anatomical structures (used internally). +- RegisterImagesICON (register_images_icon.py): + Frame-to-frame image registration (used internally). +- ContourTools (contour_tools.py): + Extracts and transforms surface meshes from segmentation masks (used internally). +- USDAnatomyTools (usd_anatomy_tools.py): + Applies clinical material colours to USD prims (used internally). + +Data Required +------------- +See data/README.md for download instructions and dataset licensing. +Dataset: DirLab-4DCT - see ``data/DirLab-4DCT/README.md``. +This script expects the HU-corrected ``Case1Pack_T??.mha`` phase volumes to +already exist under ``data/DirLab-4DCT/``. Download the DirLab-4DCT case and run +``data/DirLab-4DCT/fix_downloaded_data.py`` before running this tutorial. +""" + +# %% +# Imports +from __future__ import annotations + +import logging +from pathlib import Path + +import itk + +from physiotwin4d.register_images_icon import RegisterImagesICON +from physiotwin4d.test_tools import TestTools +from physiotwin4d.workflow_convert_image_to_usd import ( + WorkflowConvertImageToUSD, +) + +# nnUNetv2 (used by TotalSegmentator inside WorkflowConvertImageToUSD) +# spawns a multiprocessing.Pool. On Windows the spawn start method re-imports +# this script in each child; without the __name__ == "__main__" guard around +# the top-level work, that re-import fires workflow.process() again and +# Python's spawn-cascade detector raises RuntimeError. +if __name__ == "__main__": + # %% + # Data directory specification + REPO_ROOT = Path(__file__).resolve().parent.parent + TUTORIALS_DIR = Path(__file__).resolve().parent + DATA_DIR = REPO_ROOT / "data" + FULL_DATA_DIR = DATA_DIR / "DirLab-4DCT" + TEST_DATA_DIR = DATA_DIR / "test" / "DirLab-4DCT" + OUTPUT_DIR = TUTORIALS_DIR / "output" / "tutorial_01b" + LOG_LEVEL = logging.INFO + + # %% + # Data reading + test_mode = TestTools.running_as_test() + + data_dir = TEST_DATA_DIR if test_mode else FULL_DATA_DIR + output_dir = OUTPUT_DIR + log_level = LOG_LEVEL + + output_dir.mkdir(parents=True, exist_ok=True) + + if test_mode: + number_of_registration_iterations = 1 + else: + number_of_registration_iterations = 10 + + # %% + # .mha files are DirLab-4DCT data already converted to HU by + # data/DirLab-4DCT/fix_downloaded_data.py. + frame_files = sorted(data_dir.glob("Case1Pack_T??.mha")) + if test_mode: + frame_files = frame_files[:2] + + input_filenames = [str(path) for path in frame_files] + if not input_filenames: + raise FileNotFoundError( + "DirLab-4DCT data not found. Checked:\n" + + f" - {data_dir}" + + "\n" + + "See data/README.md for download instructions." + ) + + time_series_images = [itk.imread(str(path)) for path in input_filenames] + reference_image = time_series_images[int(0.7 * len(time_series_images))] + + print("Number of time-series images:", len(time_series_images)) + + # %% + # Workflow initialization + registration_method = RegisterImagesICON(log_level=log_level) + registration_method.set_number_of_iterations(number_of_registration_iterations) + + workflow = WorkflowConvertImageToUSD( + time_series_images=time_series_images, + reference_image=reference_image, + output_directory=str(output_dir), + usd_project_name="lung_model", + registration_method=registration_method, + log_level=log_level, + save_assets=True, + ) + + # %% + # Workflow execution + usd_files = workflow.process() + # if dynamic_labelmap_ids is not None, there are two USD files + if len(workflow.dynamic_labelmap_ids) > 0: + usd_file = output_dir / usd_files["dynamic"] + else: + usd_file = output_dir / usd_files["all"] + + # %% + # Result saving + tt = TestTools( + class_name="tutorial_01b_lung_gated_ct_to_usd", + results_dir=output_dir, + log_level=log_level, + ) + + screenshots: list[Path] = [] + + test_image_num = int(0.7 * len(input_filenames)) + test_image_path = output_dir / f"slice_{test_image_num:03d}_registered.mha" + if test_image_path.exists(): + test_image = itk.imread(str(test_image_path)) + screenshots.append( + tt.save_screenshot_image_slice( + test_image, + f"slice_{test_image_num:03d}_registered_test.png", + axis=0, + slice_fraction=0.5, + colormap="gray", + vmin=-200, + vmax=600, + ) + ) + + test_labelmap_path = output_dir / f"slice_{test_image_num:03d}_labelmap.mha" + if test_labelmap_path.exists(): + test_labelmap = itk.imread(str(test_labelmap_path)) + screenshots.append( + tt.save_screenshot_image_slice( + test_image, + f"slice_{test_image_num:03d}_labelmap_test.png", + axis=0, + slice_fraction=0.5, + colormap="gray", + vmin=-200, + vmax=600, + overlay_mask=test_labelmap, + ) + ) + + if usd_file.exists(): + screenshots.append( + tt.save_screenshot_openusd( + usd_file, + "lung_model_test.png", + ) + ) + + tutorial_results = {"usd_file": str(usd_file), "screenshots": screenshots} diff --git a/tutorials/tutorial_04_fit_statistical_model_to_patient.py b/tutorials/tutorial_04_fit_statistical_model_to_patient.py index 2743091..6cb78f2 100644 --- a/tutorials/tutorial_04_fit_statistical_model_to_patient.py +++ b/tutorials/tutorial_04_fit_statistical_model_to_patient.py @@ -22,106 +22,134 @@ from pathlib import Path from typing import Any, Optional, cast +import itk import pyvista as pv -from physiotwin4d.test_tools import TestTools -from physiotwin4d.workflow_fit_statistical_model_to_patient import ( +from physiotwin4d import ( + ContourTools, + SegmentChestTotalSegmentator, + # SegmentHeartSimplewareTrimmedBranches, + # SegmentChestTotalSegmentatorWithContrast, + TestTools, WorkflowFitStatisticalModelToPatient, ) +# %% +# Data directory specification + # nnUNetv2 (used by TotalSegmentator inside several workflows) spawns a # multiprocessing.Pool. On Windows the spawn start method re-imports this # script in each child; without the __name__ == "__main__" guard around # top-level work, that re-import fires the segmenter again and Python's # spawn-cascade detector raises RuntimeError. if __name__ == "__main__": - # %% - # Data directory specification REPO_ROOT = Path(__file__).resolve().parent.parent TUTORIALS_DIR = Path(__file__).resolve().parent - DATA_DIR = REPO_ROOT / "data" - FULL_DATA_DIR = DATA_DIR / "KCL-Heart-Model" - TEST_DATA_DIR = DATA_DIR / "test" / "KCL-Heart-Model" + DATA_DIR = REPO_ROOT / "data" / "DirLab-4DCT" OUTPUT_DIR = TUTORIALS_DIR / "output" / "tutorial_04" BASELINES_DIR = REPO_ROOT / "tests" / "baselines" PCA_JSON = TUTORIALS_DIR / "output" / "tutorial_03" / "pca_model.json" + PCA_MEAN_FILE = TUTORIALS_DIR / "output" / "tutorial_03" / "pca_mean_surface.vtp" + # .mha files are DirLab-4DCT data already converted to HU by + # data/DirLab-4DCT/fix_downloaded_data.py. + PATIENT_IMAGE_FILE = DATA_DIR / "Case1Pack_T70.mha" + SEGMENTATION_METHOD = SegmentChestTotalSegmentator() + SEGMENTATION_METHOD.set_has_highres_heart_license(True) + # SEGMENTATION_METHOD = SegmentHeartSimplewareTrimmedBranches() # Use when available + # and images are contrast-enhanced. + # SEGMENTATION_METHOD = SegmentChestTotalSegmentatorWithContrast() # Use when contrast-enhanced + # images and Simpleware is not available. LOG_LEVEL = logging.INFO - # %% - # Data reading - test_mode = TestTools.running_as_test() + TEST_MODE = TestTools.running_as_test() - data_dir = TEST_DATA_DIR if test_mode else FULL_DATA_DIR + data_dir = DATA_DIR output_dir = OUTPUT_DIR pca_json = PCA_JSON + pca_mean_file = PCA_MEAN_FILE + patient_image_file = PATIENT_IMAGE_FILE + segmentation_method = SEGMENTATION_METHOD log_level = LOG_LEVEL output_dir.mkdir(parents=True, exist_ok=True) - template_file = data_dir / "pca_mean.vtu" - if not template_file.exists(): + if not pca_mean_file.exists(): raise FileNotFoundError( - f"KCL-Heart-Model template not found: {template_file}\n" + f"DirLab-4DCT template not found: {pca_mean_file}\n" "See data/README.md for download instructions." ) + pca_mean = cast(pv.DataSet, pv.read(str(pca_mean_file))) - template_data = cast(pv.DataSet, pv.read(str(template_file))) - if isinstance(template_data, pv.PolyData): - template_model = template_data - else: - template_model = template_data.extract_surface(algorithm="dataset_surface") - - sample_files = sorted((data_dir / "sample_meshes").glob("*.vtu")) - if not sample_files: - sample_files = sorted(data_dir.glob("*.vtu")) - sample_files = [path for path in sample_files if path.name != "pca_mean.vtu"] - sample_files = sample_files[:3] - if not sample_files: + pca_model: Optional[dict[str, Any]] = None + if pca_json.exists(): + with pca_json.open(encoding="utf-8") as f: + pca_model = json.load(f) + + if not patient_image_file.exists(): raise FileNotFoundError( - f"No patient-like sample meshes found under {data_dir}.\n" + f"DirLab-4DCT template not found: {patient_image_file}\n" "See data/README.md for download instructions." ) + patient_image = itk.imread(str(patient_image_file)) + itk.imwrite(patient_image, output_dir / "patient_image.nii.gz") - patient_models: list[pv.DataSet] = [] - for sample_file in sample_files: - sample_data = cast(pv.DataSet, pv.read(str(sample_file))) - if isinstance(sample_data, pv.PolyData): - patient_models.append(sample_data) - else: - patient_models.append( - sample_data.extract_surface(algorithm="dataset_surface") - ) + segmentation_result = segmentation_method.segment( + patient_image, + ) + patient_labelmap = segmentation_result["labelmap"] + itk.imwrite(patient_labelmap, output_dir / "patient_labelmap.nii.gz") - pca_model: Optional[dict[str, Any]] = None - if pca_json.exists(): - with pca_json.open(encoding="utf-8") as f: - pca_model = json.load(f) + heart_labelmap = segmentation_result["heart"] + itk.imwrite(heart_labelmap, output_dir / "heart_labelmap.nii.gz") + + contour_tools = ContourTools() + heart_surface = contour_tools.extract_contours( + labelmap_image=heart_labelmap, + ) + heart_surface.save(output_dir / "heart_surface.vtp") # %% # Workflow initialization workflow = WorkflowFitStatisticalModelToPatient( - template_model=template_model, - patient_models=patient_models, + template_model=pca_mean, + patient_models=[heart_surface], + patient_image=patient_image, + patient_labelmap=heart_labelmap, log_level=log_level, + labelmap_interior_object_ids=[141, 142, 143, 144], ) if pca_model is not None: - workflow.set_use_pca_registration(True, pca_model=pca_model) + workflow.set_use_pca_registration( + use_pca_registration=True, + pca_model=pca_model, + use_surface=False, + ) # %% # Workflow execution - result = workflow.run_workflow() + workflow_results = workflow.process() # %% # Result saving - registered_surface: pv.PolyData = result["registered_template_model_surface"] - registered_file = output_dir / "registered_template.vtp" - registered_surface.save(str(registered_file)) + registered_coefficients = workflow.pca_coefficients + if registered_coefficients is not None: + registered_coefficients_path = output_dir / "registered_coefficients.json" + with registered_coefficients_path.open(mode="w", encoding="utf-8") as f: + json.dump(registered_coefficients.tolist(), f) - patient_combined = ( - pv.merge(patient_models) if len(patient_models) > 1 else patient_models[0] - ) - patient_surface = cast(pv.PolyData, patient_combined) + template_mesh = workflow.pca_template_model + template_mesh.save(str(output_dir / "template_mesh.vtp")) + + template_surface = workflow.pca_template_model_surface + template_surface.save(str(output_dir / "template_surface.vtp")) + + registered_mesh = workflow_results["registered_template_model"] + registered_mesh.save(str(output_dir / "template_mesh_registered.vtp")) + registered_surface = workflow_results["registered_template_model_surface"] + registered_surface.save(str(output_dir / "template_surface_registered.vtp")) + + # %% try: pv.start_xvfb() except Exception: @@ -131,8 +159,8 @@ before_path = output_dir / "model_before_registration.png" plotter = pv.Plotter(off_screen=True, window_size=[800, 600]) - plotter.add_mesh(template_model, color="dodgerblue", opacity=0.6) - plotter.add_mesh(patient_surface, color="tomato", opacity=0.6) + plotter.add_mesh(pca_mean, color="dodgerblue", opacity=0.6) + plotter.add_mesh(heart_surface, color="tomato", opacity=0.6) plotter.camera_position = "iso" plotter.screenshot(str(before_path)) plotter.close() @@ -141,7 +169,7 @@ after_path = output_dir / "model_after_registration.png" plotter = pv.Plotter(off_screen=True, window_size=[800, 600]) plotter.add_mesh(registered_surface, color="limegreen", opacity=0.7) - plotter.add_mesh(patient_surface, color="tomato", opacity=0.4) + plotter.add_mesh(heart_surface, color="tomato", opacity=0.4) plotter.camera_position = "iso" plotter.screenshot(str(after_path)) plotter.close() @@ -155,7 +183,7 @@ ) tutorial_results = { - "registered_model": registered_surface, - "registered_file": registered_file, + "registered_mesh": registered_mesh, + "registered_surface": registered_surface, "screenshots": screenshots, } diff --git a/tutorials/tutorial_06_reconstruct_highres_4d_ct.py b/tutorials/tutorial_06_reconstruct_highres_4d_ct.py index 1d88f65..8f33444 100644 --- a/tutorials/tutorial_06_reconstruct_highres_4d_ct.py +++ b/tutorials/tutorial_06_reconstruct_highres_4d_ct.py @@ -41,11 +41,13 @@ REPO_ROOT = Path(__file__).resolve().parent.parent TUTORIALS_DIR = Path(__file__).resolve().parent DATA_DIR = REPO_ROOT / "data" - FULL_DATA_DIR = DATA_DIR / "DirLab-4DCT" / "Case1" - TEST_DATA_DIR = DATA_DIR / "test" / "DirLab-4DCT" / "Case1" + FULL_DATA_DIR = DATA_DIR / "DirLab-4DCT" + TEST_DATA_DIR = DATA_DIR / "test" / "DirLab-4DCT" + # .mha files are DirLab-4DCT data already converted to HU by + # data/DirLab-4DCT/fix_downloaded_data.py. + CASE_GLOB = "Case1Pack_T??.mha" OUTPUT_DIR = TUTORIALS_DIR / "output" / "tutorial_06" BASELINES_DIR = REPO_ROOT / "tests" / "baselines" - MAX_FRAMES = 4 LOG_LEVEL = logging.INFO # %% @@ -57,22 +59,19 @@ log_level = LOG_LEVEL if test_mode: - max_frames = min(MAX_FRAMES, 3) number_of_iterations_Greedy = [1, 0] else: - max_frames = MAX_FRAMES number_of_iterations_Greedy = [30, 15, 7, 3] output_dir.mkdir(parents=True, exist_ok=True) - phase_files = sorted(list(data_dir.glob("*.mhd")) + list(data_dir.glob("*.mha"))) + phase_files = sorted(list(data_dir.glob(CASE_GLOB))) if not phase_files: raise FileNotFoundError( f"No DirLab phase images found under {data_dir}.\n" "See data/README.md for download instructions." ) - phase_files = phase_files[:max_frames] time_series = [itk.imread(str(path)) for path in phase_files] fixed_image = time_series[0] @@ -83,7 +82,7 @@ workflow = WorkflowReconstructHighres4DCT( time_series_images=time_series, fixed_image=fixed_image, - reference_frame=0, + reference_frame=6, registration_method=registration_method, log_level=log_level, ) @@ -91,10 +90,12 @@ # %% # Workflow execution - result = workflow.run_workflow() + result = workflow.process() # %% # Result saving + forward_transform = result["forward_transforms"] + inverse_transform = result["inverse_transforms"] reconstructed_images: list[itk.Image] = result["reconstructed_images"] reconstructed_files: list[Path] = [] for frame_index, image in enumerate(reconstructed_images): @@ -102,6 +103,12 @@ itk.imwrite(image, str(out_path), compression=True) reconstructed_files.append(out_path) + out_path = output_dir / f"reconstructed_frame_{frame_index:03d}_fwd.hdf" + itk.transformwrite(forward_transform[frame_index], str(out_path)) + + out_path = output_dir / f"reconstructed_frame_{frame_index:03d}_inv.hdf" + itk.transformwrite(inverse_transform[frame_index], str(out_path)) + tt = TestTools( class_name="tutorial_06_reconstruct_highres_4d_ct", results_dir=output_dir, diff --git a/tutorials/tutorial_08_cardiac_fit_model.py b/tutorials/tutorial_08_cardiac_fit_model.py index 4732265..c634a00 100644 --- a/tutorials/tutorial_08_cardiac_fit_model.py +++ b/tutorials/tutorial_08_cardiac_fit_model.py @@ -153,14 +153,14 @@ ssm_fit_workflow.set_use_pca_registration( use_pca_registration=True, pca_model=ssm_model, - pca_uses_surface=False, + use_surface=False, ) - ssm_fit_workflow_result = ssm_fit_workflow.run_workflow() + ssm_fit_workflow_result = ssm_fit_workflow.process() ssm_pca_coefficients = ssm_fit_workflow.pca_coefficients assert ssm_pca_coefficients is not None, ( - "pca_coefficients must be set after run_workflow() with " + "pca_coefficients must be set after process() with " "use_pca_registration=True" ) with ssm_pca_coefficients_path.open(mode="w", encoding="utf-8") as f: @@ -214,7 +214,7 @@ registration_method=icon_registration_method, ) reg_workflow.set_modality("ct") - reg_result = reg_workflow.run_workflow() + reg_result = reg_workflow.process() reconstructed_images = reg_result["reconstructed_images"] else: