{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":117682,"databundleVersionId":15062069,"isSourceIdPinned":false,"sourceType":"competition"},{"sourceId":14066012,"sourceType":"datasetVersion","datasetId":8953198},{"sourceId":14330461,"sourceType":"datasetVersion","datasetId":8944423},{"sourceId":14406731,"sourceType":"datasetVersion","datasetId":9166727},{"sourceId":14590716,"sourceType":"datasetVersion","datasetId":9319946},{"sourceId":14591140,"sourceType":"datasetVersion","datasetId":9320245},{"sourceId":14631335,"sourceType":"datasetVersion","datasetId":9321909},{"sourceId":14634817,"sourceType":"datasetVersion","datasetId":9323860},{"sourceId":14641327,"sourceType":"datasetVersion","datasetId":9353000}],"dockerImageVersionId":31193,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Vesuvius Challenge - Surface Detection with nnUNet\n\n**Purpose:** 3D semantic segmentation baseline using nnUNetv2 for detecting papyrus surfaces\nin CT scan volumes from the Vesuvius Challenge.\n\n**Competition:** https://www.kaggle.com/competitions/vesuvius-challenge-surface-detection\n\n### Problem Overview\n- Input: 3D CT volumes of ancient papyrus scrolls (.tif files)\n- Output: Binary segmentation mask identifying papyrus surface pixels\n- Challenge: Partially labeled data (some regions have no ground truth)\n\n### Solution Approach\nUses nnUNet's built-in **ignore label** mechanism:\n- Label 0 = background\n- Label 1 = papyrus surface  \n- Label 2 = unlabeled/ignore (excluded from loss computation)\n\nnnUNet automatically handles partial labels using **partial losses** - only computing\ngradients on pixels with known labels (0 or 1), while ignoring label 2 regions.\n\n### Key Features of This Baseline\n1. **Native TIFF support** - Custom SimpleTiffIO reader, no NIfTI conversion needed\n2. **Symlink-based data prep** - Fast setup without copying large files\n3. **Pre-processed data caching** - Upload to Kaggle dataset to skip preprocessing\n4. **Multi-GPU support** - DDP training with auto-detection\n5. **Configurable epochs** - Use built-in trainer variants (1-8000 epochs)\n6. **Visualization tools** - Training progress and prediction visualization\n\n### File Structure\n```\n/kaggle/input/vesuvius-challenge-surface-detection/\n├── train_images/*.tif      # Training CT volumes\n├── train_labels/*.tif      # Training labels (0=bg, 1=surface, 2=ignore)\n└── test_images/*.tif       # Test CT volumes (no labels)\n\n/kaggle/temp/nnUNet_data/    # Working directory (cleared between sessions)\n├── nnUNet_raw/              # Symlinks to competition data\n└── nnUNet_preprocessed/     # Preprocessed data (can cache as Kaggle dataset)\n\n/kaggle/working/             # Persistent outputs\n├── nnUNet_results/          # Trained models and checkpoints\n├── predictions_tiff/        # Inference outputs\n└── submission.zip           # Final submission\n```\n\n### Usage\n1. **Quick start:** Just run `full_pipeline()` with default settings\n2. **Fast experiments:** Use `epochs=50` or `config=\"2d\"`\n3. **Production:** Use `epochs=250-500` with `3d_fullres` or `3d_lowres`\n4. **Resume training:** Use `continue_training=True`\n\n### Performance Notes\n- Preprocessing: ~1-2 hours (cache as Kaggle dataset to skip)\n- Training: ~5-10 min/epoch on T4, ~2-3 min/epoch on better GPUs\n- fold=\"all\" validation: Runs on ALL training data after training (can be slow)\n- Inference: ~1-2 min per volume\n\n### Conversation Checkpoint\nThis code serves as a self-contained baseline. Key decisions made:\n- Using nnUNetPlannerResEncM (ResNet encoder, medium) - good speed/quality balance\n- Using 3d_fullres config - best quality single-stage segmentation\n- Using fold=\"all\" - trains on all data, no cross-validation (faster for submission)\n- Multi-GPU enabled by default via auto-detection\n- model_path auto-parses epochs/plans/config/fold from path structure\n\n---\n\n### Quick Command Reference\n\n```bash\n# Preprocessing\nnnUNetv2_plan_and_preprocess -d 100 -np 4 -pl nnUNetPlannerResEncM\n\n# Training (uses ignore label automatically)\nnnUNetv2_train 100 3d_fullres 0 -p nnUNetResEncUNetMPlans\n\n# Training with fewer epochs (faster)\nnnUNetv2_train 100 3d_fullres 0 -p nnUNetResEncUNetMPlans -tr nnUNetTrainer_250epochs\n\n# Training all folds in parallel (different GPUs)\nCUDA_VISIBLE_DEVICES=0 nnUNetv2_train 100 3d_fullres 0 ... &\nCUDA_VISIBLE_DEVICES=1 nnUNetv2_train 100 3d_fullres 1 ... &\n\n# Multi-GPU DDP training (single fold, multiple GPUs)\nnnUNetv2_train 100 3d_fullres 0 ... -num_gpus 2\n\n# Inference\nnnUNetv2_predict -d 100 -c 3d_fullres -f 0 -i INPUT -o OUTPUT -p nnUNetResEncUNetMPlans\n```\n\n### Configuration Reference\n\n#### Configurations\n| Config | Use Case | Speed | Quality |\n|--------|----------|-------|---------|\n| `2d` | Fast experiments, slice-by-slice | Fastest | Lower |\n| `3d_lowres` | Large volumes, good balance | Fast | Good |\n| `3d_fullres` | Best quality, smaller volumes | Slow | Best |\n| `3d_cascade_fullres` | Very large volumes (2-stage) | Slowest | Best |\n\n#### Planners\n| Planner | Plans Name | Notes |\n|---------|------------|-------|\n| `nnUNetPlanner` | `nnUNetPlans` | Default U-Net |\n| `nnUNetPlannerResEncM` | `nnUNetResEncUNetMPlans` | ResNet encoder, medium (recommended) |\n| `nnUNetPlannerResEncL` | `nnUNetResEncUNetLPlans` | ResNet encoder, large |\n\n#### Available Epochs (built-in trainer variants)\n```\n1, 5, 10, 20, 50, 100, 250, 500, 750, 1000, 2000, 4000, 8000\n```\n- `epochs=None` or `1000` = default full training\n- `epochs=50-100` = quick experiments\n- `epochs=250-500` = good balance for most datasets\n\n### Hyperparameter Tuning Guide\n\n**nnUNet auto-configures most parameters.** Only tune these if needed:\n\n#### What nnUNet Auto-Configures (DON'T TOUCH):\n- Patch size (based on GPU memory and data)\n- Batch size (based on GPU memory)\n- Network architecture depth\n- Normalization scheme\n- Data augmentation\n- Learning rate schedule (poly LR: `lr = initial_lr * (1 - epoch/max_epochs)^0.9`)\n\n#### What You CAN Tune:\n| Parameter | Default | How to Change | When to Change |\n|-----------|---------|---------------|----------------|\n| **Epochs** | 1000 | `epochs=250` | If converged early (check progress.png) |\n| **Planner** | nnUNetPlanner | `planner=\"nnUNetPlannerResEncM\"` | For potentially better accuracy |\n| **Folds** | 5-fold CV | `fold=\"all\"` | If limited time/compute |\n| **Config** | 3d_fullres | `config=\"2d\"` or `\"3d_lowres\"` | For faster experiments |\n| **GPUs** | auto | `num_gpus=2` | For multi-GPU training |\n\n#### Tuning Recommendations:\n1. **Start with defaults** - nnUNet is well-tuned out of the box\n2. **Try ResEnc planner** - often 1-2% better than default\n3. **Use 250-500 epochs** if model converges early (check progress.png)\n4. **Use 2d config** for quick experiments before 3d_fullres\n5. **Don't tune** batch size, patch size, LR - nnUNet handles this\n\n### References\n- [nnUNet GitHub](https://github.com/MIC-DKFZ/nnUNet)\n- [nnUNet Ignore Label Docs](https://github.com/MIC-DKFZ/nnUNet/blob/master/documentation/ignore_label.md)\n- [Vesuvius Challenge](https://scrollprize.org/)\n- Isensee et al. (2021) - nnU-Net: Nature Methods\n- Gotkowski et al. - Partial annotations paper (for ignore label support)","metadata":{"_uuid":"5fb700b5-f434-4bca-8f07-763a559c66d0","_cell_guid":"193fd083-4485-416d-a30c-334e721e02f9","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"markdown","source":"## 1. Configuration & Imports\n\n**Section Summary:**\nThis section defines all configurable parameters for the pipeline. Modify these\nto adapt the baseline to your needs.\n\n**Key Configuration Decisions:**\n- `FOLD=\"all\"`: Trains on all data without cross-validation (faster for submissions)\n- `CONFIGURATION=\"3d_lowres\"`: Good balance of speed and quality for large volumes\n- `NUM_GPUS=auto`: Uses all available GPUs for DDP training\n- `EPOCHS=None`: Default 1000 epochs (can reduce to 250-500 if converged)\n\n**Preprocessed Data Caching:**\nTo skip the 1-2 hour preprocessing step:\n1. Run preprocessing once\n2. Upload `nnUNet_preprocessed` folder as Kaggle dataset\n3. Set `PREPARED_PREPROCESSED_PATH` to that dataset path\n\n**Path Structure:**\n- INPUT_DIR: Competition data (read-only)\n- WORKING_DIR: Temporary files (/kaggle/temp, cleared between sessions)\n- OUTPUT_DIR: Persistent outputs (/kaggle/working)","metadata":{"_uuid":"25f7ae89-ae88-40df-a51a-660b3fdb8a2b","_cell_guid":"18f34627-74d1-4dc9-b072-760fe5a19cec","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"import os\nimport json\nimport shutil\nimport subprocess\nfrom functools import partial\nfrom multiprocessing import Pool\nfrom pathlib import Path\nfrom typing import Optional, Tuple, List, Literal, Union\n\n# =============================================================================\n# TYPE DEFINITIONS\n# =============================================================================\n\n# Available epoch counts from nnUNet's built-in trainers\n# These are pre-defined trainer classes in nnUNet - you cannot use arbitrary values\n# See: nnunetv2/training/nnUNetTrainer/variants/training_length/nnUNetTrainer_Xepochs.py\nEpochs = Literal[1, 5, 10, 20, 50, 100, 250, 500, 750, 1000, 2000, 4000, 8000]\n\n# =============================================================================\n# CONFIGURATION - MODIFY THESE FOR YOUR SETUP\n# =============================================================================\n\n# Competition data path (Kaggle dataset)\nINPUT_DIR = Path(\"/kaggle/input/vesuvius-challenge-surface-detection\")\n\n# Pre-prepared nnUNet preprocessed dataset\n# Upload your preprocessed data as a Kaggle dataset to skip the 1-2 hour preprocessing step\n# Set to non-existent path to force fresh preprocessing\nPREPARED_PREPROCESSED_PATH = Path(\"/kaggle/input/d/choudharymanas/nn-unet-preprocessed-192/nnUNet_preprocessed/Dataset100_VesuviusSurface\")\n\n# Working directories\nWORKING_DIR = Path(\"/kaggle/temp\")  # Large intermediate files (cleared between sessions)\nOUTPUT_DIR = Path(\"/kaggle/working\")  # Final outputs (persisted)\n\n# nnUNet directory structure (follows nnUNet conventions)\nNNUNET_BASE = WORKING_DIR / \"nnUNet_data\"\nNNUNET_RAW = NNUNET_BASE / \"nnUNet_raw\"  # Small - uses symlinks to competition data\nNNUNET_PREPROCESSED = NNUNET_BASE / \"nnUNet_preprocessed\"  # Large - can use pre-prepared dataset\nNNUNET_RESULTS = OUTPUT_DIR / \"nnUNet_results\"  # Trained models go here (persisted)\n\n# Dataset configuration\nDATASET_ID = 100  # nnUNet dataset ID (arbitrary, just needs to be consistent)\nDATASET_NAME = f\"Dataset{DATASET_ID:03d}_VesuviusSurface\"\n\n# =============================================================================\n# TRAINING CONFIGURATION\n# =============================================================================\n\n# Fold: 0-4 for 5-fold cross-validation, \"all\" for training on all data\n# - Use \"all\" for final submission (faster, uses all training data)\n# - Use 0-4 for model validation and ensembling\nFOLD: Union[int, str] = \"all\"\n\n# Configuration: determines network architecture and resolution\n# - \"2d\": Fastest, processes slice-by-slice (good for quick experiments)\n# - \"3d_lowres\": Fast, handles large volumes well (NOTE: creates predicted_next_stage for cascade)\n# - \"3d_fullres\": Best quality, single-stage (RECOMMENDED)\n# - \"3d_cascade_fullres\": Two-stage for very large volumes (slowest, uses 3d_lowres predictions)\nCONFIGURATION = \"3d_fullres\"\n\n# Planner: determines network architecture variant\n# - \"nnUNetPlanner\": Standard U-Net encoder\n# - \"nnUNetPlannerResEncM\": ResNet encoder, medium size (RECOMMENDED - often 1-2% better)\n# - \"nnUNetPlannerResEncL\": ResNet encoder, large (more parameters, slower)\nPLANNER = \"nnUNetPlannerResEncM\"\nPLANS_NAME = \"nnUNetResEncUNetMPlans\"  # Must match planner (see reference table in header)\n\n# Number of CPU workers for preprocessing and data preparation\nNUM_WORKERS = os.cpu_count() or 4\n\n# Epochs: Number of training epochs\n# - None or 1000: Default full training\n# - 50-100: Quick experiments\n# - 250-500: Good balance (check progress.png to see if converged)\n# Must be one of: 1, 5, 10, 20, 50, 100, 250, 500, 750, 1000, 2000, 4000, 8000\nEPOCHS: Optional[Epochs] = None\n\n# Command timeout in seconds (None for no timeout)\n# Useful for Kaggle's 9-hour limit - set to e.g., 28800 (8 hours) to leave buffer\nCOMMAND_TIMEOUT: Optional[int] = None\n\n\ndef _get_gpu_count() -> int:\n    \"\"\"Get number of available CUDA GPUs.\"\"\"\n    try:\n        import torch\n        return torch.cuda.device_count() if torch.cuda.is_available() else 1\n    except ImportError:\n        return 0\n\n\n# Number of GPUs for DDP training\n# Auto-detected by default. Set to 1 to disable multi-GPU.\n# NOTE: Multi-GPU DDP can sometimes hang in notebook environments.\n# If training hangs, try num_gpus=1\nNUM_GPUS: int = _get_gpu_count()\n\n\n# =============================================================================\n# PATH HELPER FUNCTIONS\n# =============================================================================\n\ndef _get_trainer_name_simple(epochs: Optional[int]) -> str:\n    \"\"\"Get trainer class name based on epochs (simple version for path construction).\"\"\"\n    if epochs is None or epochs == 1000:\n        return \"nnUNetTrainer\"\n    elif epochs == 1:\n        return \"nnUNetTrainer_1epoch\"  # Special case: singular form\n    else:\n        return f\"nnUNetTrainer_{epochs}epochs\"\n\n\ndef get_training_output_dir(\n    epochs: Optional[Epochs] = None,\n    plans: str = PLANS_NAME,\n    config: str = CONFIGURATION,\n    fold: Union[int, str] = FOLD\n) -> Path:\n    \"\"\"\n    Get the training output directory path based on configuration.\n    \n    nnUNet creates this folder structure:\n    NNUNET_RESULTS/Dataset100_VesuviusSurface/nnUNetTrainer__nnUNetResEncUNetMPlans__3d_lowres/fold_all/\n    \n    Use this to find checkpoints, logs, and progress.png\n    \"\"\"\n    _epochs = epochs if epochs is not None else EPOCHS\n    trainer = _get_trainer_name_simple(_epochs)\n    return NNUNET_RESULTS / DATASET_NAME / f\"{trainer}__{plans}__{config}\" / f\"fold_{fold}\"\n\n\ndef get_progress_image_path(\n    epochs: Optional[Epochs] = None,\n    plans: str = PLANS_NAME,\n    config: str = CONFIGURATION,\n    fold: Union[int, str] = FOLD\n) -> Path:\n    \"\"\"Get path to training progress image (loss curves, metrics over epochs).\"\"\"\n    return get_training_output_dir(epochs, plans, config, fold) / \"progress.png\"","metadata":{"_uuid":"62c06cdf-6d1d-4bd7-a0d9-42742bec8f93","_cell_guid":"7e0c318c-855b-4a9f-99e0-9f43eaf31d17","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 2. Environment Setup\n\n**Section Summary:**\nSets up nnUNet environment variables and directory structure.\n\n**nnUNet Environment Variables:**\n- `nnUNet_raw`: Where nnUNet looks for raw dataset\n- `nnUNet_preprocessed`: Where preprocessed data is stored\n- `nnUNet_results`: Where trained models are saved\n- `nnUNet_compile`: Disable torch.compile (can cause issues)\n\n**Pre-prepared Data Handling:**\nThe `_link_prepared_preprocessed()` function handles cached preprocessed data:\n- Copies metadata files (JSON, PKL) that nnUNet may need to modify\n- Symlinks heavy data files (NPZ, B2ND) to save space\n- Handles both direct dataset paths and parent folder paths","metadata":{"_uuid":"aa6d8ff2-81c6-49c2-848b-c3b156e21b72","_cell_guid":"ea24161f-4aee-4d0e-bd71-40401f32b400","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def setup_environment():\n    \"\"\"Set up nnUNet environment variables and directories.\"\"\"\n    for d in [NNUNET_RAW, NNUNET_PREPROCESSED, NNUNET_RESULTS, OUTPUT_DIR]:\n        d.mkdir(parents=True, exist_ok=True)\n    \n    os.environ[\"nnUNet_raw\"] = str(NNUNET_RAW)\n    os.environ[\"nnUNet_preprocessed\"] = str(NNUNET_PREPROCESSED)\n    os.environ[\"nnUNet_results\"] = str(NNUNET_RESULTS)\n    os.environ[\"nnUNet_compile\"] = \"true\"\n    \n    print(f\"nnUNet_raw: {NNUNET_RAW}\")\n    print(f\"nnUNet_preprocessed: {NNUNET_PREPROCESSED}\")\n    print(f\"nnUNet_results: {NNUNET_RESULTS}\")\n    print(f\"nnUNet_USE_BLOSC2: {os.environ.get('nnUNet_USE_BLOSC2', 'not set')} (0=NPZ, 1=blosc2)\")\n    print(f\"NUM_WORKERS: {NUM_WORKERS}\")\n\n\ndef _link_prepared_preprocessed() -> bool:\n    \"\"\"\n    Link/copy pre-prepared preprocessed data if available.\n    \n    Copies the folder structure and metadata files (which nnUNet may need to modify),\n    but symlinks the heavy .npz/.b2nd data files to save space.\n    \n    Handles two possible structures:\n    1. PREPARED_PREPROCESSED_PATH points directly to Dataset100_* folder\n    2. PREPARED_PREPROCESSED_PATH contains Dataset100_* as subfolder\n    \n    Returns True if linked/copied successfully.\n    \"\"\"\n    if not PREPARED_PREPROCESSED_PATH.exists():\n        return False\n    \n    # Determine source directory - could be the path itself or a subfolder\n    source_dir = PREPARED_PREPROCESSED_PATH\n    if not (source_dir / \"dataset.json\").exists():\n        # Look for Dataset folder inside\n        dataset_folders = list(PREPARED_PREPROCESSED_PATH.glob(f\"Dataset*_{DATASET_NAME.split('_')[1]}*\"))\n        if not dataset_folders:\n            dataset_folders = list(PREPARED_PREPROCESSED_PATH.glob(\"Dataset*\"))\n        if dataset_folders:\n            source_dir = dataset_folders[0]\n        else:\n            print(f\"No dataset folder found in {PREPARED_PREPROCESSED_PATH}\")\n            return False\n    \n    target_dir = NNUNET_PREPROCESSED / DATASET_NAME\n    \n    if target_dir.exists():\n        print(f\"Preprocessed data already exists: {target_dir}\")\n        return True\n    \n    print(f\"Linking preprocessed data from: {source_dir}\")\n    target_dir.mkdir(parents=True, exist_ok=True)\n    \n    # Files that nnUNet may need to write/modify - copy these\n    copy_patterns = ['*.json', '*.pkl', '*.txt']\n    \n    # Heavy data files - symlink these  \n    symlink_patterns = ['*.npz', '*.npy', '*.b2nd']\n    \n    copied = 0\n    linked = 0\n    \n    for src_path in source_dir.rglob('*'):\n        if src_path.is_dir():\n            continue\n            \n        # Compute relative path and create target path\n        rel_path = src_path.relative_to(source_dir)\n        dst_path = target_dir / rel_path\n        dst_path.parent.mkdir(parents=True, exist_ok=True)\n        \n        # Check if this is a heavy data file (symlink) or metadata (copy)\n        is_data_file = any(src_path.match(pat) for pat in symlink_patterns)\n        \n        if is_data_file:\n            if not dst_path.exists():\n                dst_path.symlink_to(src_path.resolve())\n                linked += 1\n        else:\n            if not dst_path.exists():\n                shutil.copy2(src_path, dst_path)\n                copied += 1\n    \n    print(f\"Prepared preprocessed data: {copied} files copied, {linked} files symlinked\")\n    print(f\"Location: {target_dir}\")\n    return True\n\n\nsetup_environment()","metadata":{"_uuid":"7a1c6bbf-63b7-4dcc-8de2-58151acb3ed6","_cell_guid":"87bff483-646f-496d-9f7d-9ee6f0de7199","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 3. Installation & Imports\n\n**Section Summary:**\nInstalls required packages and imports libraries.\n\n**Key Dependencies:**\n- `nnunetv2`: The nnUNet framework (includes PyTorch)\n- `nibabel`: For loading NIfTI files (nnUNet's default format)\n- `tifffile`: For reading/writing TIFF files (competition format)\n- `matplotlib`: For visualization\n\n**Important Notes:**\n- Set `nnUNet_USE_BLOSC2=1` BEFORE importing nnunetv2 to use blosc2 format (faster)\n- Blosc2 format can cause issues on some systems","metadata":{"_uuid":"665ab694-31f8-43dd-8780-5e98a42c1bf5","_cell_guid":"6a77e8bd-fee1-4c9f-9b5e-b21af0f52acf","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"!mkdir -p /kaggle/temp\n!mkdir predictions_tiff\n!pip install nnunetv2 nibabel tifffile tqdm -q --no-index -f \"/kaggle/input/surface-packages-offline\"\n\n# IMPORTANT: Set this BEFORE importing nnunetv2\n# Blosc2 is a newer compression format but can cause compatibility issues\nos.environ[\"nnUNet_USE_BLOSC2\"] = \"1\"  # Use blosc2 format (faster, smaller files)\n\nimport nibabel as nib\nimport numpy as np\nimport pandas as pd\nimport tifffile\nimport matplotlib.pyplot as plt\nfrom tqdm.auto import tqdm\n\n# Show GPU configuration\nprint(f\"Available GPUs: {_get_gpu_count()}\")\nprint(f\"Using NUM_GPUS={NUM_GPUS}\")","metadata":{"_uuid":"82a5abad-98d7-4e54-91bf-78afe62ddbb7","_cell_guid":"69b7c4ca-d393-41d0-ae77-e96a9e952a13","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 4. Data Utilities\n\n**Section Summary:**\nHelper functions for handling TIFF and NIfTI file formats.\n\n**TIFF Format:**\n- Competition uses 3D TIFF files for CT volumes\n- nnUNet natively supports TIFF via our custom SimpleTiffIO reader\n- No conversion to NIfTI needed (saves disk space and time)\n\n**Spacing Information:**\n- nnUNet requires voxel spacing (physical size of each voxel)\n- We use isotropic spacing (1.0, 1.0, 1.0) as competition doesn't specify\n- Stored in JSON sidecar files alongside TIFF images","metadata":{"_uuid":"4d2688e3-0654-46d7-9325-87be0c939e27","_cell_guid":"a7785be3-ccdb-4015-83dd-3fe2b43bb7eb","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def create_spacing_json(output_path: Path, shape: tuple, spacing: tuple = (1.0, 1.0, 1.0)):\n    \"\"\"Create JSON sidecar with spacing info for TIFF files.\"\"\"\n    json_data = {\"spacing\": list(spacing)}\n    with open(output_path, \"w\") as f:\n        json.dump(json_data, f)\n\n\ndef load_nifti(path: Path) -> np.ndarray:\n    \"\"\"Load NIfTI file (used for loading nnUNet predictions).\"\"\"\n    return nib.load(str(path)).get_fdata()\n\n\ndef create_dataset_json(output_dir: Path, num_training: int, file_ending: str = \".tif\") -> dict:\n    \"\"\"Create dataset.json with ignore label support and 3D TIFF reader.\"\"\"\n    \n    dataset_json = {\n        \"channel_names\": {\"0\": \"CT\"},\n        \"labels\": {\"background\": 0, \"surface\": 1, \"ignore\": 2},\n        \"numTraining\": num_training,\n        \"file_ending\": file_ending,\n        \"overwrite_image_reader_writer\": \"SimpleTiffIO\"\n    }\n    \n    json_path = output_dir / \"dataset.json\"\n    with open(json_path, \"w\") as f:\n        json.dump(dataset_json, f, indent=4)\n    \n    print(f\"Created {json_path}\")\n    print(f\"  - {num_training} training cases\")\n    print(f\"  - Labels: background(0), surface(1), ignore(2)\")\n    print(f\"  - Reader: SimpleTiffIO (3D TIFF)\")\n    \n    return dataset_json","metadata":{"_uuid":"96f7afc7-cd4d-4304-93b4-dcb7fd1b096b","_cell_guid":"d1c4bcf6-c084-4149-83ee-e5dce0c64f6c","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 5. Dataset Preparation\n\n**Section Summary:**\nConverts competition data to nnUNet's expected format.\n\n**nnUNet Dataset Structure:**\n```\nnnUNet_raw/Dataset100_VesuviusSurface/\n├── imagesTr/           # Training images\n│   ├── case001_0000.tif  # _0000 suffix = channel 0 (only one for CT)\n│   └── case001_0000.json # Spacing information\n├── labelsTr/           # Training labels\n│   ├── case001.tif\n│   └── case001.json\n└── dataset.json        # Dataset configuration\n```\n\n**Symlink Strategy:**\n- Uses symlinks instead of copying (fast, saves disk space)\n- Only JSON sidecar files are created (spacing info)\n- Original TIFF files remain in competition dataset\n\n**Parallel Processing:**\n- Uses multiprocessing for faster preparation\n- ~2 minutes for 806 training cases","metadata":{"_uuid":"9ece55f8-b8fc-457a-9563-95b732f60c1e","_cell_guid":"9cd77ab4-9f29-4203-8193-2b29d9013f0f","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def prepare_single_case(\n    src_path: Path, \n    dest_path: Path, \n    json_path: Path, \n    use_symlinks: bool = True\n) -> bool:\n    \"\"\"\n    Prepare a single TIFF file for nnUNet: create symlink/copy and JSON sidecar.\n    Returns True on success, False on failure.\n    \"\"\"\n    try:\n        # Get shape for JSON\n        with tifffile.TiffFile(src_path) as tif:\n            shape = tif.pages[0].shape if len(tif.pages) == 1 else (len(tif.pages), *tif.pages[0].shape)\n        \n        # Link or copy file\n        if use_symlinks:\n            if not dest_path.exists():\n                dest_path.symlink_to(src_path.resolve())\n        else:\n            shutil.copy2(src_path, dest_path)\n        \n        # Create JSON sidecar\n        create_spacing_json(json_path, shape)\n        return True\n        \n    except Exception as e:\n        print(f\"Error processing {src_path.name}: {e}\")\n        return False\n\n\ndef _prepare_training_case(\n    img_path: Path,\n    train_labels_dir: Path,\n    images_dir: Path,\n    labels_dir: Path,\n    use_symlinks: bool\n) -> bool:\n    \"\"\"Worker function for parallel dataset preparation.\"\"\"\n    case_id = img_path.stem\n    label_path = train_labels_dir / img_path.name\n    \n    if not label_path.exists():\n        return False\n    \n    img_ok = prepare_single_case(\n        img_path,\n        images_dir / f\"{case_id}_0000.tif\",\n        images_dir / f\"{case_id}_0000.json\",\n        use_symlinks\n    )\n    \n    label_ok = prepare_single_case(\n        label_path,\n        labels_dir / f\"{case_id}.tif\",\n        labels_dir / f\"{case_id}.json\",\n        use_symlinks\n    )\n    \n    return img_ok and label_ok\n\n\ndef prepare_dataset(input_dir: Path, max_cases: Optional[int] = None, use_symlinks: bool = True):\n    \"\"\"\n    Convert competition data to nnUNet format using TIFF directly (no NIfTI).\n    Uses multiprocessing for faster preparation.\n    \n    Competition structure:\n    - train_images/*.tif  (3D volumes)\n    - train_labels/*.tif  (3D labels: 0=bg, 1=surface, 2=ignore)\n    \"\"\"\n    dataset_dir = NNUNET_RAW / DATASET_NAME\n    images_dir = dataset_dir / \"imagesTr\"\n    labels_dir = dataset_dir / \"labelsTr\"\n    \n    images_dir.mkdir(parents=True, exist_ok=True)\n    labels_dir.mkdir(parents=True, exist_ok=True)\n    \n    train_images_dir = input_dir / \"train_images\"\n    train_labels_dir = input_dir / \"train_labels\"\n    \n    if not train_images_dir.exists():\n        print(f\"ERROR: {train_images_dir} not found!\")\n        return None\n    \n    image_files = sorted(train_images_dir.glob(\"*.tif\"))\n    if max_cases:\n        image_files = image_files[:max_cases]\n    \n    print(f\"Found {len(image_files)} training cases\")\n    print(f\"Using {'symlinks' if use_symlinks else 'copy'}\")\n    print(f\"Processing with {NUM_WORKERS} workers...\")\n    \n    # Create worker function with fixed arguments\n    worker = partial(\n        _prepare_training_case,\n        train_labels_dir=train_labels_dir,\n        images_dir=images_dir,\n        labels_dir=labels_dir,\n        use_symlinks=use_symlinks\n    )\n    \n    # Process in parallel with progress bar\n    with Pool(NUM_WORKERS) as pool:\n        results = list(tqdm(\n            pool.imap(worker, image_files),\n            total=len(image_files),\n            desc=\"Preparing dataset\"\n        ))\n    \n    num_converted = sum(results)\n    create_dataset_json(dataset_dir, num_converted, file_ending=\".tif\")\n    \n    print(f\"\\nDataset prepared: {num_converted} cases\")\n    print(f\"Location: {dataset_dir}\")\n    \n    return dataset_dir","metadata":{"_uuid":"fad5a235-b611-49cd-b8f2-5519de775dd7","_cell_guid":"0b404165-d4c2-4292-8193-99db616a0e9c","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 6. Synthetic Data (for testing)\n\n**Section Summary:**\nCreates small synthetic datasets for testing the pipeline without real data.\n\n**Use Cases:**\n- Testing pipeline setup before using real data\n- Debugging data loading issues\n- Quick iteration on code changes\n\n**Generated Data:**\n- Small 3D volumes (64x64x64 by default)\n- Spherical shell as \"surface\" (label 1)\n- Random ignore regions (label 2)\n- Random background noise","metadata":{"_uuid":"9d20116b-546f-4bd6-8b08-8c3b1f8488af","_cell_guid":"72871a71-ad2f-431d-ad48-fbdc93ce0874","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def create_single_synthetic_case(\n    case_id: str,\n    size: Tuple[int, int, int],\n    images_dir: Path,\n    labels_dir: Path\n):\n    \"\"\"Create a single synthetic training case.\"\"\"\n    # Create volume with structure\n    volume = np.random.randn(*size).astype(np.float32) * 0.1\n    \n    # Add spherical shell as \"surface\"\n    z, y, x = np.ogrid[:size[0], :size[1], :size[2]]\n    center = np.array(size) // 2\n    dist = np.sqrt((z - center[0])**2 + (y - center[1])**2 + (x - center[2])**2)\n    shell = (dist > 15) & (dist < 20)\n    volume[shell] += 1.0\n    \n    # Create labels\n    labels = np.zeros(size, dtype=np.uint8)\n    labels[shell] = 1  # Surface\n    \n    # Add ignore regions (label 2)\n    ignore_mask = np.random.random(size) < 0.15\n    labels[ignore_mask] = 2\n    \n    # Save as TIFF\n    tifffile.imwrite(images_dir / f\"{case_id}_0000.tif\", volume)\n    tifffile.imwrite(labels_dir / f\"{case_id}.tif\", labels)\n    \n    # Create JSON sidecars\n    create_spacing_json(images_dir / f\"{case_id}_0000.json\", size)\n    create_spacing_json(labels_dir / f\"{case_id}.json\", size)\n\n\ndef create_synthetic_dataset(num_cases: int = 5, size: Tuple[int, int, int] = (64, 64, 64)):\n    \"\"\"Create synthetic 3D data with ignore regions for testing (TIFF format).\"\"\"\n    \n    dataset_dir = NNUNET_RAW / DATASET_NAME\n    images_dir = dataset_dir / \"imagesTr\"\n    labels_dir = dataset_dir / \"labelsTr\"\n    \n    images_dir.mkdir(parents=True, exist_ok=True)\n    labels_dir.mkdir(parents=True, exist_ok=True)\n    \n    for i in tqdm(range(num_cases), desc=\"Creating synthetic data\"):\n        case_id = f\"case_{i:03d}\"\n        create_single_synthetic_case(case_id, size, images_dir, labels_dir)\n    \n    create_dataset_json(dataset_dir, num_cases, file_ending=\".tif\")\n    print(f\"Created synthetic dataset: {dataset_dir}\")\n    \n    return dataset_dir\n\n# To use synthetic data instead of real data, uncomment:\n# create_synthetic_dataset(num_cases=5)","metadata":{"_uuid":"ae7d78c3-e0a1-4259-9ec0-56147aca59e0","_cell_guid":"7fc46883-0414-437d-b4d0-058394c93cec","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 7. nnUNet Commands\n\n**Section Summary:**\nWrapper functions for nnUNet command-line tools.\n\n**nnUNet Pipeline Steps:**\n1. `run_preprocessing()` → `nnUNetv2_plan_and_preprocess`\n   - Analyzes dataset statistics\n   - Creates experiment plans (network architecture, patch size, etc.)\n   - Preprocesses data (resampling, normalization)\n   - Time: 1-2 hours for full dataset\n\n2. `run_training()` → `nnUNetv2_train`\n   - Trains the network\n   - Saves checkpoints every 50 epochs\n   - Generates progress.png with loss curves\n   - Time: 5-10 min/epoch on T4\n\n3. `run_inference()` → `nnUNetv2_predict`\n   - Runs trained model on test data\n   - Outputs NIfTI predictions (converted to TIFF later)\n   - Time: 1-2 min per volume\n\n**Important Training Notes:**\n- `fold=\"all\"`: After training, nnUNet validates on ALL training data\n  This can be slow (50+ min) but cannot be skipped with built-in flags\n- Multi-GPU: Use `-num_gpus X` for DDP training\n- Epochs: Use `-tr nnUNetTrainer_Xepochs` for custom epoch counts","metadata":{"_uuid":"68b6c269-e1c1-4956-9403-a23160df8c23","_cell_guid":"cda8a2f7-aa28-427e-8600-6eb9cddfe4fd","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def _run_command(\n    cmd: str, \n    name: str = \"Command\", \n    tail_lines: int = 20,\n    timeout: Optional[int] = COMMAND_TIMEOUT\n) -> bool:\n    \"\"\"\n    Execute shell command and handle output parsing.\n    \n    Args:\n        cmd: Shell command to execute\n        name: Display name for logging\n        tail_lines: Number of stdout lines to show on success\n        timeout: Timeout in seconds (None for no timeout)\n    \n    Returns:\n        True if command succeeded, False otherwise\n    \"\"\"\n    print(f\"Running: {cmd}\")\n    if timeout:\n        print(f\"Timeout: {timeout}s ({timeout/3600:.1f}h)\")\n    \n    try:\n        result = subprocess.run(\n            cmd, \n            shell=True, \n            capture_output=True, \n            text=True,\n            timeout=timeout\n        )\n    except subprocess.TimeoutExpired:\n        print(f\"{name} TIMEOUT after {timeout}s!\")\n        return False\n    \n    if result.returncode != 0:\n        print(f\"{name} FAILED!\")\n        print(f\"STDERR:\\n{result.stderr[-3000:]}\")\n        return False\n    \n    print(f\"{name} complete!\")\n    if result.stdout.strip():\n        lines = result.stdout.strip().split('\\n')\n        print('\\n'.join(lines[-tail_lines:]))\n    \n    return True\n\n\ndef run_preprocessing(\n    dataset_id: int = DATASET_ID, \n    planner: str = PLANNER,\n    num_workers: int = NUM_WORKERS,\n    configurations: Optional[List[str]] = None,\n    timeout: Optional[int] = COMMAND_TIMEOUT\n) -> bool:\n    \"\"\"\n    Run nnUNet preprocessing.\n    \n    Args:\n        dataset_id: nnUNet dataset ID\n        planner: Planner class name\n        num_workers: Number of CPU workers for parallel processing\n        configurations: List of configs to preprocess (e.g., [\"3d_fullres\"])\n        timeout: Timeout in seconds (None for no timeout)\n    \n    Returns:\n        True if preprocessing succeeded\n    \"\"\"\n    if configurations is None:\n        configurations = [CONFIGURATION]\n    \n    cmd = f\"nnUNetv2_plan_and_preprocess -d {dataset_id:03d} -np {num_workers}\"\n    cmd += f\" -pl {planner}\"\n    cmd += f\" -c {' '.join(configurations)}\"\n    \n    return _run_command(cmd, \"Preprocessing\", timeout=timeout)\n\n\ndef _get_trainer_name(epochs: Optional[Epochs]) -> str:\n    \"\"\"Get trainer class name based on epochs.\"\"\"\n    if epochs is None or epochs == 1000:\n        return \"nnUNetTrainer\"\n    elif epochs == 1:\n        return \"nnUNetTrainer_1epoch\"\n    else:\n        return f\"nnUNetTrainer_{epochs}epochs\"\n\n\ndef run_training(\n    dataset_id: int = DATASET_ID,\n    config: str = CONFIGURATION,\n    fold: Union[int, str] = FOLD,\n    plans: str = PLANS_NAME,\n    epochs: Optional[Epochs] = EPOCHS,\n    pretrained_weights: Optional[Path] = None,\n    continue_training: bool = False,\n    only_run_validation: bool = False,\n    disable_checkpointing: bool = False,\n    npz: bool = False,\n    num_gpus: int = NUM_GPUS,\n    timeout: Optional[int] = COMMAND_TIMEOUT\n) -> bool:\n    \"\"\"\n    Run nnUNet training.\n    \n    Args:\n        dataset_id: nnUNet dataset ID\n        config: Configuration name (3d_fullres, 2d, etc.)\n        fold: Fold number (0-4) or \"all\" for training on all data\n        plans: Plans name matching the planner used\n        epochs: Number of epochs. Available: 1, 5, 10, 20, 50, 100, 150, 200, 250, \n                300, 400, 500, 750, 1000, 2000, 4000, 8000. None = 1000 (default)\n        pretrained_weights: Optional path to pretrained checkpoint for fine-tuning\n        continue_training: Continue from last checkpoint (add -c flag)\n        only_run_validation: Only run validation, skip training\n        disable_checkpointing: Disable saving checkpoints (saves disk space)\n        npz: Save softmax outputs during validation (needed for ensembling)\n        num_gpus: Number of GPUs for DDP training (default: auto-detected)\n        timeout: Timeout in seconds (None for no timeout)\n    \n    Returns:\n        True if training succeeded\n    \n    Note:\n        For multi-GPU (DDP) training, batch size should be divisible by num_gpus.\n        The first run extracts preprocessed data - wait for GPU usage before \n        starting additional folds on other GPUs.\n        \n        When using fold=\"all\", nnUNet will run validation on ALL training data\n        after training completes, which can be slow. This is normal behavior.\n    \"\"\"\n    trainer = _get_trainer_name(epochs)\n    cmd = f\"nnUNetv2_train {dataset_id:03d} {config} {fold} -p {plans} -tr {trainer}\"\n    \n    if pretrained_weights:\n        cmd += f\" -pretrained_weights {pretrained_weights}\"\n    if continue_training:\n        # Find checkpoint to resume from: prefer checkpoint_final.pth, fallback to checkpoint_best.pth\n        model_dir = get_training_output_dir(epochs=epochs, plans=plans, config=config, fold=fold)\n        checkpoint_final = model_dir / \"checkpoint_final.pth\"\n        checkpoint_best = model_dir / \"checkpoint_best.pth\"\n        if checkpoint_final.exists():\n            print(f\"Resuming from: {checkpoint_final}\")\n            cmd += \" --c\"\n        elif checkpoint_best.exists():\n            print(f\"Resuming from: {checkpoint_best}\")\n            cmd += \" --c\"\n        else:\n            print(f\"WARNING: No checkpoint found in {model_dir}, starting fresh\")\n    if only_run_validation:\n        cmd += \" --val\"\n    if disable_checkpointing:\n        cmd += \" --disable_checkpointing\"\n    if npz:\n        cmd += \" --npz\"\n    if num_gpus > 1:\n        cmd += f\" -num_gpus {num_gpus}\"\n    \n    epochs_str = epochs if epochs else 1000\n    gpu_str = f\", {num_gpus} GPUs\" if num_gpus > 1 else \"\"\n    return _run_command(cmd, f\"Training ({epochs_str} epochs{gpu_str})\", tail_lines=30, timeout=timeout)\n\n\ndef run_inference(\n    input_dir: Path,\n    output_dir: Path,\n    dataset_id: int = DATASET_ID,\n    config: str = CONFIGURATION,\n    fold: Union[int, str] = FOLD,\n    plans: str = PLANS_NAME,\n    epochs: Optional[Epochs] = EPOCHS,\n    save_probabilities: bool = True,\n    num_processes_preprocessing: int = 2,\n    num_processes_segmentation: int = 2,\n    timeout: Optional[int] = COMMAND_TIMEOUT\n) -> bool:\n    \"\"\"\n    Run inference with trained model.\n    \n    Args:\n        input_dir: Directory with test images (must have _0000 suffix)\n        output_dir: Directory to save predictions\n        dataset_id: nnUNet dataset ID\n        config: Configuration name\n        fold: Fold number used for training (or \"all\" or tuple like \"0,1,2\")\n        plans: Plans name\n        epochs: Epochs used during training (must match trained model)\n        save_probabilities: Whether to save probability maps (.npz files)\n        num_processes_preprocessing: Parallel processes for preprocessing\n        num_processes_segmentation: Parallel processes for segmentation\n        timeout: Timeout in seconds (None for no timeout)\n    \n    Returns:\n        True if inference succeeded\n    \"\"\"\n    output_dir.mkdir(parents=True, exist_ok=True)\n    \n    trainer = _get_trainer_name(epochs)\n    \n    cmd = f\"nnUNetv2_predict -d {dataset_id:03d} -c {config} -f {fold}\"\n    cmd += f\" -i {input_dir} -o {output_dir} -p {plans} -tr {trainer}\"\n    cmd += f\" -npp {num_processes_preprocessing} -nps {num_processes_segmentation}\"\n    cmd += \" --verbose\"\n    \n    if save_probabilities:\n        cmd += \" --save_probabilities\"\n    \n    return _run_command(cmd, \"Inference\", timeout=timeout)","metadata":{"_uuid":"6c762a1b-4033-46bd-8518-7245c917e103","_cell_guid":"f7372f1e-dcc8-4bee-8320-037609ecd707","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 8. Test Data Preparation\n\n**Section Summary:**\nPrepares test data for nnUNet inference.\n\n**Required Format:**\n- Files must have `_0000` suffix (channel indicator)\n- JSON sidecar files with spacing info\n\n**Note:**\nTest images are in competition data, we just need to rename/link them\nto match nnUNet's expected naming convention.","metadata":{"_uuid":"0f38a3b7-0208-4b1f-9d16-96d6a61c1dc0","_cell_guid":"8d468a99-805e-476f-a555-8490b119ea0d","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def prepare_test_data(input_dir: Path, output_dir: Path, use_symlinks: bool = True) -> Path:\n    \"\"\"Prepare test TIFF images for nnUNet inference.\"\"\"\n    \n    output_dir.mkdir(parents=True, exist_ok=True)\n    \n    test_images_dir = input_dir / \"test_images\"\n    \n    if not test_images_dir.exists():\n        print(f\"ERROR: {test_images_dir} not found!\")\n        return output_dir\n    \n    test_files = sorted(test_images_dir.glob(\"*.tif\"))\n    print(f\"Found {len(test_files)} test cases\")\n    print(f\"Using {'symlinks' if use_symlinks else 'copy'}\")\n    \n    for img_path in tqdm(test_files, desc=\"Preparing test data\"):\n        case_id = img_path.stem\n        prepare_single_case(\n            img_path,\n            output_dir / f\"{case_id}_0000.tif\",\n            output_dir / f\"{case_id}_0000.json\",\n            use_symlinks\n        )\n    \n    return output_dir","metadata":{"_uuid":"a25f68d5-bb4c-42cc-b6af-777552aa5be2","_cell_guid":"f489ef92-b1b0-4dbd-b99c-9a074a90c9e7","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 9. Inference Utilities\n\n**Section Summary:**\nFunctions for converting nnUNet predictions to submission format.\n\n**nnUNet Output Formats:**\nWhen using our SimpleTiffIO reader with `save_probabilities=True`:\n- `.npz` - Probability maps (preferred for conversion)\n- `.tif` - Direct predictions (already in our format)\n- `.pkl` - Metadata\n- `.json` - Configuration info\n\n**Conversion Process:**\n1. Load probability maps from .npz files\n2. Apply argmax to get class predictions\n3. Save as uint8 TIFF (0=background, 1=surface)","metadata":{"_uuid":"9c9ffb9e-d7a5-4b75-9d24-8120da4aa8a6","_cell_guid":"e36848ed-dfb4-4d7f-9137-3a878257ad6d","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"import numpy as np\nimport scipy.ndimage as ndimage\nfrom skimage.morphology import remove_small_objects\n\n\ndef frangi_filter_3d(image, sigmas=(1, 3, 5), beta1=0.5, beta2=15, gamma=None, black_ridges=False):\n    \"\"\"\n    Frangi Filter - Enhances tubular structures in 3D images (e.g., blood vessels).\n    \n    The Frangi filter detects tubular structures based on eigenvalues of the Hessian matrix\n    by analyzing local second-order derivatives to identify regions with \"tubular\" geometric features.\n    \n    Principle:\n    1. Compute the Hessian matrix at different scales using Gaussian smoothing\n    2. Calculate eigenvalues (λ1, λ2, λ3) of the Hessian matrix\n    3. Determine if a region is tubular based on the relationship between eigenvalues\n    4. Use Frangi response function to compute enhancement values\n    \n    Parameters:\n        image: numpy array with shape (D, H, W) or (D, H, W, 1) - 3D image\n        sigmas: tuple of standard deviations for Gaussian kernels used for multi-scale detection\n                e.g., (1, 3, 5) detects tubular structures from 1 to 5 pixels in width\n        beta1: float, controls sensitivity to deviation from tubular structure (default 0.5)\n        beta2: float, controls sensitivity to background (default 15)\n        gamma: float, normalization factor; if None, calculated automatically\n        black_ridges: bool, True detects dark tubular structures (vessels), \n                      False detects bright structures\n    \n    Returns:\n        enhanced: numpy array with same shape as input, enhanced image values in [0, 1]\n    \"\"\"\n    # Ensure input is 3D array\n    if len(image.shape) == 4:\n        image = image[..., 0]\n    \n    # Convert to float type\n    image = image.astype(np.float64)\n    \n    # Invert image if detecting dark structures\n    if black_ridges:\n        image = -image\n    \n    # Initialize output\n    enhanced = np.zeros_like(image)\n    \n    # Compute Frangi response for each scale\n    for sigma in sigmas:\n        # Compute Hessian matrix components using Gaussian derivatives\n        hxx = ndimage.gaussian_filter1d(\n            ndimage.gaussian_filter1d(\n                ndimage.gaussian_filter1d(image, sigma, axis=2, order=2), \n                sigma, axis=1, order=0), \n            sigma, axis=0, order=0)\n        \n        hyy = ndimage.gaussian_filter1d(\n            ndimage.gaussian_filter1d(\n                ndimage.gaussian_filter1d(image, sigma, axis=2, order=0), \n                sigma, axis=1, order=2), \n            sigma, axis=0, order=0)\n        \n        hzz = ndimage.gaussian_filter1d(\n            ndimage.gaussian_filter1d(\n                ndimage.gaussian_filter1d(image, sigma, axis=2, order=0), \n                sigma, axis=1, order=0), \n            sigma, axis=0, order=2)\n        \n        hxy = ndimage.gaussian_filter1d(\n            ndimage.gaussian_filter1d(\n                ndimage.gaussian_filter1d(image, sigma, axis=2, order=1), \n                sigma, axis=1, order=1), \n            sigma, axis=0, order=0)\n        \n        hxz = ndimage.gaussian_filter1d(\n            ndimage.gaussian_filter1d(\n                ndimage.gaussian_filter1d(image, sigma, axis=2, order=1), \n                sigma, axis=1, order=0), \n            sigma, axis=0, order=1)\n        \n        hyz = ndimage.gaussian_filter1d(\n            ndimage.gaussian_filter1d(\n                ndimage.gaussian_filter1d(image, sigma, axis=2, order=0), \n                sigma, axis=1, order=1), \n            sigma, axis=0, order=1)\n        \n        # Simplified Frangi response calculation based on Hessian trace and determinant\n        # For tubular structures: one eigenvalue close to 0, two larger eigenvalues with same sign\n        trace = hxx + hyy + hzz\n        det = (hxx * hyy * hzz + \n               2 * hxy * hxz * hyz - \n               hxx * hyz * hyz - \n               hyy * hxz * hxz - \n               hzz * hxy * hxy)\n        \n        # Normalization factor\n        if gamma is None:\n            gamma = np.max(np.abs(trace)) if np.max(np.abs(trace)) > 0 else 1.0\n        \n        # Compute Frobenius norm of Hessian matrix\n        hessian_norm = np.sqrt(hxx**2 + hyy**2 + hzz**2 + \n                              2 * (hxy**2 + hxz**2 + hyz**2))\n        \n        # Frangi response function\n        response = np.exp(-beta1 * (trace**2) / (hessian_norm + 1e-10)) * \\\n                   (1 - np.exp(-beta2 * hessian_norm**2 / (gamma**2 + 1e-10)))\n        \n        # Keep only positive values (tubular structures)\n        response = np.maximum(response, 0)\n        \n        # Take maximum response across all scales\n        enhanced = np.maximum(enhanced, response)\n    \n    # Normalize to [0, 1]\n    if enhanced.max() > 0:\n        enhanced = enhanced / enhanced.max()\n    \n    return enhanced.astype(np.float32)\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import numpy as np\nfrom pathlib import Path\nfrom PIL import Image\nfrom multiprocessing import Pool\nfrom skimage.morphology import remove_small_objects\nfrom scipy.ndimage import distance_transform_edt\n\n# assumes these already exist:\n# load_tif_volume\n# apply_hysteresis_3d\n# apply_anisotropic_closing_3d\n# apply_meijering_slice_safe\n# frangi_filter_3d\n\n# ==========================================\n# FRANGI PIPELINE HELPER\n# ==========================================\ndef apply_frangi_pipeline(\n    vol,\n   dt_thresh=0.05,\n    frangi_alpha=0.3,\n    final_thresh=0.3,\n    frangi_sigmas=[2],\n    frangi_beta1=0.5,\n    frangi_beta2=15\n):\n     \n\n    \"\"\"\n    Frangi-based processing pipeline:\n    1. Threshold probability volume\n    2. Compute distance transform\n    3. Apply Frangi filter\n    4. Blend original probability with Frangi enhancement\n    5. Apply threshold\n    \n    Parameters:\n        vol: probability volume (D, H, W)\n        dt_thresh: threshold for initial binarization before distance transform\n        frangi_alpha: weight for Frangi in blend (combined = alpha*frangi + (1-alpha)*prob)\n        final_thresh: final binarization threshold after blending\n        frangi_sigmas: Frangi filter scale parameters\n        frangi_beta1, frangi_beta2: Frangi filter parameters\n    \n    Returns:\n        mask3d: processed binary mask\n    \"\"\"\n    # Step 1: Initial threshold to create binary mask\n    binary_mask = vol > dt_thresh\n    \n    # Step 2: Compute distance transform from the binary mask\n    dt = distance_transform_edt(binary_mask)\n    \n    # Step 3: Apply Frangi filter to distance transform\n    frangi_enhanced = frangi_filter_3d(\n        dt,\n        sigmas=frangi_sigmas,\n        beta1=frangi_beta1,\n        beta2=frangi_beta2,\n        black_ridges=True\n    )\n    \n    # Step 4: Blend original probability with Frangi enhancement\n    blended = frangi_alpha * frangi_enhanced + (1 - frangi_alpha) * vol\n    \n    # Step 5: Apply final threshold\n    mask3d = blended > final_thresh\n    \n    return mask3d\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"######### HELPER FUNCTIONS \n\nimport tifffile\nimport scipy.ndimage as ndi\nfrom skimage.morphology import remove_small_objects\n\n\ndef load_volume(path):\n    vol = tifffile.imread(path)          # (D, H, W)\n    vol = vol.astype(np.float32)\n    vol = vol[None, ..., None]           # (1, D, H, W, 1)\n    return vol\n\n\n\n\n# ==========================================\n# HELPER: Anisotropic Structure Builder\n# ==========================================\ndef build_anisotropic_struct(z_radius: int, xy_radius: int):\n    z, r = z_radius, xy_radius\n\n    if z == 0 and r == 0:\n        return None\n\n    if z == 0 and r > 0:\n        size = 2 * r + 1\n        struct = np.zeros((1, size, size), dtype=bool)\n        cy, cx = r, r\n        for dy in range(-r, r + 1):\n            for dx in range(-r, r + 1):\n                if dy * dy + dx * dx <= r * r:\n                    struct[0, cy + dy, cx + dx] = True\n        return struct\n\n    if z > 0 and r == 0:\n        struct = np.zeros((2 * z + 1, 1, 1), dtype=bool)\n        struct[:, 0, 0] = True\n        return struct\n\n    depth = 2 * z + 1\n    size = 2 * r + 1\n    struct = np.zeros((depth, size, size), dtype=bool)\n    cz, cy, cx = z, r, r\n    for dz in range(-z, z + 1):\n        for dy in range(-r, r + 1):\n            for dx in range(-r, r + 1):\n                if dy * dy + dx * dx <= r * r:\n                    struct[cz + dz, cy + dy, cx + dx] = True\n    return struct\n\n### Fill 3D SHEETS\nimport numpy as np\nfrom skimage.segmentation import watershed\nfrom skimage.measure import label\nfrom scipy import ndimage as ndi\nfrom scipy.ndimage import binary_dilation\n\ndef fill_3d_sheets(prob_volume, seed_volume, prob_threshold=0.10):\n    \"\"\"\n    Fills holes in 3D sheets using a probability map and seeds.\n    \n    Parameters:\n    - prob_volume: 3D numpy array (float 0.0 to 1.0). High values = sheet.\n    - seed_volume: 3D numpy array (binary or labeled). Your \"Hysteresis\" output.\n    - prob_threshold: Only allow sheets to grow where probability is above this.\n                      Prevents sheets from expanding into \"dead air\" between layers.\n    \"\"\"\n    \n    # 1. Ensure seeds are labeled\n    # If your seed_volume is binary (0 and 1), we need to give each \n    # disconnected segment a unique ID so the watershed knows they are distinct.\n    if len(np.unique(seed_volume)) <= 2:\n        print(\"Labeling seeds...\")\n        labeled_seeds = label(seed_volume)\n    else:\n        labeled_seeds = seed_volume\n\n\n    # 2. Create the \"Basin\" (The Inverse Probability Map)\n    # Watershed floods into 'low' areas. \n    # We want it to flood into 'high' probability areas, so we invert it.\n    cost_surface = 1 - prob_volume\n\n    # 3. Create a Mask\n    # We don't want the watershed to fill the entire 3D box.\n    # We only want it to grow where there is a reasonable chance of a sheet existing.\n    mask = prob_volume > prob_threshold\n\n    # 4. Run the Watershed\n    # - cost_surface:  sThe landscape the 'water' flows through.\n    # - labeled_seeds: Thetarting points (your broken segments).\n    # - mask: Restricts the growth to the sheet-like regions.\n    # - compactness: Increasing this (e.g. 0.01) makes the edges smoother if data is noisy.\n    print(\"Running 3D Watershed (this may take a moment for large volumes)...\")\n    filled_labels = watershed(cost_surface, labeled_seeds, mask=mask,watershed_line=True)\n    # ws_line = (filled_labels == 0) & mask\n    # barrier = binary_dilation(ws_line, iterations=1)\n    # filled_labels[barrier] = 0\n\n    # 5. Convert back to a binary sheet mask (optional)\n    # If you just want a single \"White=FG\" mask like your GT:\n    filled_sheets_binary = (filled_labels > 0).astype(np.uint8)\n\n    return filled_sheets_binary, filled_labels\n\n# --- Example Usage ---\n# prob_vol = load_your_probability_map()\n# hysteresis_seeds = load_your_thresholded_lines()\n# \n# result_binary, result_labels = fill_3d_sheets(prob_vol, hysteresis_seeds)\n\n\n\n\n\n# ==========================================\n# MAIN POST-PROCESSING LOGIC\n# ==========================================\ndef topo_postprocess(\n    probs,          # (D, H, W)\n    T_low=0.60,\n    T_high=0.80,\n    z_radius=0,\n    xy_radius=0,\n    dust_min_size=100,\n):\n    \n    # --- Step 1: 3D Hysteresis ---\n    strong = probs >= T_high\n    weak   = probs >= T_low\n\n    if not strong.any():\n        return np.zeros_like(probs, dtype=np.uint8)\n    struct_hyst = ndi.generate_binary_structure(3, 3)\n    mask = ndi.binary_propagation(strong, mask=weak, structure=struct_hyst)\n    if not mask.any():\n        return np.zeros_like(probs, dtype=np.uint8)\n    # --- Step 2: 3D Anisotropic Closing ---\n    if z_radius > 0 or xy_radius > 0:\n        struct_close = build_anisotropic_struct(z_radius, xy_radius)\n        if struct_close is not None:\n            mask = ndi.binary_closing(mask, structure=struct_close)\n\n    # --- Step 3: Dust Removal ---\n    if dust_min_size > 0:\n        mask = remove_small_objects(mask.astype(bool), min_size=dust_min_size)\n        \n    mask,_ = fill_3d_sheets(probs,mask)\n\n    \n    return mask.astype(np.uint8)\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"from skimage.filters import meijering\nimport numpy as np\n\ndef apply_meijering_slice_safe(img_slice, sigma_range, filter_padding=20):\n    padded = np.pad(img_slice, filter_padding, mode='reflect')\n    filtered = meijering(padded, sigmas=sigma_range, black_ridges=False)\n    if filter_padding > 0:\n        return filtered[filter_padding:-filter_padding, filter_padding:-filter_padding]\n    return filtered\n\n\ndef hybrid_meijering_threshold_slice(\n    prob_slice,\n    margin,\n    sigma_min=1,\n    sigma_max=4,\n    meijering_thresh=0.1,\n    simple_thresh=0.8,\n    pad=20,\n):\n    H, W = prob_slice.shape\n\n    center_mask = np.zeros_like(prob_slice, dtype=bool)\n    if margin == 0:\n        center_mask[:] = True\n    elif margin < H//2 and margin < W//2:\n        center_mask[margin:H-margin, margin:W-margin] = True\n    else:\n        center_mask[:] = True\n\n    tubeness = apply_meijering_slice_safe(\n        prob_slice,\n        range(sigma_min, max(sigma_min+1,sigma_max)),\n        filter_padding=pad\n    )\n\n    if tubeness.max() > 0:\n        tubeness /= tubeness.max()\n\n    mask_inner = tubeness > meijering_thresh\n    mask_outer = prob_slice > simple_thresh\n\n    return np.where(center_mask, mask_inner, mask_outer)\n\nimport scipy.ndimage as ndi\nfrom skimage.morphology import remove_small_objects\nfrom joblib import Parallel, delayed\n\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def load_probabilities(npz_path: Path) -> np.ndarray:\n    \"\"\"\n    Load probability maps from nnUNet inference.\n    \n    Only available if inference was run with save_probabilities=True.\n    Shape: (num_classes, D, H, W) with float32 values in [0, 1].\n    \"\"\"\n    data = np.load(npz_path)\n    return data['probabilities']\n\n\n\n\ndef predictions_to_tiff(\n    pred_dir: Path,\n    output_dir: Path,\n    *,\n    T_low=0.25,\n    T_high=0.250,##CHANGING THIS DOESN MEAN SHIT\n    z_radius=0,\n    xy_radius=0,\n    dust_min_size=100,\n):\n    \"\"\"\n    Convert nnUNet probability outputs to TIFF using topology-aware postprocessing.\n    \"\"\"\n\n    output_dir.mkdir(parents=True, exist_ok=True)\n\n    npz_files = sorted(pred_dir.glob(\"*.npz\"))\n    tif_files = sorted(pred_dir.glob(\"*.tif\"))\n    nii_files = sorted(pred_dir.glob(\"*.nii.gz\"))\n\n    if npz_files:\n        print(f\"Post-processing {len(npz_files)} NPZ probability files...\")\n\n        for npz_path in tqdm(npz_files, desc=\"Topo postprocess\"):\n            case_id = npz_path.stem\n\n            # ---- load nnUNet probabilities ----\n            probs = load_probabilities(npz_path)   # (C, D, H, W)\n            prob_surface = probs[1]                # class-1 = surface\n\n            # ---- YOUR helper is used here ----\n            # pred = apply_frangi_pipeline(prob_surface)\n            pred = topo_postprocess(\n                prob_surface\n            )\n\n            tifffile.imwrite(\n                output_dir / f\"{case_id}.tif\",\n                pred.astype(np.uint8),\n                compression=\"zlib\",\n            )\n\n    elif tif_files:\n        print(\"No probabilities found — copying TIFF predictions\")\n        for tif_path in tqdm(tif_files):\n            pred = tifffile.imread(tif_path).astype(np.uint8)\n            tifffile.imwrite(output_dir / tif_path.name, pred)\n\n    elif nii_files:\n        print(\"Converting NIfTI predictions\")\n        for nii_path in tqdm(nii_files):\n            case_id = nii_path.stem.replace(\".nii\", \"\")\n            pred = load_nifti(nii_path).astype(np.uint8)\n            tifffile.imwrite(output_dir / f\"{case_id}.tif\", pred)\n\n    else:\n        raise RuntimeError(f\"No predictions found in {pred_dir}\")\n","metadata":{"_uuid":"d5edfb02-bf4f-4782-8c10-9b6a535f9def","_cell_guid":"77f89718-b429-4287-bcdb-c2977dc93fb4","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 10. Full Pipeline\n\n**Section Summary:**\nMain entry point that orchestrates the entire nnUNet pipeline.\n\n**Pipeline Steps:**\n1. Environment setup (create directories, set env vars)\n2. Raw data preparation (symlinks + JSON sidecars)\n3. Preprocessing (or link pre-prepared data)\n4. Training\n5. Inference + TIFF conversion\n6. Visualization (progress + predictions)\n\n**Usage Patterns:**\n```python\n# Full pipeline with defaults\nfull_pipeline()\n\n# Quick experiment\nfull_pipeline(epochs=50, config=\"2d\")\n\n# Resume training\nfull_pipeline(continue_training=True, epochs=250)\n\n# Inference only (trained model exists)\nfull_pipeline(do_preprocess=False, do_train=False)\n```\n\n**Important Notes:**\n- When `fold=\"all\"`, validation runs on ALL training data after training\n- This validation step can take 50+ minutes but cannot be skipped\n- Progress image and predictions are auto-displayed after completion","metadata":{"_uuid":"5462cd2d-43da-4208-aaea-d978d66ff6e8","_cell_guid":"cbcac86d-dfc8-4dc2-b7de-a00ddf092c58","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def _parse_model_path(model_path: Union[str, Path, None]) -> Tuple[Optional[Path], Optional[int], Optional[str], Optional[str], Optional[str]]:\n    \"\"\"\n    Parse model path to extract configuration parameters.\n    \n    Model path format: .../DatasetXXX_Name/TrainerName__PlansName__Config/fold_X/checkpoint.pth\n    \n    Returns:\n        (model_dir, epochs, plans, config, fold) - extracted from path, None if not found\n    \n    Note: Parses path structure even if file doesn't exist (for path validation).\n    \"\"\"\n    if model_path is None:\n        return None, None, None, None, None\n    \n    model_path = Path(model_path)\n    \n    # Determine model directory from path (even if doesn't exist)\n    # If path ends with .pth, use parent directory\n    if model_path.suffix == \".pth\" or (model_path.exists() and model_path.is_file()):\n        model_dir = model_path.parent\n    else:\n        model_dir = model_path\n    \n    # Try to parse folder structure\n    # Expected: .../fold_X or .../TrainerName__PlansName__Config/fold_X\n    try:\n        parts = model_dir.parts\n        \n        # Find fold\n        fold = None\n        for part in reversed(parts):\n            if part.startswith(\"fold_\"):\n                fold = part.replace(\"fold_\", \"\")\n                break\n        \n        # Find trainer__plans__config\n        epochs = None\n        plans = None\n        config = None\n        for part in parts:\n            if \"__\" not in part:\n                continue\n            segments = part.split(\"__\")\n            if len(segments) >= 3:\n                trainer_name = segments[0]\n                plans = segments[1]\n                config = segments[2]\n                # Extract epochs from trainer name\n                if \"epochs\" in trainer_name:\n                    import re\n                    match = re.search(r'(\\d+)epochs?', trainer_name)\n                    if match:\n                        epochs = int(match.group(1))\n                elif trainer_name == \"nnUNetTrainer\":\n                    epochs = 1000  # Default\n                break\n        \n        return model_dir, epochs, plans, config, fold\n    except Exception:\n        return model_dir, None, None, None, None","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def full_pipeline(\n    # Data options\n    max_cases: Optional[int] = None,\n    # Stage control\n    do_preprocess: bool = True,\n    do_train: bool = True,\n    do_inference: bool = True,\n    # Training options (all tunable parameters)\n    config: str = CONFIGURATION,\n    fold: Union[int, str] = FOLD,\n    planner: str = PLANNER,\n    plans: str = PLANS_NAME,\n    epochs: Optional[Epochs] = EPOCHS,\n    pretrained_weights: Optional[Path] = None,\n    continue_training: bool = False,\n    num_gpus: int = NUM_GPUS,\n    # Inference options\n    save_probabilities: bool = True,\n    # External model (for inference without training)\n    model_path: Optional[Union[str, Path]] = None,\n    # Timeout\n    timeout: Optional[int] = COMMAND_TIMEOUT,\n):\n    \"\"\"\n    Run complete pipeline: setup -> data prep -> preprocess -> train -> predict.\n    \n    Args:\n        max_cases: Limit number of training cases (None = use all)\n        do_preprocess: Run preprocessing step\n        do_train: Run training step\n        do_inference: Run inference step\n        \n        config: nnUNet configuration (3d_fullres, 2d, 3d_lowres, 3d_cascade_fullres)\n        fold: Fold number (0-4) or \"all\" for training on all data\n        planner: Planner class name\n        plans: Plans name matching the planner\n        epochs: Number of training epochs (1, 5, 10, 20, 50, 100, 150, 200, 250, \n                300, 400, 500, 750, 1000, 2000, 4000, 8000). None = 1000\n        pretrained_weights: Path to pretrained weights for fine-tuning\n        continue_training: Continue from last checkpoint\n        num_gpus: Number of GPUs for DDP training\n        \n        save_probabilities: Save probability maps during inference\n        \n        model_path: Path to trained model checkpoint or directory (str or Path).\n                    When provided with do_train=False, parameters (epochs, plans, config, fold)\n                    are auto-extracted from the path if possible.\n                    Example: \"/path/to/nnUNetTrainer_5epochs__nnUNetResEncUNetMPlans__3d_fullres/fold_all/checkpoint_best.pth\"\n        timeout: Command timeout in seconds (None for no timeout)\n    \n    Returns:\n        True if pipeline completed successfully\n    \"\"\"\n    \n    print(\"=\" * 60)\n    print(\"Vesuvius Surface Detection - nnUNet Pipeline\")\n    print(\"=\" * 60)\n    print(f\"Stages: preprocess={do_preprocess}, train={do_train}, inference={do_inference}\")\n    print(f\"Config: {config}, Fold: {fold}, Epochs: {epochs or 1000}, GPUs: {num_gpus}\")\n    \n    # 1. Setup\n    print(\"\\n[1/5] Environment setup...\")\n    setup_environment()\n    \n    # 2. Prepare raw data (always - uses symlinks, fast)\n    print(\"\\n[2/5] Preparing raw dataset (symlinks)...\")\n    raw_dataset_dir = NNUNET_RAW / DATASET_NAME\n    if not raw_dataset_dir.exists():\n        prepare_dataset(INPUT_DIR, max_cases=max_cases)\n    else:\n        print(f\"Raw dataset already exists: {raw_dataset_dir}\")\n    \n    # 3. Preprocessing\n    if do_preprocess:\n        # Check if pre-prepared preprocessed data exists\n        if _link_prepared_preprocessed():\n            print(\"\\n[3/5] Using pre-prepared preprocessed data...\")\n        else:\n            print(\"\\n[3/5] Preprocessing...\")\n            success = run_preprocessing(planner=planner, configurations=[config], timeout=timeout)\n            if not success:\n                print(\"Preprocessing failed!\")\n                return False\n    else:\n        print(\"\\n[3/5] Skipping preprocessing...\")\n        _link_prepared_preprocessed()  # Still link if available\n    \n    # 4. Training\n    if do_train:\n        print(\"\\n[4/5] Training...\")\n        success = run_training(\n            config=config,\n            fold=fold,\n            plans=plans,\n            epochs=epochs,\n            pretrained_weights=pretrained_weights,\n            continue_training=continue_training,\n            num_gpus=num_gpus,\n            timeout=timeout\n        )\n        if not success:\n            print(\"Training failed!\")\n            return False\n    else:\n        print(\"\\n[4/5] Skipping training...\")\n        \n        # Parse model_path to extract configuration if provided\n        if model_path is not None:\n            model_path = Path(model_path) if isinstance(model_path, str) else model_path\n            _, parsed_epochs, parsed_plans, parsed_config, parsed_fold = _parse_model_path(model_path)\n            \n            # Override parameters with parsed values\n            if parsed_epochs is not None:\n                epochs = parsed_epochs\n                print(f\"  Detected epochs={epochs} from model path\")\n            if parsed_plans is not None:\n                plans = parsed_plans\n                print(f\"  Detected plans={plans} from model path\")\n            if parsed_config is not None:\n                config = parsed_config\n                print(f\"  Detected config={config} from model path\")\n            if parsed_fold is not None:\n                fold = parsed_fold\n                print(f\"  Detected fold={fold} from model path\")\n            \n            if model_path.exists():\n                print(f\"Using model: {model_path}\")\n            else:\n                print(f\"WARNING: model_path does not exist: {model_path}\")\n        \n        # Verify model exists\n        expected_model_dir = get_training_output_dir(epochs=epochs, plans=plans, config=config, fold=fold)\n        expected_model_dir.mkdir(exist_ok=True, parents=True)\n        checkpoint_final = expected_model_dir / \"checkpoint_final.pth\"\n        checkpoint_best = expected_model_dir / \"checkpoint_best.pth\"\n\n        if model_path and Path(model_path).exists():\n            # Symlink external -> final so nnUNet finds it\n            checkpoint_final.symlink_to(model_path)\n            print(f\"Map model: {checkpoint_final}\")\n        elif checkpoint_final.exists():\n            print(f\"Found model: {checkpoint_final}\")\n        elif checkpoint_best.exists():\n            # Symlink best -> final so nnUNet finds it\n            checkpoint_final.symlink_to(checkpoint_best)\n            print(f\"Found model: {checkpoint_best} (symlinked to checkpoint_final.pth)\")\n        elif do_inference:\n            print(f\"WARNING: No model found at {expected_model_dir}\")\n            print(\"  Provide model_path to a valid nnUNet checkpoint file\")\n    \n    # 5. Inference\n    if do_inference:\n        print(\"\\n[5/5] Running inference...\")\n        \n        # Prepare test data in temp location\n        test_input_dir = WORKING_DIR / \"test_input\"\n        prepare_test_data(INPUT_DIR, test_input_dir)\n        \n        # Run inference\n        predictions_dir = WORKING_DIR / \"predictions\"\n        success = run_inference(\n            test_input_dir, \n            predictions_dir,\n            config=config,\n            fold=fold,\n            plans=plans,\n            epochs=epochs,\n            save_probabilities=save_probabilities,\n            timeout=timeout\n        )\n        if not success:\n            print(\"Inference failed!\")\n            return False\n        \n        # Convert predictions to TIFF\n        print(\"\\nConverting predictions to TIFF...\")\n        tiff_output_dir = OUTPUT_DIR / \"predictions_tiff\"\n        predictions_to_tiff(predictions_dir, tiff_output_dir)\n        \n        print(f\"\\nPredictions saved to: {tiff_output_dir}\")\n    else:\n        print(\"\\n[5/5] Skipping inference...\")\n    \n    print(\"\\n\" + \"=\" * 60)\n    print(\"Pipeline complete!\")\n    print(\"=\" * 60)\n    \n    # Show training progress if training was done\n    if do_train:\n        print(\"\\n[Visualization] Training progress:\")\n        show_progress(epochs=epochs, plans=plans, config=config, fold=fold)\n    \n    # Visualize predictions if inference was done\n    if do_inference:\n        print(\"\\n[Visualization] Sample prediction:\")\n        visualize_predictions(num_samples=1)\n    \n    return True","metadata":{"_uuid":"9079d68c-f6f7-4363-9090-849178090320","_cell_guid":"42cebe6a-6001-4cad-9755-3d79e31e51ad","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 11. Create Submission\n\n**Section Summary:**\nCreates a ZIP file for Kaggle submission from prediction TIFF files.\n\n**Expected Format:**\n- ZIP containing .tif prediction files\n- Filenames should match test image names\n- Predictions should be uint8 (0=background, 1=surface)\n\n**Output Location:**\n- `/kaggle/working/submission.zip`","metadata":{"_uuid":"5c82059b-90f5-4cee-8478-c4f74d1ff880","_cell_guid":"19eba4e0-b5fe-4471-a7ac-575a9e869df5","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def generate_submission(\n    predictions_tiff_dir: Path = OUTPUT_DIR / \"predictions_tiff\",\n    output_zip: Path = OUTPUT_DIR / \"submission.zip\",\n    delete_after_zip: bool = True  # Default True to save space on Kaggle\n) -> Optional[Path]:\n    \"\"\"\n    Create submission ZIP from TIFF predictions.\n    \n    Args:\n        predictions_tiff_dir: Directory containing predicted TIFF files\n        output_zip: Output ZIP file path\n        delete_after_zip: Delete TIFF files after adding to ZIP (saves space)\n    \n    Returns:\n        Path to submission ZIP if successful, None otherwise\n    \"\"\"\n    import zipfile\n    \n    if not predictions_tiff_dir.exists():\n        print(f\"ERROR: Predictions directory not found: {predictions_tiff_dir}\")\n        print(\"Run inference first!\")\n        return None\n    \n    tiff_files = sorted(predictions_tiff_dir.glob(\"*.tif\"))\n    \n    if not tiff_files:\n        print(f\"No TIFF files found in {predictions_tiff_dir}\")\n        return None\n    \n    print(f\"Creating submission ZIP with {len(tiff_files)} files...\")\n    \n    with zipfile.ZipFile(output_zip, 'w', zipfile.ZIP_DEFLATED) as zipf:\n        for tiff_path in tqdm(tiff_files, desc=\"Zipping predictions\"):\n            # Add file with just the filename (no directory structure)\n            zipf.write(tiff_path, tiff_path.name)\n            \n            if delete_after_zip:\n                tiff_path.unlink()\n    \n    zip_size_mb = output_zip.stat().st_size / (1024 * 1024)\n    print(f\"Submission saved: {output_zip} ({zip_size_mb:.1f} MB)\")\n    \n    return output_zip\n\n\n# Example usage:\n# generate_submission()\n# generate_submission(delete_after_zip=True)  # Delete TIFFs after zipping to save space","metadata":{"_uuid":"7defcaab-59f8-4f13-baa4-b1de03be7bd4","_cell_guid":"540caf8b-e834-4fc6-9ec9-7ca77ca1ffc0","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 12. Visualization\n\n**Section Summary:**\nFunctions for visualizing training progress and predictions.\n\n**`show_progress()`:**\n- Displays nnUNet's training progress image (loss curves, metrics)\n- Auto-detects correct path based on configuration\n\n**`plot_three_axis_cuts()`:**\n- Shows middle slices of XY, XZ, YZ planes\n- Displays image and mask side-by-side\n\n**`visualize_predictions()`:**\n- Auto-finds predictions and corresponding input images\n- Displays multiple samples if requested","metadata":{"_uuid":"a8d74280-4aad-486e-8f1f-a125033c334d","_cell_guid":"1267b34d-c872-413f-a830-157ec2492bf8","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def show_progress(\n    epochs: Optional[Epochs] = None,\n    plans: str = PLANS_NAME,\n    config: str = CONFIGURATION,\n    fold: Union[int, str] = FOLD\n):\n    \"\"\"\n    Display training progress image from nnUNet.\n    \n    Args:\n        epochs: Number of epochs used for training (to find correct folder)\n        plans: Plans name\n        config: Configuration name\n        fold: Fold number or \"all\"\n    \"\"\"\n    progress_path = get_progress_image_path(epochs, plans, config, fold)\n    \n    if not progress_path.exists():\n        print(f\"Progress image not found: {progress_path}\")\n        print(\"Training may not have started or completed yet.\")\n        return\n    \n    from IPython.display import Image, display\n    print(f\"Training progress: {progress_path}\")\n    display(Image(filename=str(progress_path)))\n\n\ndef plot_three_axis_cuts(\n    image_vol_path: Path,\n    mask_vol_path: Path,\n    figsize: tuple = (12, 15)\n):\n    \"\"\"\n    Plot middle slices of XY, XZ, and YZ planes for image volume and predicted mask.\n    \n    Args:\n        image_vol_path: Path to input image TIFF\n        mask_vol_path: Path to prediction mask TIFF\n        figsize: Figure size (width, height)\n    \"\"\"\n    print(f\"Visualizing: {image_vol_path.name}\")\n    \n    # Load volumes\n    image_vol = tifffile.imread(str(image_vol_path))\n    mask_vol = tifffile.imread(str(mask_vol_path)).astype(np.uint8)\n    \n    # Get dimensions\n    d, h, w = image_vol.shape\n    z_mid, y_mid, x_mid = d // 2, h // 2, w // 2\n    \n    # Extract slices\n    slices = {\n        'XY Plane (Z-axis)': (image_vol[z_mid, :, :], mask_vol[z_mid, :, :]),\n        'XZ Plane (Y-axis)': (image_vol[:, y_mid, :], mask_vol[:, y_mid, :]),\n        'YZ Plane (X-axis)': (image_vol[:, :, x_mid], mask_vol[:, :, x_mid])\n    }\n    \n    fig, axes = plt.subplots(3, 2, figsize=figsize)\n    for i, (plane_name, (img_slice, mask_slice)) in enumerate(slices.items()):\n        # Image Volume\n        axes[i, 0].imshow(img_slice, cmap='gray')\n        axes[i, 0].set_title(f\"{plane_name} - Image Volume\")\n        axes[i, 0].axis('off')\n        \n        # Mask\n        axes[i, 1].imshow(mask_slice, cmap='gray')\n        axes[i, 1].set_title(f\"{plane_name} - Predicted Mask\")\n        axes[i, 1].axis('off')\n    \n    plt.tight_layout()\n    plt.show()\n\n\ndef visualize_predictions(\n    predictions_dir: Path = OUTPUT_DIR / \"predictions_tiff\",\n    test_images_dir: Path = INPUT_DIR / \"test_images\",\n    num_samples: int = 1\n):\n    \"\"\"\n    Visualize prediction results by showing image/mask pairs.\n    \n    Args:\n        predictions_dir: Directory containing prediction TIFFs\n        test_images_dir: Directory containing input test images\n        num_samples: Number of samples to visualize\n    \"\"\"\n    if not predictions_dir.exists():\n        print(f\"Predictions directory not found: {predictions_dir}\")\n        return\n    \n    predictions = sorted(predictions_dir.glob(\"*.tif\"))\n    if not predictions:\n        print(f\"No TIFF predictions found in {predictions_dir}\")\n        return\n    \n    print(f\"Found {len(predictions)} predictions\")\n    \n    for pred_path in predictions[:num_samples]:\n        image_path = test_images_dir / pred_path.name\n        if image_path.exists():\n            plot_three_axis_cuts(image_path, pred_path)\n        else:\n            print(f\"Warning: Input image not found: {image_path}\")","metadata":{"_uuid":"16d62e5e-eb9a-4fea-bfa0-60912fb0380a","_cell_guid":"9e7cd5bf-abbf-4c5e-995e-4bd79ee13c0c","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 13. Run Pipeline\n\n**Section Summary:**\nExecute the main pipeline and generate submission.\n\n**Default Configuration:**\n- `fold=\"all\"` - Train on all data (no cross-validation)\n- `config=\"3d_fullres\"` - Best quality single-stage\n- `epochs=1000` - Full training (can reduce to 250-500)\n- `num_gpus=auto` - Uses all available GPUs\n\n**After Running:**\n1. Training progress displayed (progress.png)\n2. Sample prediction visualized\n3. Submission ZIP created at `/kaggle/working/submission.zip`","metadata":{"_uuid":"f58def5b-2797-4470-8310-e82d9eb07f97","_cell_guid":"a87af08c-9922-46c3-9a62-94e2fb7064d8","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"!rm -rf nnUNet_results/\n!cp -r /kaggle/input/nnunet-192models/ nnUNet_results/\n# dont leak the used models for inference\n# !rm -rf /kaggle/working/nnUNet_results/Dataset100_VesuviusSurface/nnUNetTrainer_750epochs__nnUNetResEncUNetMPlans__3d_fullres/fold_1\n# !rm -f /kaggle/working/nnUNet_results/Dataset100_VesuviusSurface/nnUNetTrainer_750epochs__nnUNetResEncUNetMPlans__3d_fullres/fold_1/checkpoint_best.pth\n!rm -f /kaggle/working/nnUNet_results/Dataset100_VesuviusSurface/nnUNetTrainer_500epochs__nnUNetResEncUNetLPlans__3d_fullres/fold_all/checkpoint_final.pth","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import json\n\nsplits_path = \"/kaggle/input/split-final/splits_final.json\"\nfold = 0  # change this to the fold you want\n\nwith open(splits_path, \"r\") as f:\n    splits = json.load(f)\n\nval_cases = splits[fold][\"val\"]\n\nprint(f\"Fold {fold} validation cases ({len(val_cases)}):\")\nprint(val_cases)\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# # =============================================================================\n# # FAST ONE-SHOT INFERENCE ON SELECTED TRAINING CASES\n# # =============================================================================\n\n# import shutil\n# import numpy as np\n# from pathlib import Path\n\n# # -----------------------------------------------------------------------------\n# # CONFIG\n# # -----------------------------------------------------------------------------\n\n# MODEL_CHECKPOINT = Path(\n# \"/kaggle/input/nnunet-192models/Dataset100_VesuviusSurface/nnUNetTrainer__nnUNetResEncUNetMPlans__3d_fullres/fold_all/checkpoint_best.pth\"\n# )\n\n# TRAIN_CASES = ['1006462223']\n\n# DTYPE = np.float16\n# SURFACE_CLASS = 1\n\n# # -----------------------------------------------------------------------------\n# # PATHS\n# # -----------------------------------------------------------------------------\n\n# FAST_INPUT_DIR = WORKING_DIR / \"fast_selected_imagesTs\"\n# PRED_DIR = WORKING_DIR / \"fast_predictions\"\n# OUT_NPY = OUTPUT_DIR / \"probabilities_npy_fast\"\n\n# train_images_dir = INPUT_DIR / \"train_images\"\n\n# # -----------------------------------------------------------------------------\n# # 1. MAP MODEL → nnUNet expected structure\n# # -----------------------------------------------------------------------------\n\n# assert MODEL_CHECKPOINT.exists()\n\n# model_dir, epochs, plans, config, fold = _parse_model_path(MODEL_CHECKPOINT)\n# assert all([epochs, plans, config, fold])\n\n# expected_model_dir = get_training_output_dir(\n#     epochs=epochs, plans=plans, config=config, fold=fold\n# )\n# expected_model_dir.mkdir(parents=True, exist_ok=True)\n\n# checkpoint_final = expected_model_dir / \"checkpoint_final.pth\"\n# if not checkpoint_final.exists():\n#     checkpoint_final.symlink_to(MODEL_CHECKPOINT)\n\n# print(\"Using model:\", checkpoint_final)\n\n# # -----------------------------------------------------------------------------\n# # 2. Build proper batched nnUNet test folder (FAST PATH)\n# # -----------------------------------------------------------------------------\n\n# if FAST_INPUT_DIR.exists():\n#     shutil.rmtree(FAST_INPUT_DIR)\n# FAST_INPUT_DIR.mkdir(parents=True)\n\n# print(f\"Preparing {len(TRAIN_CASES)} cases (fast layout)...\")\n\n# for cid in TRAIN_CASES:\n#     src = train_images_dir / f\"{cid}.tif\"\n#     dst = FAST_INPUT_DIR / f\"{cid}_0000.tif\"\n#     dst.symlink_to(src)\n\n# print(\"Done.\")\n\n# # -----------------------------------------------------------------------------\n# # 3. Run inference ONCE (batched, efficient)\n# # -----------------------------------------------------------------------------\n\n# print(\"\\nRunning fast nnUNet inference...\")\n\n# success = run_inference(\n#     input_dir=FAST_INPUT_DIR,\n#     output_dir=PRED_DIR,\n#     config=config,\n#     fold=fold,\n#     plans=plans,\n#     epochs=epochs,\n#     save_probabilities=True,  # still needed to extract class probs\n# )\n\n# if not success:\n#     raise RuntimeError(\"Inference failed\")\n\n# # -----------------------------------------------------------------------------\n# # 4. Extract and save surface probabilities efficiently\n# # -----------------------------------------------------------------------------\n\n# OUT_NPY.mkdir(parents=True, exist_ok=True)\n\n# print(\"\\nExtracting probabilities...\")\n\n# saved = 0\n# for npz_path in PRED_DIR.glob(\"*.npz\"):\n#     case_id = npz_path.stem\n#     if case_id not in TRAIN_CASES:\n#         continue\n\n#     with np.load(npz_path, mmap_mode=\"r\") as data:\n#         probs = data[\"probabilities\"]  # (C, Z, Y, X)\n#         surface = probs[SURFACE_CLASS].astype(DTYPE)\n\n#     np.save(OUT_NPY / f\"{case_id}.npy\", surface)\n\n#     print(f\"[OK] {case_id} | shape={surface.shape}\")\n#     saved += 1\n\n# assert saved > 0\n# print(f\"\\nDONE. Saved {saved} files → {OUT_NPY}\")\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# !rm -rf nnUNet_results/","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# # Sample pipeline with 50 epochs and time limit\n# # full_pipeline(timeout=10*60*60)\n\n# # Then resume when training finishes and run inference with the best model\nfull_pipeline(\n    do_preprocess=False, do_train=False,\n    model_path=\"/kaggle/input/nnunet-192models/Dataset100_VesuviusSurface/nnUNetTrainer_500epochs__nnUNetResEncUNetLPlans__3d_fullres/fold_all/checkpoint_best.pth\",\n)\n\n# Generate submission ZIP\ngenerate_submission()","metadata":{"_uuid":"463f2b49-d26e-45c9-9472-d0779a99322e","_cell_guid":"07aca3f2-b73c-47c9-916c-0683c19971d1","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}