{"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":15062069,"sourceType":"competition"}],"dockerImageVersionId":31192,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import os\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nfrom PIL import Image\nimport tifffile\nimport cv2\nimport warnings\nwarnings.filterwarnings('ignore')\n\nBASE_PATH = \"/kaggle/input/vesuvius-challenge-surface-detection\"\nTRAIN_CSV = os.path.join(BASE_PATH, \"train.csv\")\nTEST_CSV = os.path.join(BASE_PATH, \"test.csv\")\nTRAIN_IMAGES_DIR = os.path.join(BASE_PATH, \"train_images\")\nTEST_IMAGES_DIR = os.path.join(BASE_PATH, \"test_images\")\nTRAIN_LABELS_DIR = os.path.join(BASE_PATH, \"train_labels\")\n\ndef load_tiff_volume(path):\n    try:\n        if not os.path.exists(path):\n            return None\n        \n        with Image.open(path) as img:\n            if img.n_frames == 1:\n                volume = np.array(img)\n                if volume.ndim == 2:\n                    volume = volume[np.newaxis, ...]\n            else:\n                slices = []\n                for i in range(img.n_frames):\n                    img.seek(i)\n                    slices.append(np.array(img))\n                volume = np.stack(slices, axis=0)\n        \n        return volume\n    except Exception:\n        return None\n\ndef get_available_files():\n    train_files = {}\n    for f in os.listdir(TRAIN_IMAGES_DIR):\n        if f.endswith('.tif'):\n            vol_id = f.replace('.tif', '')\n            train_files[vol_id] = {\n                'image_path': os.path.join(TRAIN_IMAGES_DIR, f),\n                'label_path': os.path.join(TRAIN_LABELS_DIR, f)\n            }\n    return train_files\n\ndef analyze_volume(volume, vol_id):\n    if volume is None:\n        return None\n    \n    depth, height, width = volume.shape\n    \n    fig = plt.figure(figsize=(20, 12))\n    \n    gs = plt.GridSpec(3, 5, figure=fig, hspace=0.3, wspace=0.3)\n    \n    slice_indices = [0, depth//4, depth//2, 3*depth//4, depth-1]\n    \n    for i, idx in enumerate(slice_indices):\n        ax = fig.add_subplot(gs[0, i])\n        slice_data = volume[idx]\n        im = ax.imshow(slice_data, cmap='gray', aspect='equal')\n        ax.set_title(f'Slice {idx}\\nMin: {slice_data.min():.1f}\\nMax: {slice_data.max():.1f}', fontsize=9)\n        ax.axis('off')\n        plt.colorbar(im, ax=ax, fraction=0.046, pad=0.04)\n    \n    ax_hist = fig.add_subplot(gs[1, :3])\n    ax_hist.hist(volume.flatten(), bins=100, alpha=0.7, color='steelblue', edgecolor='black')\n    ax_hist.set_title('Intensity Distribution', fontsize=12)\n    ax_hist.set_xlabel('Intensity')\n    ax_hist.set_ylabel('Frequency')\n    ax_hist.grid(True, alpha=0.3)\n    \n    ax_stats = fig.add_subplot(gs[1, 3:])\n    ax_stats.axis('off')\n    stats_text = f\"\"\"\n    Volume Statistics:\n    Shape: {volume.shape}\n    Data Type: {volume.dtype}\n    Min Value: {volume.min():.2f}\n    Max Value: {volume.max():.2f}\n    Mean: {volume.mean():.2f}\n    Std: {volume.std():.2f}\n    Size: {volume.nbytes / (1024**3):.2f} GB\n    \"\"\"\n    ax_stats.text(0.1, 0.5, stats_text, fontsize=11, family='monospace',\n                 verticalalignment='center', bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))\n    \n    projections = [\n        ('XY (max)', np.max(volume, axis=0), gs[2, 0]),\n        ('XZ (max)', np.max(volume, axis=1), gs[2, 1]),\n        ('YZ (max)', np.max(volume, axis=2), gs[2, 2]),\n        ('XY (mean)', np.mean(volume, axis=0), gs[2, 3]),\n        ('YZ (std)', np.std(volume, axis=2), gs[2, 4])\n    ]\n    \n    for i, (title, proj, pos) in enumerate(projections):\n        ax = fig.add_subplot(pos)\n        im = ax.imshow(proj, cmap='viridis', aspect='auto')\n        ax.set_title(title, fontsize=10)\n        ax.axis('off')\n        plt.colorbar(im, ax=ax, fraction=0.046, pad=0.04)\n    \n    plt.suptitle(f'3D Volume Analysis: {vol_id}', fontsize=16, y=1.02)\n    plt.tight_layout()\n    plt.show()\n    \n    return {\n        'shape': volume.shape,\n        'dtype': str(volume.dtype),\n        'min': float(volume.min()),\n        'max': float(volume.max()),\n        'mean': float(volume.mean()),\n        'std': float(volume.std())\n    }\n\ndef analyze_mask(mask, vol_id):\n    if mask is None:\n        return None\n    \n    unique_values = np.unique(mask)\n    \n    fig, axes = plt.subplots(2, 3, figsize=(15, 10))\n    \n    slice_idx = mask.shape[0] // 2\n    mask_slice = mask[slice_idx]\n    \n    axes[0, 0].imshow(mask_slice, cmap='tab10', vmin=0, vmax=2)\n    axes[0, 0].set_title(f'Mask Slice {slice_idx}')\n    axes[0, 0].axis('off')\n    \n    class_names = {0: 'Background', 1: 'Foreground', 2: 'Unlabeled'}\n    colors = ['lightgray', 'red', 'orange']\n    \n    total_voxels = mask.size\n    class_counts = []\n    class_percents = []\n    \n    for class_id in unique_values:\n        count = np.sum(mask == class_id)\n        percent = (count / total_voxels) * 100\n        class_counts.append(count)\n        class_percents.append(percent)\n    \n    axes[0, 1].bar(range(len(unique_values)), class_counts, color=colors[:len(unique_values)])\n    axes[0, 1].set_title('Class Distribution (Count)')\n    axes[0, 1].set_xticks(range(len(unique_values)))\n    axes[0, 1].set_xticklabels([class_names.get(i, f'Class {i}') for i in unique_values])\n    axes[0, 1].tick_params(axis='x', rotation=45)\n    \n    axes[0, 2].bar(range(len(unique_values)), class_percents, color=colors[:len(unique_values)])\n    axes[0, 2].set_title('Class Distribution (%)')\n    axes[0, 2].set_ylabel('Percentage')\n    axes[0, 2].set_xticks(range(len(unique_values)))\n    axes[0, 2].set_xticklabels([class_names.get(i, f'Class {i}') for i in unique_values])\n    axes[0, 2].tick_params(axis='x', rotation=45)\n    \n    foreground = (mask == 1).astype(np.uint8)\n    \n    projections = [\n        ('XY Foreground', np.max(foreground, axis=0), axes[1, 0]),\n        ('XZ Foreground', np.max(foreground, axis=1), axes[1, 1]),\n        ('YZ Foreground', np.max(foreground, axis=2), axes[1, 2])\n    ]\n    \n    for title, proj, ax in projections:\n        ax.imshow(proj, cmap='hot')\n        ax.set_title(title)\n        ax.axis('off')\n    \n    plt.suptitle(f'Mask Analysis: {vol_id}', fontsize=16)\n    plt.tight_layout()\n    plt.show()\n    \n    stats = {}\n    for class_id in unique_values:\n        count = np.sum(mask == class_id)\n        percent = (count / total_voxels) * 100\n        stats[f'class_{class_id}'] = {'count': count, 'percent': percent}\n    \n    return stats\n\ndef analyze_texture_features(volume, vol_id):\n    if volume is None:\n        return None\n    \n    slice_idx = volume.shape[0] // 2\n    slice_data = volume[slice_idx]\n    \n    slice_norm = cv2.normalize(slice_data, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n    \n    sobel_x = cv2.Sobel(slice_norm, cv2.CV_64F, 1, 0, ksize=3)\n    sobel_y = cv2.Sobel(slice_norm, cv2.CV_64F, 0, 1, ksize=3)\n    sobel_mag = np.sqrt(sobel_x**2 + sobel_y**2)\n    \n    laplacian = cv2.Laplacian(slice_norm, cv2.CV_64F)\n    \n    edges = cv2.Canny(slice_norm, 50, 150)\n    \n    _, binary = cv2.threshold(slice_norm, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)\n    \n    fig, axes = plt.subplots(3, 3, figsize=(15, 12))\n    \n    visualizations = [\n        ('Original Slice', slice_data, 'gray', axes[0, 0]),\n        ('Sobel Gradient', sobel_mag, 'hot', axes[0, 1]),\n        ('Laplacian', np.abs(laplacian), 'coolwarm', axes[0, 2]),\n        ('Canny Edges', edges, 'gray', axes[1, 0]),\n        ('Otsu Binary', binary, 'gray', axes[1, 1]),\n        ('Intensity Profile', slice_data[slice_data.shape[0]//2, :], 'line', axes[1, 2]),\n        ('Local Std Dev', cv2.blur(slice_norm.astype(np.float32), (5, 5)).std(), 'viridis', axes[2, 0]),\n        ('Gradient Magnitude', np.sqrt(cv2.Sobel(slice_norm, cv2.CV_64F, 1, 1, ksize=3)**2), 'plasma', axes[2, 1]),\n        ('Histogram Equalized', cv2.equalizeHist(slice_norm), 'gray', axes[2, 2])\n    ]\n    \n    for title, data, cmap, ax in visualizations:\n        ax.clear()\n        if title == 'Intensity Profile':\n            ax.plot(data, color='darkblue', linewidth=2)\n            ax.set_title(title)\n            ax.set_xlabel('X Position')\n            ax.set_ylabel('Intensity')\n            ax.grid(True, alpha=0.3)\n        elif title == 'Local Std Dev':\n            ax.imshow(np.ones((10, 10)) * data, cmap=cmap)\n            ax.set_title(f'{title}: {data:.2f}')\n            ax.axis('off')\n        else:\n            ax.imshow(data, cmap=cmap)\n            ax.set_title(title)\n            ax.axis('off')\n    \n    plt.suptitle(f'Texture Analysis: {vol_id}', fontsize=16)\n    plt.tight_layout()\n    plt.show()\n\ndef compare_volumes(vol_ids, train_files, n_volumes=3):\n    fig, axes = plt.subplots(n_volumes, 5, figsize=(20, 4*n_volumes))\n    \n    if n_volumes == 1:\n        axes = axes.reshape(1, -1)\n    \n    for i, vol_id in enumerate(vol_ids[:n_volumes]):\n        if vol_id not in train_files:\n            continue\n        \n        volume = load_tiff_volume(train_files[vol_id]['image_path'])\n        mask = load_tiff_volume(train_files[vol_id]['label_path'])\n        \n        if volume is None or mask is None:\n            continue\n        \n        mid_slice = volume.shape[0] // 2\n        \n        axes[i, 0].imshow(volume[mid_slice], cmap='gray')\n        axes[i, 0].set_title(f'{vol_id}\\nVolume')\n        axes[i, 0].axis('off')\n        \n        axes[i, 1].hist(volume.flatten(), bins=50, alpha=0.7, color='steelblue')\n        axes[i, 1].set_title('Intensity Hist')\n        axes[i, 1].grid(True, alpha=0.3)\n        \n        axes[i, 2].imshow(mask[mid_slice], cmap='tab10', vmin=0, vmax=2)\n        axes[i, 2].set_title('Mask')\n        axes[i, 2].axis('off')\n        \n        unique_mask_vals = np.unique(mask)\n        mask_colors = ['lightgray', 'red', 'orange']\n        \n        for idx, class_id in enumerate(unique_mask_vals):\n            count = np.sum(mask == class_id)\n            axes[i, 3].bar(idx, count, color=mask_colors[class_id], alpha=0.7)\n        \n        axes[i, 3].set_title('Class Counts')\n        axes[i, 3].set_xticks(range(len(unique_mask_vals)))\n        axes[i, 3].set_xticklabels([str(v) for v in unique_mask_vals])\n        \n        foreground_mask = (mask == 1).astype(np.uint8)\n        if foreground_mask.any():\n            foreground_voxels = volume[foreground_mask == 1]\n            axes[i, 4].hist(foreground_voxels.flatten(), bins=50, alpha=0.7, color='red')\n            axes[i, 4].set_title('Foreground Intensity')\n            axes[i, 4].grid(True, alpha=0.3)\n        else:\n            axes[i, 4].text(0.5, 0.5, 'No\\nForeground', \n                           ha='center', va='center', transform=axes[i, 4].transAxes)\n            axes[i, 4].axis('off')\n    \n    plt.suptitle('Volume Comparison', fontsize=16, y=1.02)\n    plt.tight_layout()\n    plt.show()\n\ndef create_summary_report(train_df, train_files):\n    available_ids = list(train_files.keys())\n    \n    print(f\"\\n📊 Dataset Overview:\")\n    print(f\"   • Total entries in train.csv: {len(train_df)}\")\n    print(f\"   • Available 3D volumes: {len(available_ids)}\")\n    print(f\"   • Unique scrolls: {train_df['scroll_id'].nunique()}\")\n    \n    scroll_counts = train_df[train_df['id'].astype(str).isin(available_ids)]['scroll_id'].value_counts()\n    print(f\"\\n📜 Scroll Distribution:\")\n    for scroll_id, count in scroll_counts.items():\n        print(f\"   • Scroll {scroll_id}: {count} volumes\")\n    \n    print(f\"\\n🔍 Data Quality Check:\")\n    all_stats = []\n    \n    for i, vol_id in enumerate(available_ids[:5]):\n        volume = load_tiff_volume(train_files[vol_id]['image_path'])\n        mask = load_tiff_volume(train_files[vol_id]['label_path'])\n        \n        if volume is not None:\n            vol_stats = analyze_volume(volume, vol_id)\n            if vol_stats:\n                all_stats.append(vol_stats)\n            \n            if mask is not None:\n                analyze_mask(mask, vol_id)\n            \n            analyze_texture_features(volume, vol_id)\n    \n    if all_stats:\n        print(f\"\\n📈 Volume Statistics (sample):\")\n        stats_df = pd.DataFrame(all_stats)\n        print(stats_df.describe().to_string())\n    \n    compare_volumes(available_ids[:3], train_files, n_volumes=3)\n\ndef main():\n    train_df = pd.read_csv(TRAIN_CSV)\n    test_df = pd.read_csv(TEST_CSV)\n    \n    train_files = get_available_files()\n    \n    print(f\"Found {len(train_files)} available volume files\")\n    \n    if train_files:\n        create_summary_report(train_df, train_files)\n    else:\n        print(\"No train files found. Check directory structure.\")\n\nif __name__ == \"__main__\":\n    main()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-26T18:08:48.657559Z","iopub.execute_input":"2026-01-26T18:08:48.657964Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}