{"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":"none","dataSources":[{"sourceId":117682,"databundleVersionId":14443416,"sourceType":"competition"},{"sourceId":13784990,"sourceType":"datasetVersion","datasetId":8775318}],"dockerImageVersionId":31192,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"!pip install /kaggle/input/imagedecodecs-whl/imagecodecs-2025.11.11-cp311-abi3-manylinux_2_27_x86_64.manylinux_2_28_x86_64.whl","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-11-19T09:00:10.267494Z","iopub.execute_input":"2025-11-19T09:00:10.267808Z","iopub.status.idle":"2025-11-19T09:00:16.012745Z","shell.execute_reply.started":"2025-11-19T09:00:10.267784Z","shell.execute_reply":"2025-11-19T09:00:16.011887Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!pip install imageio","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-19T09:00:18.761170Z","iopub.execute_input":"2025-11-19T09:00:18.761734Z","iopub.status.idle":"2025-11-19T09:00:22.421595Z","shell.execute_reply.started":"2025-11-19T09:00:18.761704Z","shell.execute_reply":"2025-11-19T09:00:22.420534Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import imageio.v3 as iio\nimport numpy as np\nfrom scipy.ndimage import gaussian_filter\nfrom skimage.morphology import disk, binary_dilation\nimport matplotlib.pyplot as plt\nimport os\nimport glob\nfrom scipy.ndimage import distance_transform_edt # Nécessaire si vous recalculez E_LABEL\n\n# --- CONSTANTES GLOBALES (Basées sur votre environnement) ---\nPATH_VOLUME_ROOT = \"/kaggle/input/vesuvius-challenge-surface-detection/train_images/\"\nPATH_LABEL_ROOT = \"/kaggle/input/vesuvius-challenge-surface-detection/train_labels/\" \nTOTAL_SLICES = 320\nREF_Z = TOTAL_SLICES // 2 # 160\n\n# Projection 3D élargie\nZ_SHIFTS = list(range(-5, 6))      # [-5, -4, ..., 0, ..., +5]\n\n# --- PARAMÈTRES DU DERNIER TEST GLOBAL (IoU 0.0611) ---\nQUANTILE_THRESHOLD = 0.9995         # Seuil qui a donné des 0 extrêmes\nGAUSSIAN_SIGMA = 1.0                # Lissage\n\n# --- ÉPAISSEUR MOYENNE (Basée sur votre calcul) ---\nE_LABEL = 12.39\nBASE_RADIUS = int(E_LABEL / 2) # 6 px\n\n# --- IDs où le score IoU était 0.0000 ---\n# Ces IDs sont les cas limites que nous devons analyser.\nproblem_ids = [\n    \"3060150865\", \"4024699648\", \"3020371188\", \"693501383\", \"3137156884\", \n    \"3686818985\", \"2945808105\", \"601001728\", \"2116132949\", \"1967300661\",\n    \"3664792107\", \"2051981635\" \n] \n# =========================================================\n# FONCTIONS DE BASE (Re-définies pour garantir l'exécution)\n# =========================================================\n\ndef detect_bright_dots(ct_slice, quantile_thr):    \n    \"\"\"Applique le seuil quantile ultra-sélectif.\"\"\"\n    thr = np.quantile(ct_slice, quantile_thr)\n    mask = (ct_slice >= thr).astype(np.uint8)\n    return mask, thr\n\ndef dynamic_dilation(mask, base_radius):\n    \"\"\"Dilatation dynamique basée sur le nombre de pixels détectés.\"\"\"\n    n = mask.sum()\n\n    if n > 10000:\n        r = max(1, base_radius - 4)\n    elif n > 2000:\n        r = max(1, base_radius - 2)\n    elif n > 200:\n        r = base_radius\n    else:\n        r = base_radius + 2 \n\n    se = disk(r)\n    dil = binary_dilation(mask, se).astype(np.uint8)\n\n    return dil, r\n\n# =========================================================\n# FONCTION DE VISUALISATION (Adaptation de votre run_projection_3d)\n# =========================================================\n\ndef visualize_zero_iou(image_id):\n    \"\"\"Calcule, affiche les masques et les détails de seuillage pour un ID donné.\"\"\"\n    \n    print(f\"\\n===== ID {image_id} | SIGMA = {GAUSSIAN_SIGMA:.1f} | QUANTILE = {QUANTILE_THRESHOLD:.5f} =====\")\n\n    # Load\n    try:\n        ct_vol = iio.imread(f\"{PATH_VOLUME_ROOT}{image_id}.tif\")\n        label_vol = iio.imread(f\"{PATH_LABEL_ROOT}{image_id}.tif\")\n    except Exception as e:\n        print(f\"[ERREUR] Impossible de charger {image_id}: {e}\")\n        return\n\n    true_mask = (label_vol[REF_Z] == 1).astype(np.uint8)\n    final_mask = np.zeros_like(true_mask)\n\n    for dz in Z_SHIFTS:\n        z = REF_Z + dz\n        if z < 0 or z >= TOTAL_SLICES:\n            continue\n\n        ct_slice = ct_vol[z]\n        \n        # 1. Filtre Gaussien\n        ct_slice_filtered = gaussian_filter(ct_slice, sigma=GAUSSIAN_SIGMA)\n        \n        # 2. Seuil adaptatif (utilisation du seuil global défini)\n        dots, thr = detect_bright_dots(ct_slice_filtered, QUANTILE_THRESHOLD)\n\n        # 3. Dilatation dynamique\n        dilated, used_r = dynamic_dilation(dots, BASE_RADIUS)\n\n        final_mask = (final_mask | dilated).astype(np.uint8)\n\n        print(f\"  • Z {dz:+d} | seuil={thr:.1f} | radius={used_r} | pixels={dots.sum()}\")\n\n    # Score IoU\n    union = np.sum((final_mask | true_mask) > 0)\n    inter = np.sum((final_mask & true_mask) > 0)\n    iou = inter / union if union > 0 else 0\n\n    print(f\"  → IoU final pour {image_id}: {iou:.4f}\")\n    \n    # --- VISUALISATION ---\n    fig, ax = plt.subplots(1, 3, figsize=(28, 9))\n    ref_ct = ct_vol[REF_Z]\n\n    ax[0].imshow(ref_ct, cmap='gray')\n    ax[0].set_title(\"Image CT (Z=ref)\")\n    ax[0].axis(\"off\")\n\n    ax[1].imshow(ref_ct, cmap='gray')\n    ax[1].imshow(final_mask, cmap='Greens', alpha=0.6) # Prédiction (Vert)\n    ax[1].set_title(f\"Détection Bright Dots (IoU={iou:.4f})\")\n    ax[1].axis(\"off\")\n\n    ax[2].imshow(ref_ct, cmap='gray')\n    ax[2].imshow(true_mask, cmap='Reds', alpha=0.6) # Vérité Terrain (Rouge)\n    ax[2].set_title(\"Label (Vérité Terrain - Rouge)\")\n    ax[2].axis(\"off\")\n\n    plt.show()\n\n# =========================================================\n# EXÉCUTION DE LA VISUALISATION CIBLÉE\n# =========================================================\nprint(f\"--- DÉBUT DE LA VISUALISATION DES {len(problem_ids)} CAS D'ÉCHEC (IoU=0) ---\")\n\nfor img_id in problem_ids:\n    visualize_zero_iou(img_id)\n\nprint(\"--- ANALYSE TERMINÉE : PRÊT POUR LE GRID SEARCH ---\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-19T09:00:28.411549Z","iopub.execute_input":"2025-11-19T09:00:28.411902Z","iopub.status.idle":"2025-11-19T09:00:28.805546Z","shell.execute_reply.started":"2025-11-19T09:00:28.411865Z","shell.execute_reply":"2025-11-19T09:00:28.804813Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import imageio.v3 as iio\nimport numpy as np\nimport random\nfrom scipy.ndimage import distance_transform_edt, gaussian_filter\nfrom skimage.morphology import disk, binary_dilation\nimport os\nimport glob\nimport matplotlib.pyplot as plt # Importé mais utilisé uniquement pour run_projection_3d s'il était réactivé\n\n# --- CONSTANTES ---\nPATH_VOLUME_ROOT = \"/kaggle/input/vesuvius-challenge-surface-detection/train_images/\"\nPATH_LABEL_ROOT = \"/kaggle/input/vesuvius-challenge-surface-detection/train_labels/\"\nTOTAL_SLICES = 320\nREF_Z = TOTAL_SLICES // 2 # 160\n\n# Projection 3D élargie\nZ_SHIFTS = list(range(-5, 6))      # [-5, -4, ..., 0, ..., +5]\n\n# --- PARAMÈTRES GLOBALS MODIFIABLES PAR LA GRID SEARCH ---\nQUANTILE_THRESHOLD = 0.9995         \nGAUSSIAN_SIGMA = 1.0                \nE_LABEL = 1.0 # Sera recalculé\n\n# =========================================================\n# A. Charger TOUS les IDs d'entraînement\n# =========================================================\ndef load_all_train_ids(path_root):\n    files = glob.glob(os.path.join(path_root, \"*.tif\"))\n    ids = [os.path.basename(f).replace(\".tif\", \"\") for f in files]\n    print(f\"[INFO] {len(ids)} IDs d'entraînement chargés depuis le répertoire.\")\n    return ids\ntrain_ids_all = load_all_train_ids(PATH_VOLUME_ROOT)\n\n# =========================================================\n# 1. Calcul propre de l’épaisseur moyenne E_label\n# =========================================================\ndef compute_average_ink_thickness(num_samples=10):\n    thicknesses = []\n    samples = [\n        (random.choice(train_ids_all), random.randint(20, TOTAL_SLICES - 20))\n        for _ in range(num_samples)\n    ]\n    print(f\"\\n[INFO] Échantillonnage de {num_samples} tranches pour estimer E_label...\")\n\n    for image_id, z in samples:\n        label_path = f\"{PATH_LABEL_ROOT}{image_id}.tif\"\n        try:\n            volume = iio.imread(label_path)\n            mask = (volume[z] == 1)\n            if mask.sum() == 0:\n                continue\n            dist = distance_transform_edt(mask)\n            thickness = dist[mask].max() * 2\n            thicknesses.append(thickness)\n        except Exception:\n            continue\n\n    if len(thicknesses) == 0:\n        print(\"[WARN] Aucun label détecté → E_label = 1 par défaut\")\n        return 1.0\n\n    E_label_mean = np.mean(thicknesses)\n    print(f\"[INFO] Épaisseur moyenne estimée = {E_label_mean:.2f} px\")\n    return E_label_mean\n\nE_LABEL = compute_average_ink_thickness()\nBASE_RADIUS = int(E_LABEL / 2)\nprint(f\"[INFO] Rayon de base (Base Radius) = {BASE_RADIUS} px\")\n\n# =========================================================\n# 2. Détection bright dots avec seuil ADAPTATIF\n# =========================================================\ndef detect_bright_dots(ct_slice):    \n    \"\"\"Applique le seuil quantile ultra-sélectif en utilisant le global QUANTILE_THRESHOLD.\"\"\"\n    thr = np.quantile(ct_slice, QUANTILE_THRESHOLD)\n    mask = (ct_slice >= thr).astype(np.uint8)\n    return mask, thr\n\n# =========================================================\n# 3. Dilatation dynamique selon le nombre de pixels détectés\n# =========================================================\ndef dynamic_dilation(mask, base_radius):\n    n = mask.sum()\n\n    if n > 10000:\n        r = max(1, base_radius - 4)\n    elif n > 2000:\n        r = max(1, base_radius - 2)\n    elif n > 200:\n        r = base_radius\n    else:\n        r = base_radius + 2 \n\n    se = disk(r)\n    dil = binary_dilation(mask, se).astype(np.uint8)\n\n    return dil, r\n\n# =========================================================\n# 4. Projection 3D et calcul IoU (OPTIMISÉ SANS VISUALISATION)\n# =========================================================\ndef run_projection_3d(image_id):\n    # Charge les volumes, retourne 0.0 si erreur de chargement.\n    try:\n        ct_vol = iio.imread(f\"{PATH_VOLUME_ROOT}{image_id}.tif\")\n        label_vol = iio.imread(f\"{PATH_LABEL_ROOT}{image_id}.tif\")\n    except Exception:\n        return 0.0\n\n    true_mask = (label_vol[REF_Z] == 1).astype(np.uint8)\n    final_mask = np.zeros_like(true_mask)\n    \n    # Projection 3D (OR)\n    for dz in Z_SHIFTS:\n        z = REF_Z + dz\n        if 0 <= z < TOTAL_SLICES:\n            ct_slice = ct_vol[z]\n            \n            # 1. Filtre Gaussien (utilise GAUSSIAN_SIGMA)\n            ct_slice_filtered = gaussian_filter(ct_slice, sigma=GAUSSIAN_SIGMA)\n            \n            # 2. Seuil adaptatif (utilise QUANTILE_THRESHOLD)\n            dots, thr = detect_bright_dots(ct_slice_filtered)\n\n            # 3. Dilatation dynamique\n            dilated, used_r = dynamic_dilation(dots, BASE_RADIUS)\n\n            # OR des masques\n            final_mask = (final_mask | dilated).astype(np.uint8)\n    \n    # SCORE IoU\n    union = np.sum((final_mask | true_mask) > 0)\n    inter = np.sum((final_mask & true_mask) > 0)\n    iou = inter / union if union > 0 else 0\n\n    return iou\n\n# =========================================================\n# 5. EXECUTION DE LA GRID SEARCH (Phase Finale - 16 Tests ciblés)\n# =========================================================\n# Paramètres de la grille\nQUANTILE_TESTS = [0.9980, 0.9985, 0.9990, 0.9993]  \nSIGMA_TESTS = [0.7, 0.8, 0.9, 1.0] # NOUVEAU: Exploration sous 0.8\n\n# Variables de suivi\nbest_iou = 0.0\nbest_params = {}\nall_results = []\n\nprint(\"\\n--- DÉBUT DE LA GRILLE DE RECHERCHE (PHASE FINALE : CIBLAGE OPTIMAL) ---\")\nprint(f\"Test de {len(QUANTILE_TESTS) * len(SIGMA_TESTS)} combinaisons sur {len(train_ids_all)} images.\")\n\nfor q_test in QUANTILE_TESTS:\n    # Mise à jour de la variable globale\n    QUANTILE_THRESHOLD = q_test \n    \n    for s_test in SIGMA_TESTS:\n        GAUSSIAN_SIGMA = s_test\n        \n        print(f\"\\n[TEST] Quantile={q_test:.4f}, Sigma={s_test:.1f}\")\n\n        # Exécution sur tous les IDs\n        all_ious = []\n        for img_id in train_ids_all:\n             iou_score = run_projection_3d(img_id) \n             all_ious.append(iou_score)\n        \n        # Calcul du score moyen\n        mean_iou = np.mean(all_ious)\n        \n        # Sauvegarde et comparaison\n        current_result = {\n            'Quantile': q_test,\n            'Sigma': s_test,\n            'Mean_IoU': mean_iou\n        }\n        all_results.append(current_result)\n        \n        print(f\"  → IOU MOYEN pour ce test : {mean_iou:.4f}\")\n        \n        if mean_iou > best_iou:\n            best_iou = mean_iou\n            best_params = current_result\n\n# =========================================================\n# Affichage du Meilleur Résultat\n# =========================================================\nprint(\"\\n--- RÉSULTATS FINAUX DU GRID SEARCH ---\")\nfor res in all_results:\n    print(f\"Q={res['Quantile']:.4f}, S={res['Sigma']:.1f} : IoU={res['Mean_IoU']:.4f}\")\n\nprint(\"\\n=======================================================\")\nprint(f\"🥇 MEILLEUR SCORE PROTOTYPE (Final Baseline): {best_iou:.4f}\")\nprint(f\"  avec Paramètres: Quantile={best_params['Quantile']:.4f}, Sigma={best_params['Sigma']:.1f}\")\nprint(\"=======================================================\")","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}