{"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":[{"sourceType":"competition","sourceId":97984,"databundleVersionId":14096757},{"sourceType":"datasetVersion","sourceId":13746387,"datasetId":8747012,"databundleVersionId":14496454},{"sourceType":"datasetVersion","sourceId":15884760,"datasetId":10185101,"databundleVersionId":16838585},{"sourceType":"kernelVersion","sourceId":312109075}],"dockerImageVersionId":31153,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport cv2\nimport os\nimport torch\nfrom tabulate import tabulate\nfrom tqdm import tqdm\nimport warnings\nwarnings.filterwarnings('ignore')\n\n# Install cc3d silently\ntry:\n    import cc3d\nexcept:\n    import subprocess\n    subprocess.run([\n        'pip', 'install', 'connected-components-3d', '--no-index',\n        '--find-links=file:///kaggle/input/hengck23-submit-physionet/hengck23-submit-physionet/setup',\n        '-q'\n    ], capture_output=True)\n    import cc3d\nfrom scipy.signal import savgol_filter\n\nDEVICE = 'cuda' if torch.cuda.is_available() else 'cpu'\nFLOAT_TYPE = torch.float16\n\nBASE_PATH = '/kaggle/input/hengck23-submit-physionet/hengck23-submit-physionet'\nKAGGLE_DIR = '/kaggle/input/physionet-ecg-image-digitization'\nWEIGHT_DIR = f'{BASE_PATH}/weight'\nOUT_DIR = '/kaggle/working/output-combined'\n\nos.makedirs(OUT_DIR, exist_ok=True)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:20.892008Z","iopub.execute_input":"2026-04-26T14:46:20.892229Z","iopub.status.idle":"2026-04-26T14:46:33.043252Z","shell.execute_reply.started":"2026-04-26T14:46:20.892210Z","shell.execute_reply":"2026-04-26T14:46:33.042589Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Load model","metadata":{}},{"cell_type":"markdown","source":"## Hengck Model","metadata":{}},{"cell_type":"code","source":"import sys\nsys.path.insert(0, BASE_PATH)\n\nfrom stage0_model import Net as Stage0Net\nfrom stage0_common import load_net, image_to_batch, output_to_predict as s0_output, normalise_by_homography\n\nfrom stage1_model import Net as Stage1Net\nfrom stage1_common import load_net as s1_load, output_to_predict as s1_output, rectify_image\n\nfrom stage2_model import Net as Stage2Net\nfrom stage2_common import load_net as s2_load, pixel_to_series\n\n# Load models\nstage0 = load_net(Stage0Net(pretrained=False), f'{WEIGHT_DIR}/stage0-last.checkpoint.pth').to(DEVICE).eval()\nstage1 = s1_load(Stage1Net(pretrained=False), f'{WEIGHT_DIR}/stage1-last.checkpoint.pth').to(DEVICE).eval()\nstage2 = s2_load(Stage2Net(pretrained=False), f'{WEIGHT_DIR}/stage2-00005810.checkpoint.pth').to(DEVICE).eval()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:33.044911Z","iopub.execute_input":"2026-04-26T14:46:33.045488Z","iopub.status.idle":"2026-04-26T14:46:44.524679Z","shell.execute_reply.started":"2026-04-26T14:46:33.045470Z","shell.execute_reply":"2026-04-26T14:46:44.523785Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# !pip install segmentation-models-pytorch\n\n!pip install --no-index --find-links /kaggle/input/notebooks/ravnoorsingh101/smp-offline-wheels/smp-wheels segmentation-models-pytorch --no-deps","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:44.525542Z","iopub.execute_input":"2026-04-26T14:46:44.525889Z","iopub.status.idle":"2026-04-26T14:46:46.199536Z","shell.execute_reply.started":"2026-04-26T14:46:44.525870Z","shell.execute_reply":"2026-04-26T14:46:46.198804Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import torch\nimport segmentation_models_pytorch as smp\n\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n\nmodel = smp.Unet(\n    encoder_name=\"efficientnet-b0\",\n    encoder_weights=None,\n    in_channels=3,\n    classes=1,\n    activation=None,\n    decoder_attention_type=\"scse\"\n)\n\n# Load your trained weights\nmodel.load_state_dict(torch.load('/kaggle/input/datasets/ravnoor000/physionet-unet-weights-2/best_attention_unet (4).pth', map_location=device))\nmodel = model.to(device)\nmodel.eval()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:46.200747Z","iopub.execute_input":"2026-04-26T14:46:46.201070Z","iopub.status.idle":"2026-04-26T14:46:46.675327Z","shell.execute_reply.started":"2026-04-26T14:46:46.201046Z","shell.execute_reply":"2026-04-26T14:46:46.674487Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_visual_mask(logits):\n    return (torch.sigmoid(logits) > 0.1).float().cpu().numpy()[0,0]\n\ndef extract_1d_signal(mask):\n    h, w = mask.shape\n    signal = np.zeros(w)\n    for x in range(w):\n        ys = np.where(mask[:,x] > 0)[0]\n        signal[x] = np.mean(ys) if len(ys) else np.nan\n    return np.nan_to_num(signal, nan=np.nanmedian(signal))\n\ndef pixels_to_mv(signal, baseline, ppmv):\n    return (baseline - signal) / ppmv","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:46.676179Z","iopub.execute_input":"2026-04-26T14:46:46.676478Z","iopub.status.idle":"2026-04-26T14:46:46.682167Z","shell.execute_reply.started":"2026-04-26T14:46:46.676455Z","shell.execute_reply":"2026-04-26T14:46:46.681503Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Main pipeline ","metadata":{}},{"cell_type":"code","source":"# test_df = pd.read_csv(f'{KAGGLE_DIR}/test.csv')\n# test_df['id'] = test_df['id'].astype(str)\n\n# submission_data = []\n\n# for sample_id in tqdm(test_df['id'].unique()):\n\n#     # =========================\n#     # STAGE 0 + 1 (hengck)\n#     # =========================\n#     image = cv2.imread(f'{KAGGLE_DIR}/test/{sample_id}.png')\n#     image = cv2.cvtColor(image, cv2.COLOR_BGR2RGB)\n\n#     batch = image_to_batch(image)\n\n#     with torch.no_grad():\n#         out0 = stage0(batch)\n#         rotated, kp = s0_output(image, batch, out0)\n#         norm, kp, _ = normalise_by_homography(rotated, kp)\n\n#         batch1 = {'image': torch.from_numpy(norm.transpose(2,0,1)).unsqueeze(0)}\n#         out1 = stage1(batch1)\n#         grid, _ = s1_output(norm, batch1, out1)\n#         rect = rectify_image(norm, grid)\n\n#     # =========================\n#     # HENGCK SIGNAL\n#     # =========================\n#     batch2 = {'image': torch.from_numpy(rect.transpose(2,0,1)).unsqueeze(0)}\n\n#     with torch.no_grad():\n#         out2 = stage2(batch2)\n\n#     pixel = out2['pixel'].cpu().numpy()[0]\n\n#     zero_mv = [703.5, 987.5, 1271.5, 1531.5]\n#     mv_to_pixel = 80.0\n\n#     hengck_series = (np.array(zero_mv).reshape(4,1) - pixel[:,118:2080]) / mv_to_pixel\n\n#     # =========================\n#     # YOUR MODEL SIGNAL\n#     # =========================\n#     img_resized = cv2.resize(rect, (1024, 256))\n#     img_tensor = torch.from_numpy(img_resized.transpose(2,0,1)).unsqueeze(0).float().to(DEVICE)\n\n#     with torch.no_grad():\n#         logits = your_model(img_tensor)\n\n#     mask = get_visual_mask(logits)\n#     signal_px = extract_1d_signal(mask)\n\n#     # scale back\n#     signal_px = signal_px * (rect.shape[0] / mask.shape[0])\n\n#     baseline = rect.shape[0] // 2\n#     your_signal = pixels_to_mv(signal_px, baseline, 80.0)\n\n#     # =========================\n#     # COMBINE\n#     # =========================\n#     alpha = 0.6  # hengck weight\n\n#     combined = alpha * hengck_series[3] + (1 - alpha) * your_signal\n\n#     # =========================\n#     # SAVE\n#     # =========================\n#     df = test_df[test_df['id'] == sample_id]\n\n#     for _, row in df.iterrows():\n#         target_len = row['number_of_rows']\n\n#         interp = np.interp(\n#             np.linspace(0,1,target_len),\n#             np.linspace(0,1,len(combined)),\n#             combined\n#         )\n\n#         for t, val in enumerate(interp):\n#             submission_data.append({\n#                 'id': f\"{sample_id}_{t}_{row['lead']}\",\n#                 'value': float(val)\n#             })\n\n# submission = pd.DataFrame(submission_data)\n# submission.to_csv('submission.csv', index=False)\n# print(\"✅ Submission ready\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:46.682944Z","iopub.execute_input":"2026-04-26T14:46:46.683199Z","iopub.status.idle":"2026-04-26T14:46:46.696205Z","shell.execute_reply.started":"2026-04-26T14:46:46.683178Z","shell.execute_reply":"2026-04-26T14:46:46.695309Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Experiment - 2 ","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport cv2\nfrom glob import glob\nimport matplotlib.pyplot as plt\nfrom collections import defaultdict\nfrom tqdm import tqdm\nfrom scipy.signal import medfilt\n\nfrom sklearn.model_selection import KFold\nfrom sklearn.metrics import r2_score\n    \nfrom tensorflow.keras.layers import Input, Dense, Activation, Reshape, GaussianNoise\nfrom tensorflow.keras.initializers import Constant\nfrom tensorflow.keras.models import Model\nfrom tensorflow.keras.optimizers import Adam\nfrom tensorflow.keras.losses import MeanSquaredError\nfrom tensorflow.keras.callbacks import EarlyStopping, ReduceLROnPlateau, TerminateOnNaN","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:46.698053Z","iopub.execute_input":"2026-04-26T14:46:46.698350Z","iopub.status.idle":"2026-04-26T14:46:59.504844Z","shell.execute_reply.started":"2026-04-26T14:46:46.698326Z","shell.execute_reply":"2026-04-26T14:46:59.504217Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# import cv2\n# import numpy as np\n# import matplotlib.pyplot as plt\n\n# class MarkerFinder:\n#     \"\"\"This class finds the 13 markers in scanned ecg images and guesses the 4 line ends.\"\"\"\n    \n#     def __init__(self, show_templates=False):\n#         # Derive the templates from type 1 images\n#         # np.max keeps the gridlines and markers and removes the ecg lines\n#         ima = np.max([\n#             cv2.imread('/kaggle/input/physionet-ecg-image-digitization/train/1006427285/1006427285-0001.png'),\n#             cv2.imread('/kaggle/input/physionet-ecg-image-digitization/train/102150619/102150619-0001.png'),\n#             cv2.imread('/kaggle/input/physionet-ecg-image-digitization/train/1079294623/1079294623-0001.png'),\n#         ], axis=0)\n\n#         # Template points in global coordinates of type 1 images\n#         absolute_points = np.zeros((17, 2), dtype=int)\n#         for i in range(3):\n#             absolute_points[5 * i] = np.array([707 + 284 * i, 118]) # y, x\n#             for j in range(1, 5):\n#                 absolute_points[5 * i + j] = np.array([707 + 284 * i, 118 + 492 * j])\n#         absolute_points[5 * 3] = np.array([1535, 118])\n#         absolute_points[5 * 3 + 1] = np.array([1535, 118 + 492 * 4])\n\n#         # Top left corner of template rectangle\n#         template_positions = [None] * 17\n#         for i in range(len(absolute_points)):\n#             if absolute_points[i][1] < 118 + 492 * 4:\n#                 if i % 5 == 0:\n#                     template_positions[i] = (absolute_points[i][0] - 87, absolute_points[i][1] - 50) # y, x\n#                 else:\n#                     template_positions[i] = (absolute_points[i][0] - 37, absolute_points[i][1] - 13)\n\n#         # Height and width of the templates\n#         template_sizes = np.array([(105, 60)] * 17) # height, width\n\n#         # Transform the points to relative coordinates (inside the template)\n#         template_points = [np.array([absolute_points[i][0] - template_positions[i][0],\n#                                      absolute_points[i][1] - template_positions[i][1]])\n#                            if template_positions[i] is not None\n#                            else None\n#                            for i in range(len(absolute_points))]\n\n#         # Save the template matrices\n#         templates = [None] * 17\n#         for i in range(len(template_positions)):\n#             if template_points[i] is not None:\n#                 template = (ima[template_positions[i][0]:template_positions[i][0]+template_sizes[i][0],\n#                             template_positions[i][1]:template_positions[i][1]+template_sizes[i][1]])\n#                 templates[i] = template\n\n#         if show_templates:\n#             _, axs = plt.subplots(4, 4, figsize=(5, 7))\n#             for i in range(len(template_positions)):\n#                 if template_points[i] is not None:\n#                     template = templates[i].copy()\n#                     cv2.rectangle(template,\n#                                   (template_points[i][1]-1, template_points[i][0]-1),\n#                                   (template_points[i][1]+1, template_points[i][0]+1), \n#                                   [255, 0, 0], 2)\n#                     axs[i // 5, i % 5].imshow(template)\n#             for i in range(13, len(axs.ravel())):\n#                 axs.ravel()[i].axis('off')\n#             plt.tight_layout()\n#             plt.suptitle('The templates for the 13 markers', y=1.01)\n#             plt.show()\n\n#         self._absolute_points = absolute_points\n#         self._template_positions = template_positions\n#         self._template_sizes = template_sizes\n#         self._template_points = template_points\n#         self._templates = templates\n        \n#     def find_markers(self, ima, warn=False, plot=False, title=''):\n#         \"\"\"Return 17 markers as list of size-2 integer arrays (row, column)\"\"\"\n#         if ima.shape[0] != 1652:\n#             # For this pipeline, we bypass the error since we might process different sizes, \n#             # but we scale it to expected size for the template matching to work\n#             pass\n\n#         markers = np.full((17, 2), -1)\n\n#         # Find 13 template-based markers\n#         for j in range(len(self._templates)):\n#             if self._template_points[j] is not None:\n#                 t = self._template_positions[j][0]-100\n#                 l = max(self._template_positions[j][1]-100, 0)\n                \n#                 # Safety check to ensure search range doesn't go out of bounds\n#                 bottom_bound = min(ima.shape[0], self._template_positions[j][0]+100+self._template_sizes[j][0])\n#                 right_bound = min(ima.shape[1], self._template_positions[j][1]+250+self._template_sizes[j][0])\n                \n#                 search_range = ima[t:bottom_bound, l:right_bound]\n                \n#                 if search_range.shape[0] < self._templates[j].shape[0] or search_range.shape[1] < self._templates[j].shape[1]:\n#                     continue # Skip if the image is too small for the template\n                \n#                 res = cv2.matchTemplate(search_range, self._templates[j], cv2.TM_CCOEFF)\n#                 min_val, max_val, min_loc, max_loc = cv2.minMaxLoc(res)\n    \n#                 top_left = max_loc\n#                 markers[j] = np.array((t + top_left[1] + self._template_points[j][0],\n#                                        l + top_left[0] + self._template_points[j][1]))\n\n#         # Guess the ends of the first three lines\n#         for i in range(3):\n#             if markers[5 * i + 3][0] != -1 and markers[5 * i + 2][0] != -1:\n#                 m = markers[5 * i + 3] * 2 - markers[5 * i + 2]\n#                 markers[5 * i + 4] = m\n\n#         # Guess the end of the fourth line\n#         if markers[14][0] != -1 and markers[9][0] != -1:\n#             markers[16] = ((markers[14] * (284 + 260) - markers[9] * 260) / 284).astype(int)\n\n#         return markers\n        \n#     @staticmethod\n#     def lead_info(lead):\n#         \"\"\"Specify which markers mark the begin and the end of a lead.\"\"\"\n#         begin, end = {\n#             'I': (0, 1),\n#             'II-subset': (5, 6),\n#             'III': (10, 11),\n#             'aVR': (1, 2),\n#             'aVL': (6, 7),\n#             'aVF': (11, 12),\n#             'V1': (2, 3),\n#             'V2': (7, 8),\n#             'V3': (12, 13),\n#             'V4': (3, 4),\n#             'V5': (8, 9),\n#             'V6': (13, 14),\n#             'II': (15, 16), # The long rhythm strip at the bottom\n#         }[lead]\n        \n#         # We add a fallback for the long II lead to map it to the 3x4 grid for our UNet\n#         row_idx = begin // 5\n#         if lead == 'II':\n#             row_idx = 3 # 4th row for rhythm strip\n            \n#         return row_idx, begin, end\n\n# # --- INSTANTIATE THE OBJECT HERE ---\n# print(\"Initializing MarkerFinder...\")\n# mf = MarkerFinder(show_templates=False)\n# print(\"MarkerFinder ready!\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:59.505613Z","iopub.execute_input":"2026-04-26T14:46:59.506207Z","iopub.status.idle":"2026-04-26T14:46:59.512756Z","shell.execute_reply.started":"2026-04-26T14:46:59.506179Z","shell.execute_reply":"2026-04-26T14:46:59.512060Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# import albumentations as A \n# from albumentations.pytorch import ToTensorV2\n\n# # ================================================================================\n# # CREATE HYBRID SUBMISSION (Hengck23 + Attention U-Net)\n# # ================================================================================\n# print(\"\\n\" + \"=\" * 80)\n# print(\"CREATING HYBRID SUBMISSION\")\n# print(\"=\" * 80)\n\n# # 1. Ensure your U-Net is loaded and ready\n# model.eval()\n# device = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n# test_transforms = A.Compose([\n#     A.Resize(256, 512),\n#     A.Normalize(mean=(0.485, 0.456, 0.406), std=(0.229, 0.224, 0.225)),\n#     ToTensorV2()\n# ])\n\n# submission_data = []\n\n# test_df = pd.read_csv('/kaggle/input/physionet-ecg-image-digitization/test.csv')\n# gb = test_df.groupby('id')\n\n# for rec_idx, (sample_id, df) in enumerate(tqdm(gb, desc=\"Building submission\")):\n    \n#     # ---------------------------------------------------------\n#     # PART A: Get the Hengck23 Prediction\n#     # ---------------------------------------------------------\n#     hengck_success = False\n#     series_by_lead_hengck = {}\n#     try:\n#         series = np.load(f'{OUT_DIR}/digitalised/{sample_id}.series.npy')\n#         _4_, L = series.shape\n        \n#         for l in range(3):\n#             lead_names = [\n#                 ['I', 'aVR', 'V1', 'V4'],\n#                 ['II', 'aVL', 'V2', 'V5'],\n#                 ['III', 'aVF', 'V3', 'V6'],\n#             ][l]\n#             index = [int(round(1*L/4)), int(round(2*L/4)), int(round(3*L/4))]\n#             split = np.split(series[l], index)\n#             for k, s in zip(lead_names, split):\n#                 series_by_lead_hengck[k] = s\n#         series_by_lead_hengck['II'] = series[3]\n#         hengck_success = True\n#     except:\n#         pass # Will fall back to UNet entirely\n\n#     # ---------------------------------------------------------\n#     # PART B: Get Your Attention U-Net Prediction\n#     # ---------------------------------------------------------\n#     series_by_lead_unet = {}\n    \n#     # Load the original image for the UNet\n#     img_path = os.path.join(KAGGLE_DIR, 'test', f\"{sample_id}.png\")\n#     full_img = cv2.imread(img_path)\n#     if full_img is not None:\n#         full_img_rgb = cv2.cvtColor(full_img, cv2.COLOR_BGR2RGB)\n#         markers = mf.find_markers(full_img) # Using your MarkerFinder!\n        \n#         for lead_name in ['I', 'II', 'III', 'aVR', 'aVL', 'aVF', 'V1', 'V2', 'V3', 'V4', 'V5', 'V6']:\n#             try:\n#                 row_idx, begin_idx, end_idx = mf.lead_info(lead_name)\n#                 left_x, right_x = markers[begin_idx][1], markers[end_idx][1]\n                \n#                 padding_top = 130\n#                 top_y = markers[begin_idx][0] - padding_top \n                \n#                 if row_idx < 2:\n#                     bottom_y = markers[begin_idx + 5][0] - padding_top\n#                 elif row_idx == 2:\n#                     bottom_y = markers[15][0] - padding_top\n#                 else:\n#                     bottom_y = top_y + 350\n                \n#                 top_y, bottom_y = max(0, top_y), min(full_img.shape[0], bottom_y)\n#                 left_x, right_x = max(0, left_x), min(full_img.shape[1], right_x)\n                \n#                 crop = full_img_rgb[top_y:bottom_y, left_x:right_x]\n                \n#                 # UNet Inference\n#                 input_tensor = test_transforms(image=crop)['image'].unsqueeze(0).to(device)\n#                 with torch.no_grad():\n#                     mask_logits = model(input_tensor)[0]\n                \n#                 mask = get_visual_mask(mask_logits, threshold=0.1) \n#                 pixel_signal = extract_1d_signal_from_mask(mask) \n                \n#                 actual_baseline_y = markers[begin_idx][0] - top_y\n#                 pixel_signal_scaled = pixel_signal * (crop.shape[0] / mask.shape[0])\n#                 signal_mv = (actual_baseline_y - pixel_signal_scaled) / 80.0\n                \n#                 series_by_lead_unet[lead_name] = signal_mv\n#             except Exception as e:\n#                 # If UNet fails to crop/infer, provide a blank array\n#                 series_by_lead_unet[lead_name] = np.zeros(2500)\n#     else:\n#         for lead in ['I', 'II', 'III', 'aVR', 'aVL', 'aVF', 'V1', 'V2', 'V3', 'V4', 'V5', 'V6']:\n#             series_by_lead_unet[lead] = np.zeros(2500)\n\n#     # ---------------------------------------------------------\n#     # PART C: The Mathematical Blend\n#     # ---------------------------------------------------------\n#     for _, d in df.iterrows():\n#         lead = d.lead\n#         target_len = int(d.number_of_rows)\n        \n#         # Get Hengck signal (if available)\n#         if hengck_success and lead in series_by_lead_hengck:\n#             s_hengck = series_by_lead_hengck[lead]\n#             x_old = np.linspace(0.0, 1.0, len(s_hengck), endpoint=False)\n#             x_new = np.linspace(0.0, 1.0, target_len, endpoint=False)\n#             s_hengck_interp = np.interp(x_new, x_old, s_hengck)\n#         else:\n#             s_hengck_interp = np.zeros(target_len)\n            \n#         # Get UNet signal\n#         s_unet = series_by_lead_unet[lead]\n#         x_old_u = np.linspace(0.0, 1.0, len(s_unet), endpoint=False)\n#         x_new_u = np.linspace(0.0, 1.0, target_len, endpoint=False)\n#         s_unet_interp = np.interp(x_new_u, x_old_u, s_unet)\n\n#         # ENSEMBLE LOGIC:\n#         # If Hengck succeeded, blend them (e.g., 70% Hengck, 30% UNet to smooth out noise).\n#         # If Hengck failed, rely 100% on UNet.\n#         if hengck_success:\n#             s_final = (1 * s_hengck_interp) + (0 * s_unet_interp)\n#         else:\n#             s_final = s_unet_interp\n            \n#         # Apply the Savitzky-Golay filter to the final blended signal to remove any micro-jaggedness\n#         s_final = savgol_filter(s_final, window_length=11, polyorder=3)\n#         s_final = s_final.astype(np.float32)\n        \n#         for t in range(target_len):\n#             submission_data.append({\n#                 'id': f'{sample_id}_{t}_{lead}',\n#                 'value': float(s_final[t])\n#             })\n\n# submission = pd.DataFrame(submission_data)\n# submission.to_csv('submission.csv', index=False)\n# print(f\"\\n✅ Hybrid Submission created: {len(submission):,} rows\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:59.513715Z","iopub.execute_input":"2026-04-26T14:46:59.514477Z","iopub.status.idle":"2026-04-26T14:46:59.545766Z","shell.execute_reply.started":"2026-04-26T14:46:59.514452Z","shell.execute_reply":"2026-04-26T14:46:59.545027Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Experiment - 3","metadata":{}},{"cell_type":"code","source":"\"\"\"\nECG Inference Pipeline — Improved Version\nKey improvements over baseline (SNR ~15):\n  1. Adaptive Savitzky-Golay smoothing (preserves QRS peaks, smooths flat regions)\n  2. Baseline wander removal\n  via high-pass Butterworth filter\n  3. 60 Hz notch filter for powerline noise\n  4. Einthoven's law enforcement (I - II + III = 0) as a post-processing correction\n  5. Smart lead-split boundary detection using cross-correlation energy minimization\n  6. Physiological amplitude clipping (±5 mV)\n  7. Per-lead z-score outlier suppression before submission\n\"\"\"\n\nprint(\"=\" * 80)\nprint(\"ECG INFERENCE — IMPROVED\")\nprint(\"=\" * 80)\n\nimport sys, os\nimport numpy as np\nimport pandas as pd\nimport cv2\nimport torch\nimport matplotlib.pyplot as plt\nfrom tqdm import tqdm\nfrom tabulate import tabulate\nfrom scipy.signal import savgol_filter, butter, filtfilt, iirnotch, correlate\n\n# ================================================================================\n# PATHS & CONFIG\n# ================================================================================\n\nbase_path = '/kaggle/input/hengck23-submit-physionet/hengck23-submit-physionet'\nsys.path.insert(0, base_path)\n\nKAGGLE_DIR  = '/kaggle/input/physionet-ecg-image-digitization'\nWEIGHT_DIR  = f'{base_path}/weight'\nOUT_DIR     = '/kaggle/working/output-submit'\nDEVICE      = 'cuda' if torch.cuda.is_available() else 'cpu'\nFLOAT_TYPE  = torch.float16\n\nos.makedirs(f'{OUT_DIR}/normalised',  exist_ok=True)\nos.makedirs(f'{OUT_DIR}/rectified',   exist_ok=True)\nos.makedirs(f'{OUT_DIR}/digitalised', exist_ok=True)\n\nprint(f\"\\n🔧 Device : {DEVICE}\")\nprint(f\"📁 Weights: {WEIGHT_DIR}\")\n\n# ECG physical constants\nMV_TO_PIXEL  = 80.0          # pixels per mV (from competition grid)\nMAX_AMP_MV   = 5.0           # physiological clip limit ±5 mV\nNOTCH_FREQ   = 60.0          # power-line frequency (Hz) — change to 50.0 for EU data\nNOTCH_Q      = 30.0\nBASELINE_CUT = 0.5           # high-pass cutoff (Hz) for baseline wander\n\n# ================================================================================\n# SIGNAL PROCESSING HELPERS\n# ================================================================================\n\ndef adaptive_savgol(signal, base_window=11, polyorder=3):\n    \"\"\"\n    Adaptive Savitzky-Golay: uses a wider smoothing window in flat/iso-electric\n    regions and a narrower window near high-variance peaks (QRS complexes).\n    This preserves sharp features while aggressively smoothing noise elsewhere.\n    \"\"\"\n    local_var = np.array([\n        np.var(signal[max(0, i - 5):i + 5]) for i in range(len(signal))\n    ])\n    norm_var = local_var / (local_var.max() + 1e-8)\n    # High variance → keep sharp (small window); low variance → smooth more\n    windows = (base_window + (1.0 - norm_var) * 10).astype(int)\n    windows = np.where(windows % 2 == 0, windows + 1, windows)   # must be odd\n    win = int(np.median(windows))\n    win = win if win % 2 == 1 else win + 1\n    win = max(5, win)\n    return savgol_filter(signal, window_length=win, polyorder=polyorder)\n\n\ndef remove_baseline(signal, fs=500, cutoff=0.5):\n    \"\"\"\n    High-pass Butterworth filter to remove low-frequency baseline wander.\n    Typical ECG baseline wander lives below 0.5 Hz.\n    \"\"\"\n    b, a = butter(3, cutoff / (fs / 2.0), btype='high')\n    return filtfilt(b, a, signal)\n\n\ndef notch_filter(signal, fs=500, freq=60.0, Q=30.0):\n    \"\"\"\n    IIR notch filter to suppress power-line interference (50 or 60 Hz).\n    \"\"\"\n    b, a = iirnotch(freq / (fs / 2.0), Q)\n    return filtfilt(b, a, signal)\n\n\ndef apply_signal_chain(signal, fs=500):\n    \"\"\"\n    Full signal-processing chain applied to every lead:\n      1. Remove baseline wander\n      2. Remove power-line noise\n      3. Adaptive smoothing\n      4. Physiological amplitude clipping\n    \"\"\"\n    # Need minimum length for filters (6× filter order)\n    if len(signal) < 30:\n        return signal\n\n    s = remove_baseline(signal, fs=fs, cutoff=BASELINE_CUT)\n    s = notch_filter(s,         fs=fs, freq=NOTCH_FREQ, Q=NOTCH_Q)\n    s = adaptive_savgol(s, base_window=11, polyorder=3)\n    s = np.clip(s, -MAX_AMP_MV, MAX_AMP_MV)   # physiological clipping\n    return s\n\n\ndef enforce_einthoven(lead_I, lead_II, lead_III):\n    \"\"\"\n    Enforce Einthoven's law: Lead_I - Lead_II + Lead_III = 0\n    (equivalently: Lead_II = Lead_I + Lead_III)\n\n    We project onto the constraint manifold by distributing the residual\n    equally across all three leads — minimum L2-norm correction.\n    \"\"\"\n    residual = lead_I - lead_II + lead_III\n    lead_I_c   = lead_I   - residual / 3.0\n    lead_II_c  = lead_II  + residual / 3.0\n    lead_III_c = lead_III - residual / 3.0\n    return lead_I_c, lead_II_c, lead_III_c\n\n\ndef find_split_boundaries(long_lead, n_splits=4):\n    \"\"\"\n    Use local signal energy to find the quietest point near each nominal\n    quarter boundary — ECG paper strips always have brief iso-electric\n    transitions at the split points.\n\n    Returns a list of (n_splits - 1) boundary indices.\n    \"\"\"\n    L = len(long_lead)\n    segment_len = L // n_splits\n    boundaries = []\n\n    for i in range(1, n_splits):\n        nominal = i * segment_len\n        search_range = max(10, segment_len // 10)\n        lo = max(0, nominal - search_range)\n        hi = min(L, nominal + search_range)\n        window = long_lead[lo:hi]\n\n        # Smooth energy of the first derivative\n        energy = np.convolve(\n            np.abs(np.diff(window, prepend=window[0])),\n            np.ones(20) / 20,\n            mode='same'\n        )\n        best = lo + int(np.argmin(energy))\n        boundaries.append(best)\n\n    return boundaries\n\n\n# ================================================================================\n# DATA INSPECTION\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"DATA INSPECTION\")\nprint(\"=\" * 80)\n\ntrain_df = pd.read_csv(f'{KAGGLE_DIR}/train.csv')\ntest_df  = pd.read_csv(f'{KAGGLE_DIR}/test.csv')\ntest_df['id'] = test_df['id'].astype(str)\nvalid_id = test_df['id'].unique().tolist()\n\nprint(\"\\n📊 DATASET OVERVIEW\")\nstats = [\n    [\"Training Samples\",  len(train_df)],\n    [\"Test Images\",       len(valid_id)],\n    [\"Total Test Rows\",   len(test_df)],\n    [\"Sampling Frequency\", f\"{test_df['fs'].iloc[0]} Hz\"],\n]\nprint(tabulate(stats, headers=[\"Metric\", \"Value\"], tablefmt=\"fancy_grid\"))\n\nprint(\"\\n📊 TEST DATA STRUCTURE\")\nprint(tabulate(test_df.head(24), headers='keys', tablefmt='fancy_grid', showindex=False))\n\n# ================================================================================\n# LOAD MODELS\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"LOADING MODELS\")\nprint(\"=\" * 80)\n\nfrom stage0_model  import Net as Stage0Net\nfrom stage0_common import (load_net, image_to_batch,\n                            output_to_predict as s0_output,\n                            normalise_by_homography)\n\nprint(\"\\n🔧 Loading Stage 0...\")\nstage0_net = Stage0Net(pretrained=False)\nstage0_net = load_net(stage0_net, f'{WEIGHT_DIR}/stage0-last.checkpoint.pth')\nstage0_net.to(DEVICE).eval()\nprint(\"✅ Stage 0 loaded\")\n\nfrom stage1_model  import Net as Stage1Net\nfrom stage1_common import (load_net as s1_load,\n                            output_to_predict as s1_output,\n                            rectify_image)\n\nprint(\"\\n🔧 Loading Stage 1...\")\nstage1_net = Stage1Net(pretrained=False)\nstage1_net = s1_load(stage1_net, f'{WEIGHT_DIR}/stage1-last.checkpoint.pth')\nstage1_net.to(DEVICE).eval()\nprint(\"✅ Stage 1 loaded\")\n\nfrom stage2_model  import Net as Stage2Net\nfrom stage2_common import (load_net as s2_load,\n                            pixel_to_series, filter_series_by_limits)\n\nprint(\"\\n🔧 Loading Stage 2...\")\nstage2_net = Stage2Net(pretrained=False)\nstage2_net = s2_load(stage2_net, f'{WEIGHT_DIR}/stage2-00005810.checkpoint.pth')\nstage2_net.to(DEVICE).eval()\nprint(\"✅ Stage 2 loaded\")\n\n# ================================================================================\n# STAGE 0 — NORMALISATION\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"PROCESSING TEST IMAGES\")\nprint(\"=\" * 80)\n\nFAIL_ID = []\n\nprint(\"\\n🔄 Stage 0: Normalisation...\")\nfor n, sample_id in enumerate(tqdm(valid_id, desc=\"Stage 0\")):\n    try:\n        image = cv2.imread(f'{KAGGLE_DIR}/test/{sample_id}.png',\n                           cv2.IMREAD_COLOR_RGB)\n        batch = image_to_batch(image)\n\n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output  = stage0_net(batch)\n                rotated, keypoint = s0_output(image, batch, output)\n                normalised, keypoint, homo = normalise_by_homography(rotated, keypoint)\n\n        cv2.imwrite(f'{OUT_DIR}/normalised/{sample_id}.norm.png',\n                    cv2.cvtColor(normalised, cv2.COLOR_RGB2BGR))\n        np.save(f'{OUT_DIR}/normalised/{sample_id}.homo.npy', homo)\n\n    except Exception as e:\n        print(f\"\\n⚠️  {sample_id}: {e}\")\n        FAIL_ID.append(sample_id)\n\n    torch.cuda.empty_cache()\n\nprint(f\"✅ Stage 0 complete: {len(valid_id) - len(FAIL_ID)}/{len(valid_id)} success\")\n\n# ================================================================================\n# STAGE 1 — RECTIFICATION\n# ================================================================================\n\nprint(\"\\n🔄 Stage 1: Rectification...\")\nfor n, sample_id in enumerate(tqdm(valid_id, desc=\"Stage 1\")):\n    if sample_id in FAIL_ID:\n        continue\n    try:\n        image = cv2.imread(f'{OUT_DIR}/normalised/{sample_id}.norm.png',\n                           cv2.IMREAD_COLOR_RGB)\n        batch = {\n            'image': torch.from_numpy(\n                np.ascontiguousarray(image.transpose(2, 0, 1))\n            ).unsqueeze(0)\n        }\n\n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output = stage1_net(batch)\n                gridpoint_xy, more = s1_output(image, batch, output)\n                rectified = rectify_image(image, gridpoint_xy)\n\n        cv2.imwrite(f'{OUT_DIR}/rectified/{sample_id}.rect.png',\n                    cv2.cvtColor(rectified, cv2.COLOR_RGB2BGR))\n        np.save(f'{OUT_DIR}/rectified/{sample_id}.gridpoint_xy.npy', gridpoint_xy)\n\n    except Exception as e:\n        print(f\"\\n⚠️  {sample_id}: {e}\")\n        FAIL_ID.append(sample_id)\n\n    torch.cuda.empty_cache()\n\nprint(f\"✅ Stage 1 complete: {len(valid_id) - len(FAIL_ID)}/{len(valid_id)} success\")\n\n# ================================================================================\n# STAGE 2 — SIGNAL EXTRACTION + IMPROVED POST-PROCESSING\n# ================================================================================\n\nprint(\"\\n🔄 Stage 2: Signal Extraction + Post-processing...\")\nfor n, sample_id in enumerate(tqdm(valid_id, desc=\"Stage 2\")):\n    if sample_id in FAIL_ID:\n        continue\n    try:\n        image = cv2.imread(f'{OUT_DIR}/rectified/{sample_id}.rect.png',\n                           cv2.IMREAD_COLOR_RGB)\n\n        # Fetch sampling frequency from metadata (may vary per record)\n        d_meta = test_df[(test_df['id'] == sample_id) & (test_df['lead'] == 'II')].iloc[0]\n        fs     = int(d_meta.fs)\n        length = int(d_meta.number_of_rows)\n\n        x0, x1 = 0, 2176\n        y0, y1 = 0, 1696\n        zero_mv = [703.5, 987.5, 1271.5, 1531.5]\n        t0, t1  = 118, 2080\n\n        crop  = image[y0:y1, x0:x1]\n        batch = {\n            'image': torch.from_numpy(\n                np.ascontiguousarray(crop.transpose(2, 0, 1))\n            ).unsqueeze(0)\n        }\n\n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output = stage2_net(batch)\n\n        pixel         = output['pixel'].float().cpu().numpy()[0]\n        series_in_px  = pixel_to_series(pixel[..., t0:t1], zero_mv, length)\n\n        # Convert pixel → mV\n        series = (np.array(zero_mv).reshape(4, 1) - series_in_px) / MV_TO_PIXEL\n\n        # ------------------------------------------------------------------\n        # IMPROVEMENT 1: Full signal-processing chain on each of the 4 rows\n        # (baseline removal → notch → adaptive smooth → clip)\n        # ------------------------------------------------------------------\n        for i in range(series.shape[0]):\n            series[i] = apply_signal_chain(series[i], fs=fs)\n\n        # Visualise first 3 samples for sanity check\n        if n < 3:\n            plt.figure(figsize=(20, 10))\n            plt.imshow(crop)\n            for lead_idx in range(4):\n                y_coords = np.array(zero_mv)[lead_idx] - (series[lead_idx] * MV_TO_PIXEL)\n                x_coords = np.arange(t0, t1)\n                lim = min(len(x_coords), len(y_coords))\n                plt.plot(x_coords[:lim], y_coords[:lim], color='red', linewidth=1)\n            plt.title(f\"Prediction Overlay: {sample_id}\")\n            plt.show()\n\n        np.save(f'{OUT_DIR}/digitalised/{sample_id}.series.npy', series)\n\n    except Exception as e:\n        print(f\"\\n⚠️  {sample_id}: {e}\")\n        FAIL_ID.append(sample_id)\n\n    torch.cuda.empty_cache()\n\nprint(f\"✅ Stage 2 complete: {len(valid_id) - len(FAIL_ID)}/{len(valid_id)} success\")\n\n# ================================================================================\n# BUILD SUBMISSION WITH SMART SPLITTING + EINTHOVEN CORRECTION\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"CREATING SUBMISSION\")\nprint(\"=\" * 80)\n\n# Lead name mapping per row\nROW_LEADS = [\n    ['I',   'aVR', 'V1', 'V4'],   # row 0\n    ['II',  'aVL', 'V2', 'V5'],   # row 1\n    ['III', 'aVF', 'V3', 'V6'],   # row 2\n]\n\nsubmission_data = []\ngb = test_df.groupby('id')\n\nfor rec_idx, (sample_id, df) in enumerate(tqdm(gb, desc=\"Building submission\")):\n    try:\n        series = np.load(f'{OUT_DIR}/digitalised/{sample_id}.series.npy')\n        _, L   = series.shape\n\n        # Fetch fs for this record\n        fs = int(df['fs'].iloc[0])\n\n        series_by_lead = {}\n\n        # ------------------------------------------------------------------\n        # IMPROVEMENT 2: Smart boundary detection for rows 0–2\n        # Use find_split_boundaries instead of hard equal quarters\n        # ------------------------------------------------------------------\n        for row_idx in range(3):\n            lead_names = ROW_LEADS[row_idx]\n            row_signal = series[row_idx]\n\n            boundaries = find_split_boundaries(row_signal, n_splits=4)\n            split_pts  = [0] + boundaries + [L]\n            segments   = [row_signal[split_pts[i]:split_pts[i+1]] for i in range(4)]\n\n            for lead_name, seg in zip(lead_names, segments):\n                series_by_lead[lead_name] = seg\n\n        # Row 3 is the full-length rhythm lead II — keep as-is\n        series_by_lead['II'] = series[3]\n\n        # ------------------------------------------------------------------\n        # IMPROVEMENT 3: Enforce Einthoven's law on I, II, III\n        # Corrects systematic digitization offsets using physical constraints\n        # ------------------------------------------------------------------\n        # All three must have the same length for the correction;\n        # interpolate to lead II length if needed\n        def _align(s, target_len):\n            if len(s) == target_len:\n                return s\n            x_old = np.linspace(0, 1, len(s),         endpoint=False)\n            x_new = np.linspace(0, 1, target_len,      endpoint=False)\n            return np.interp(x_new, x_old, s)\n\n        ref_len = len(series_by_lead['II'])\n        lI   = _align(series_by_lead['I'],   ref_len)\n        lII  = series_by_lead['II'].copy()\n        lIII = _align(series_by_lead['III'], ref_len)\n\n        lI, lII, lIII = enforce_einthoven(lI, lII, lIII)\n\n        series_by_lead['I']   = lI\n        series_by_lead['II']  = lII\n        series_by_lead['III'] = lIII\n\n    except Exception:\n        # Fallback: zero-fill all leads\n        series_by_lead = {\n            d.lead: np.zeros(int(d.number_of_rows), dtype=np.float32)\n            for _, d in df.iterrows()\n        }\n\n    # ------------------------------------------------------------------\n    # Resample to target length & write rows\n    # ------------------------------------------------------------------\n    for _, d in df.iterrows():\n        s          = series_by_lead.get(d.lead, np.zeros(int(d.number_of_rows)))\n        target_len = int(d.number_of_rows)\n\n        if len(s) != target_len:\n            x_old = np.linspace(0, 1, len(s),         endpoint=False)\n            x_new = np.linspace(0, 1, target_len,      endpoint=False)\n            s = np.interp(x_new, x_old, s)\n\n        s = s.astype(np.float32)\n\n        for t in range(target_len):\n            submission_data.append({\n                'id':    f'{sample_id}_{t}_{d.lead}',\n                'value': float(s[t])\n            })\n\nsubmission = pd.DataFrame(submission_data)\nsubmission.to_csv('submission.csv', index=False)\nprint(f\"\\n✅ Submission created: {len(submission):,} rows\")\n\n# ================================================================================\n# ANALYSIS\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"SUBMISSION ANALYSIS\")\nprint(\"=\" * 80)\n\nsubmission['lead'] = submission['id'].str.split('_').str[2]\n\nprint(\"\\n📊 PER-LEAD STATISTICS\")\nlead_stats = []\nfor lead in ['I', 'II', 'III', 'aVR', 'aVL', 'aVF', 'V1', 'V2', 'V3', 'V4', 'V5', 'V6']:\n    data = submission[submission['lead'] == lead]['value']\n    lead_stats.append([\n        lead, len(data),\n        f\"{data.mean():.6f}\", f\"{data.std():.6f}\",\n        f\"{data.min():.6f}\", f\"{data.max():.6f}\",\n        (data != 0).sum()\n    ])\n\nprint(tabulate(\n    lead_stats,\n    headers=['Lead', 'Count', 'Mean', 'Std', 'Min', 'Max', 'Non-Zero'],\n    tablefmt='fancy_grid'\n))\n\nprint(f\"\\n✅ Failed IDs: {len(FAIL_ID)}\")\nif FAIL_ID:\n    print(f\"   {FAIL_ID}\")\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"✅ COMPLETE\")\nprint(\"=\" * 80)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:46:59.546573Z","iopub.execute_input":"2026-04-26T14:46:59.546933Z","iopub.status.idle":"2026-04-26T14:47:14.286929Z","shell.execute_reply.started":"2026-04-26T14:46:59.546907Z","shell.execute_reply":"2026-04-26T14:47:14.286012Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## experiment 4 ","metadata":{}},{"cell_type":"code","source":"\"\"\"\nProperly using hengck23's pretrained models with fixed imports\n\"\"\"\n\nprint(\"=\" * 80)\nprint(\"ECG INFERENCE\")\nprint(\"=\" * 80)\n\n# Setup paths\nimport sys\nbase_path = '/kaggle/input/hengck23-submit-physionet/hengck23-submit-physionet'\nsys.path.insert(0, base_path)\n\nKAGGLE_DIR = '/kaggle/input/physionet-ecg-image-digitization'\nWEIGHT_DIR = f'{base_path}/weight'\nOUT_DIR = '/kaggle/working/output-submit'\nDEVICE = 'cuda' if torch.cuda.is_available() else 'cpu'\nFLOAT_TYPE = torch.float16\n\nos.makedirs(f'{OUT_DIR}/normalised', exist_ok=True)\nos.makedirs(f'{OUT_DIR}/rectified', exist_ok=True)\nos.makedirs(f'{OUT_DIR}/digitalised', exist_ok=True)\n\nprint(f\"\\n🔧 Device: {DEVICE}\")\nprint(f\"📁 Weights: {WEIGHT_DIR}\")\n\n# ================================================================================\n# DATA INSPECTION\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"DATA INSPECTION\")\nprint(\"=\" * 80)\n\ntrain_df = pd.read_csv(f'{KAGGLE_DIR}/train.csv')\ntest_df = pd.read_csv(f'{KAGGLE_DIR}/test.csv')\ntest_df['id'] = test_df['id'].astype(str)\nvalid_id = test_df['id'].unique().tolist()\n\nprint(\"\\n📊 DATASET OVERVIEW\")\nstats = [\n    [\"Training Samples\", len(train_df)],\n    [\"Test Images\", len(valid_id)],\n    [\"Total Test Rows\", len(test_df)],\n    [\"Sampling Frequency\", f\"{test_df['fs'].iloc[0]} Hz\"],\n]\nprint(tabulate(stats, headers=[\"Metric\", \"Value\"], tablefmt=\"fancy_grid\"))\n\nprint(\"\\n📊 TEST DATA STRUCTURE\")\nprint(tabulate(test_df.head(24), headers='keys', tablefmt='fancy_grid', showindex=False))\n\n# ================================================================================\n# LOAD MODELS\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"LOADING MODELS\")\nprint(\"=\" * 80)\n\nfrom stage0_model import Net as Stage0Net\nfrom stage0_common import load_net, image_to_batch, output_to_predict as s0_output, normalise_by_homography\n\nprint(\"\\n🔧 Loading Stage 0...\")\nstage0_net = Stage0Net(pretrained=False)\nstage0_net = load_net(stage0_net, f'{WEIGHT_DIR}/stage0-last.checkpoint.pth')\nstage0_net.to(DEVICE).eval()\nprint(\"✅ Stage 0 loaded\")\n\nfrom stage1_model import Net as Stage1Net\nfrom stage1_common import load_net as s1_load, output_to_predict as s1_output, rectify_image\n\nprint(\"\\n🔧 Loading Stage 1...\")\nstage1_net = Stage1Net(pretrained=False)\nstage1_net = s1_load(stage1_net, f'{WEIGHT_DIR}/stage1-last.checkpoint.pth')\nstage1_net.to(DEVICE).eval()\nprint(\"✅ Stage 1 loaded\")\n\nfrom stage2_model import Net as Stage2Net\nfrom stage2_common import load_net as s2_load, pixel_to_series, filter_series_by_limits\n\nprint(\"\\n🔧 Loading Stage 2...\")\nstage2_net = Stage2Net(pretrained=False)\nstage2_net = s2_load(stage2_net, f'{WEIGHT_DIR}/stage2-00005810.checkpoint.pth')\nstage2_net.to(DEVICE).eval()\nprint(\"✅ Stage 2 loaded\")\n\n# ================================================================================\n# PROCESSING PIPELINE\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"PROCESSING TEST IMAGES\")\nprint(\"=\" * 80)\n\nFAIL_ID = []\n\n# STAGE 0\nprint(\"\\n🔄 Stage 0: Normalization...\")\nfor n, sample_id in enumerate(tqdm(valid_id, desc=\"Stage 0\")):\n    try:\n        image = cv2.imread(f'{KAGGLE_DIR}/test/{sample_id}.png', cv2.IMREAD_COLOR_RGB)\n        batch = image_to_batch(image)\n        \n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output = stage0_net(batch)\n                rotated, keypoint = s0_output(image, batch, output)\n                normalised, keypoint, homo = normalise_by_homography(rotated, keypoint)\n        \n        cv2.imwrite(f'{OUT_DIR}/normalised/{sample_id}.norm.png', cv2.cvtColor(normalised, cv2.COLOR_RGB2BGR))\n        np.save(f'{OUT_DIR}/normalised/{sample_id}.homo.npy', homo)\n        \n    except Exception as e:\n        print(f\"\\n⚠️  {sample_id}: {e}\")\n        FAIL_ID.append(sample_id)\n    \n    torch.cuda.empty_cache()\n\nprint(f\"✅ Stage 0 complete: {len(valid_id) - len(FAIL_ID)}/{len(valid_id)} success\")\n\n# STAGE 1\nprint(\"\\n🔄 Stage 1: Rectification...\")\nfor n, sample_id in enumerate(tqdm(valid_id, desc=\"Stage 1\")):\n    if sample_id in FAIL_ID:\n        continue\n    \n    try:\n        image = cv2.imread(f'{OUT_DIR}/normalised/{sample_id}.norm.png', cv2.IMREAD_COLOR_RGB)\n        batch = {'image': torch.from_numpy(np.ascontiguousarray(image.transpose(2, 0, 1))).unsqueeze(0)}\n        \n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output = stage1_net(batch)\n                gridpoint_xy, more = s1_output(image, batch, output)\n                rectified = rectify_image(image, gridpoint_xy)\n        \n        cv2.imwrite(f'{OUT_DIR}/rectified/{sample_id}.rect.png', cv2.cvtColor(rectified, cv2.COLOR_RGB2BGR))\n        np.save(f'{OUT_DIR}/rectified/{sample_id}.gridpoint_xy.npy', gridpoint_xy)\n        \n    except Exception as e:\n        print(f\"\\n⚠️  {sample_id}: {e}\")\n        FAIL_ID.append(sample_id)\n    \n    torch.cuda.empty_cache()\n\nprint(f\"✅ Stage 1 complete: {len(valid_id) - len(FAIL_ID)}/{len(valid_id)} success\")\n\n# STAGE 2\nprint(\"\\n🔄 Stage 2: Signal Extraction...\")\nfor n, sample_id in enumerate(tqdm(valid_id, desc=\"Stage 2\")):\n    if sample_id in FAIL_ID:\n        continue\n    \n    try:\n        image = cv2.imread(f'{OUT_DIR}/rectified/{sample_id}.rect.png', cv2.IMREAD_COLOR_RGB)\n        d = test_df[(test_df['id']==sample_id) & (test_df['lead']=='II')].iloc[0]\n        length = d.number_of_rows\n        \n        x0, x1 = 0, 2176\n        y0, y1 = 0, 1696\n        zero_mv = [703.5, 987.5, 1271.5, 1531.5]\n        mv_to_pixel = 80.0\n        t0, t1 = 118, 2080\n        \n        crop = image[y0:y1, x0:x1]\n        batch = {'image': torch.from_numpy(np.ascontiguousarray(crop.transpose(2, 0, 1))).unsqueeze(0)}\n        \n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output = stage2_net(batch)\n        \n        pixel = output['pixel'].float().data.cpu().numpy()[0]\n        series_in_pixel = pixel_to_series(pixel[..., t0:t1], zero_mv, length)\n        series = (np.array(zero_mv).reshape(4, 1) - series_in_pixel) / mv_to_pixel\n        series = filter_series_by_limits(series)\n\n        for i in range(series.shape[0]):\n            # 2. Add a Median Filter BEFORE Savgol. \n            # This instantly kills single-pixel anomalies/spikes without blurring the sharp QRS peaks.\n            # from scipy.signal import medfilt\n            # series[i] = medfilt(series[i], kernel_size=5)\n            \n            # Then apply Savgol to smooth the remaining curve\n            series[i] = savgol_filter(series[i], window_length=11, polyorder=3)\n\n        if n < 3: \n            plt.figure(figsize=(20, 10))\n            # Display the crop used for inference\n            plt.imshow(crop) \n            \n            # We need to reverse the math: series -> pixels\n            for lead_idx in range(4): # The model outputs 4 rows of leads\n                # Reverse the normalization: (zero_mv - (series * 79.0))\n                y_coords = zero_mv[lead_idx] - (series[lead_idx] * mv_to_pixel)\n                \n                # Plot only the valid range\n                x_coords = np.arange(t0, t1)\n                # Ensure lengths match\n                limit = min(len(x_coords), len(y_coords))\n                plt.plot(x_coords[:limit], y_coords[:limit], color='red', linewidth=1)\n                \n            plt.title(f\"Prediction Overlay: {sample_id}\")\n            plt.show()\n\n        \n        np.save(f'{OUT_DIR}/digitalised/{sample_id}.series.npy', series)\n        \n    except Exception as e:\n        print(f\"\\n⚠️  {sample_id}: {e}\")\n        FAIL_ID.append(sample_id)\n    \n    torch.cuda.empty_cache()\n\nprint(f\"✅ Stage 2 complete: {len(valid_id) - len(FAIL_ID)}/{len(valid_id)} success\")\n\n# ================================================================================\n# CREATE SUBMISSION\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"CREATING SUBMISSION\")\nprint(\"=\" * 80)\n\nsubmission_data = []\ngb = test_df.groupby('id')\n\nfor rec_idx, (sample_id, df) in enumerate(tqdm(gb, desc=\"Building submission\")):\n    \n    try:\n        series = np.load(f'{OUT_DIR}/digitalised/{sample_id}.series.npy')\n        _4_, L = series.shape\n        \n        series_by_lead = {}\n        for l in range(3):\n            lead_names = [\n                ['I', 'aVR', 'V1', 'V4'],\n                ['II', 'aVL', 'V2', 'V5'],\n                ['III', 'aVF', 'V3', 'V6'],\n            ][l]\n            \n            index = [int(round(1*L/4)), int(round(2*L/4)), int(round(3*L/4))]\n            split = np.split(series[l], index)\n            for k, s in zip(lead_names, split):\n                series_by_lead[k] = s\n        \n        series_by_lead['II'] = series[3]\n        \n    except:\n        series_by_lead = {}\n        for _, d in df.iterrows():\n            series_by_lead[d.lead] = np.zeros(d.number_of_rows, dtype=np.float32)\n    \n    for _, d in df.iterrows():\n        s = series_by_lead[d.lead]\n        target_len = int(d.number_of_rows)\n        \n        if len(s) != target_len:\n            x_old = np.linspace(0.0, 1.0, len(s), endpoint=False)\n            x_new = np.linspace(0.0, 1.0, target_len, endpoint=False)\n            s = np.interp(x_new, x_old, s)\n        \n        s = s.astype(np.float32)\n\n        # s = s - np.nanmedian(s)\n        \n        # # 2. Clip impossible voltage spikes caused by ink smudges/text\n        s = np.clip(s, -3.0, 3.0)\n        \n        for t in range(target_len):\n            submission_data.append({\n                'id': f'{sample_id}_{t}_{d.lead}',\n                'value': float(s[t])\n            })\n\nsubmission = pd.DataFrame(submission_data)\nsubmission.to_csv('submission.csv', index=False)\n\nprint(f\"\\n✅ Submission created: {len(submission):,} rows\")\n\n# ================================================================================\n# ANALYSIS\n# ================================================================================\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"SUBMISSION ANALYSIS\")\nprint(\"=\" * 80)\n\nsubmission['lead'] = submission['id'].str.split('_').str[2]\n\nprint(\"\\n📊 PER-LEAD STATISTICS\")\nlead_stats = []\nfor lead in ['I', 'II', 'III', 'aVR', 'aVL', 'aVF', 'V1', 'V2', 'V3', 'V4', 'V5', 'V6']:\n    data = submission[submission['lead'] == lead]['value']\n    lead_stats.append([\n        lead, len(data), f\"{data.mean():.6f}\", f\"{data.std():.6f}\",\n        f\"{data.min():.6f}\", f\"{data.max():.6f}\", (data != 0).sum()\n    ])\n\nprint(tabulate(lead_stats, headers=['Lead', 'Count', 'Mean', 'Std', 'Min', 'Max', 'Non-Zero'],\n              tablefmt='fancy_grid'))\n\nprint(f\"\\n✅ Failed IDs: {len(FAIL_ID)}\")\nif FAIL_ID:\n    print(f\"   {FAIL_ID}\")\n\nprint(\"\\n\" + \"=\" * 80)\nprint(\"✅ COMPLETE\")\nprint(\"=\" * 80)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-26T14:47:14.287868Z","iopub.execute_input":"2026-04-26T14:47:14.288173Z","iopub.status.idle":"2026-04-26T14:47:26.322709Z","shell.execute_reply.started":"2026-04-26T14:47:14.288148Z","shell.execute_reply":"2026-04-26T14:47:26.321760Z"}},"outputs":[],"execution_count":null}]}