{"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":16005256,"datasetId":10265254,"databundleVersionId":16968762},{"sourceType":"datasetVersion","sourceId":13731160,"datasetId":8733970,"databundleVersionId":14479231},{"sourceType":"datasetVersion","sourceId":16010811,"datasetId":9997958,"databundleVersionId":16974783},{"sourceType":"datasetVersion","sourceId":14695416,"datasetId":9387663,"databundleVersionId":15539655},{"sourceType":"kernelVersion","sourceId":295259486}],"dockerImageVersionId":31192,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"!ls /kaggle/input/notebooks/packagemanager/pm-108440548-at-02-01-2026-09-40-45\n!pip install --no-index --no-deps \\\n  /kaggle/input/notebooks/packagemanager/pm-108440548-at-02-01-2026-09-40-45/segmentation_models_pytorch-0.5.0-py3-none-any.whl","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:22.661047Z","iopub.execute_input":"2026-05-02T16:22:22.661818Z","iopub.status.idle":"2026-05-02T16:22:25.662829Z","shell.execute_reply.started":"2026-05-02T16:22:22.661789Z","shell.execute_reply":"2026-05-02T16:22:25.662115Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import segmentation_models_pytorch as smp\nprint(smp.__version__)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:25.664362Z","iopub.execute_input":"2026-05-02T16:22:25.664622Z","iopub.status.idle":"2026-05-02T16:22:36.755497Z","shell.execute_reply.started":"2026-05-02T16:22:25.664589Z","shell.execute_reply":"2026-05-02T16:22:36.754873Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"try:\n    import cc3d  # connected components in 3D (used in Stage 1 grid detection)\nexcept:         \n    #https://pypi.org/project/connected-components-3d/\n    #!pip install connected-components-3d\n\n    !ls /kaggle/input/hengck23-demo-submit-physionet/setup\n    !pip install connected-components-3d --no-index --find-links=file:///kaggle/input/hengck23-demo-submit-physionet/setup/\n\nimport cc3d\nimport cv2          # OpenCV: all image processing operations\nimport pandas as pd\nimport numpy as np\nimport torch\nimport matplotlib.pyplot as plt\nimport matplotlib\n#matplotlib.use('TkAgg')\nimport shutil\n\nimport sys\nsys.path.insert(0, '/kaggle/input/datasets/zahouaniyacine/my-stage2-lead-model') #my stage2 attention model added \nsys.path.append('/kaggle/input/hengck23-demo-submit-physionet')\nsys.path.append('/kaggle/input/datasets/takashisomeya/physionet-final-submission-models')\n# This adds hengck23's pre-packaged code to Python's module search path\n# Inside this package:\n# - stage0model.py     : the Stage 0 neural network definition\n# - stage0common.py    : Stage 0 helper functions (image→batch, output→prediction)\n# - stage1model.py     : the Stage 1 neural network definition\n# - stage1common.py    : Stage 1 helper functions (output→grid points, rectification)\n# - stage2model.py     : Stage 2 model + pixel_to_series conversion\n# - stage2common.py    : Stage 2 helper functions\n\nprint('import ok!!!')","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:36.756218Z","iopub.execute_input":"2026-05-02T16:22:36.756397Z","iopub.status.idle":"2026-05-02T16:22:40.745282Z","shell.execute_reply.started":"2026-05-02T16:22:36.756382Z","shell.execute_reply":"2026-05-02T16:22:40.744606Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import stage2_lead_model\nprint(stage2_lead_model.__file__)\n# Should print: /kaggle/input/my-physionet-stage2/stage2_lead_model.py\n# NOT: /kaggle/input/physionet-final-submission-models/stage2_lead_model.py","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:40.746311Z","iopub.execute_input":"2026-05-02T16:22:40.746949Z","iopub.status.idle":"2026-05-02T16:22:40.816523Z","shell.execute_reply.started":"2026-05-02T16:22:40.746914Z","shell.execute_reply":"2026-05-02T16:22:40.816000Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"import os\nimport numpy as np\nimport pandas as pd\nimport scipy.signal\nimport scipy.optimize\n\nLEADS = ['I', 'II', 'III', 'aVR', 'aVL', 'aVF', 'V1', 'V2', 'V3', 'V4', 'V5', 'V6']\nMAX_TIME_SHIFT = 0.2\nPERFECT_SCORE = 384.0\n\nclass ParticipantVisibleError(Exception):\n    pass\n\ndef compute_power(label: np.ndarray, prediction: np.ndarray):\n    if label.ndim != 1 or prediction.ndim != 1:\n        raise ParticipantVisibleError(\"Inputs must be 1-dimensional arrays.\")\n    finite_mask = np.isfinite(prediction)\n    if not np.any(finite_mask):\n        raise ParticipantVisibleError(\"Prediction contains no finite values.\")\n    prediction = prediction.copy()\n    prediction[~finite_mask] = 0\n    noise = label - prediction\n    p_signal = np.sum(label ** 2)\n    p_noise = np.sum(noise ** 2)\n    return p_signal, p_noise\n\ndef compute_snr(signal_power: float, noise_power: float):\n    if noise_power == 0:\n        return PERFECT_SCORE\n    elif signal_power == 0:\n        return 0.0\n    else:\n        return min(signal_power / noise_power, PERFECT_SCORE)\n\ndef align_signals(label: np.ndarray, pred: np.ndarray, max_shift=float(\"inf\")):\n    if np.any(~np.isfinite(label)):\n        raise ParticipantVisibleError(\"Label contains non-finite values.\")\n    if np.sum(np.isfinite(pred)) == 0:\n        raise ParticipantVisibleError(\"Prediction cannot be all NaN/Inf.\")\n\n    label_arr = np.asarray(label, dtype=np.float64)\n    pred_arr  = np.asarray(pred,  dtype=np.float64)\n\n    label_centered = label_arr - np.mean(label_arr)\n    pred_centered  = pred_arr  - np.mean(pred_arr)\n\n    correlation = scipy.signal.correlate(label_centered, pred_centered, mode=\"full\")\n    n_label = label_arr.size\n    n_pred  = pred_arr.size\n    lags = scipy.signal.correlation_lags(n_label, n_pred, mode=\"full\")\n\n    valid_lags_mask = (lags >= -max_shift) & (lags <= max_shift)\n    valid_corr = correlation[valid_lags_mask]\n    valid_lags = lags[valid_lags_mask]\n\n    max_corr = np.nanmax(valid_corr)\n    best_candidates = np.flatnonzero(valid_corr == max_corr)\n    best_idx_local = min(best_candidates, key=lambda i: abs(valid_lags[i]))\n    time_shift = int(valid_lags[best_idx_local])\n\n    start_padding_len = max(time_shift, 0)\n    pred_slice_start  = max(-time_shift, 0)\n    pred_slice_end    = min(n_label - time_shift, n_pred)\n    end_padding_len   = max(n_label - n_pred - time_shift, 0)\n\n    aligned_pred = np.concatenate([\n        np.full(start_padding_len, np.nan),\n        pred_arr[pred_slice_start:pred_slice_end],\n        np.full(end_padding_len, np.nan)\n    ])\n\n    def objective_func(vshift):\n        return np.nansum((label_arr - aligned_pred - vshift) ** 2)\n\n    if np.any(np.isfinite(label_arr) & np.isfinite(aligned_pred)):\n        result = scipy.optimize.minimize_scalar(objective_func, method=\"Brent\")\n        vertical_shift = result.x\n        aligned_pred = aligned_pred + vertical_shift\n\n    return aligned_pred","metadata":{"execution":{"iopub.status.busy":"2026-05-02T16:22:40.818254Z","iopub.execute_input":"2026-05-02T16:22:40.818779Z","iopub.status.idle":"2026-05-02T16:22:41.358416Z","shell.execute_reply.started":"2026-05-02T16:22:40.818761Z","shell.execute_reply":"2026-05-02T16:22:41.357862Z"},"_kg_hide-input":true}},{"cell_type":"code","source":"import importlib \nimport torch.nn as nn \nimport stage2_lead_model as _m\n\nimportlib.reload(_m)\n\nclass CrossLeadAttentionFusion(nn.Module):\n    def __init__(self, channels, num_leads=4, num_heads=4, dropout=0.1):\n        super().__init__()\n        self.channels = channels\n        self.num_leads = num_leads\n        self.num_heads = next((h for h in [8, 4, 2, 1] if channels % h == 0), 1)\n\n        self.norm1 = nn.LayerNorm(channels)\n        self.norm2 = nn.LayerNorm(channels)\n\n        self.attn = nn.MultiheadAttention(\n            embed_dim=channels,\n            num_heads=self.num_heads,\n            dropout=dropout,\n            batch_first=True,\n        )\n\n        self.ffn = nn.Sequential(\n            nn.Linear(channels, channels * 2),\n            nn.GELU(),\n            nn.Dropout(dropout),\n            nn.Linear(channels * 2, channels),\n        )\n\n        self.attn_scale = nn.Parameter(torch.tensor(0.1))\n        self.ffn_scale  = nn.Parameter(torch.tensor(0.1))\n\n    def forward(self, x, batch_size=None, batchsize=None):\n        if batch_size is None:\n            batch_size = batchsize\n        if batch_size is None:\n            raise ValueError(\"batch_size must be provided\")\n\n        B = batch_size\n        _, C, H, W = x.shape\n\n        x_leads = x.view(B, self.num_leads, C, H, W)\n        x_hw = x_leads.permute(0, 3, 4, 1, 2)\n        x_seq = x_hw.reshape(B * H * W, self.num_leads, C)\n\n        attn_in = self.norm1(x_seq)\n        attn_out, _ = self.attn(attn_in, attn_in, attn_in)\n        x_seq = x_seq + self.attn_scale * attn_out\n\n        ffn_in = self.norm2(x_seq)\n        x_seq = x_seq + self.ffn_scale * self.ffn(ffn_in)\n\n        x_hw = x_seq.reshape(B, H, W, self.num_leads, C)\n        x_out = x_hw.permute(0, 3, 4, 1, 2)\n        return x_out.reshape(B * self.num_leads, C, H, W)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:41.359111Z","iopub.execute_input":"2026-05-02T16:22:41.359432Z","iopub.status.idle":"2026-05-02T16:22:41.372659Z","shell.execute_reply.started":"2026-05-02T16:22:41.359414Z","shell.execute_reply":"2026-05-02T16:22:41.372030Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"MODE   = 'submit'  # submit  local fake\nDEVICE = 'cuda'\nFLOAT_TYPE = torch.float16 #torch.bfloat16\nFAIL_ID = []\n\nKAGGLE_DIR = \\\n\t'/kaggle/input/physionet-ecg-image-digitization'\nWEIGHT_DIR = \\\n\t'/kaggle/input/hengck23-demo-submit-physionet/weight'\nOUT_DIR = \\\n    f'/kaggle/working/output-{MODE}'\n\ndef make_test_fake_df(): \n    valid_df = pd.read_csv(f'{KAGGLE_DIR}/train.csv')\n    valid_df.loc[:,'id']=valid_df['id'].astype(str) \n    fake_test_df=[]\n    for i,d in valid_df.iterrows():\n        #if i==4: break\n        image_id = d['id']\n    \n        truth_df = pd.read_csv(f'{KAGGLE_DIR}/train/{image_id}/{image_id}.csv')\n        non_nan_count = truth_df.count()\n        #print(i,image_id,non_nan_count)\n        #print(non_nan_count.index)\n    \n        #lead\tfs\tnumber_of_rows \n        this_df = pd.DataFrame({\n            'id':image_id ,\n            'lead':non_nan_count.index,\n            'fs': d['fs'],\n            'number_of_rows':non_nan_count.values \n        })\n        fake_test_df.append(this_df)\n        if i==0: print(this_df)\n    fake_test_df = pd.concat(fake_test_df)\n    return fake_test_df","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:41.373265Z","iopub.execute_input":"2026-05-02T16:22:41.373503Z","iopub.status.idle":"2026-05-02T16:22:41.380246Z","shell.execute_reply.started":"2026-05-02T16:22:41.373483Z","shell.execute_reply":"2026-05-02T16:22:41.379607Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# set valid/test data\nif MODE == 'local':\n    valid_df = pd.read_csv(f'{KAGGLE_DIR}/train.csv')\n    valid_df['id'] = valid_df['id'].astype(str)\n\n    all_ids = valid_df['id'].unique().tolist()\n    rng = np.random.default_rng(42)\n    rng.shuffle(all_ids)\n\n    N_VAL_IDS = 50\n    val_image_ids = all_ids[:N_VAL_IDS]\n\n    TYPE_IDS = ['0001', '0003', '0004', '0005', '0006', '0009', '0010', '0011', '0012']\n    valid_id = [f'{image_id}-{type_id}' for image_id in val_image_ids for type_id in TYPE_IDS]\n\n","metadata":{"execution":{"iopub.status.busy":"2026-05-02T16:22:41.380895Z","iopub.execute_input":"2026-05-02T16:22:41.381169Z","iopub.status.idle":"2026-05-02T16:22:41.417779Z","shell.execute_reply.started":"2026-05-02T16:22:41.381152Z","shell.execute_reply":"2026-05-02T16:22:41.417253Z"},"_kg_hide-input":true,"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"    \nif MODE == 'submit':\n\tvalid_df = pd.read_csv(f'{KAGGLE_DIR}/test.csv')\n\tvalid_df['id']=valid_df['id'].astype(str) \n\tvalid_id = valid_df['id'].unique().tolist()\n\nif MODE == 'fake':\n\tvalid_df = make_test_fake_df()\n\tvalid_df['id']=valid_df['id'].astype(str) \n\tvalid_id = valid_df['id'].unique().tolist()\n\n#--------------------------------------\n\ndef read_image(sample_id):\n    if MODE == 'local':\n        image_id, type_id = sample_id.split('-')\n        image = cv2.imread(f'{KAGGLE_DIR}/train/{image_id}/{image_id}-{type_id}.png', cv2.IMREAD_COLOR_RGB)\n        return image\n    if MODE == 'submit':\n        image_id = sample_id\n        image = cv2.imread(f'{KAGGLE_DIR}/test/{image_id}.png', cv2.IMREAD_COLOR_RGB)\n        return image\n    if MODE == 'fake':\n        image_id = sample_id \n        type_id = ['0001', '0003', '0004', '0005', '0006', '0009', '0010', '0011', '0012'][\n            int(image_id)%9\n        ] \n        image = cv2.imread(f'{KAGGLE_DIR}/train/{image_id}/{image_id}-{type_id}.png', cv2.IMREAD_COLOR_RGB)\n        return image\n\ndef read_sampling_length(sample_id):\n\tif MODE == 'local':\n\t\timage_id, type_id = sample_id.split('-')\n\t\td = valid_df[valid_df['id']==image_id].iloc[0]\n\t\tlength = d.sig_len\n\t\treturn length\n\tif MODE == 'submit':\n\t\timage_id = sample_id\n\t\td = valid_df[\n\t\t\t(valid_df['id']==image_id) & (valid_df['lead']=='II')\n\t\t].iloc[0]\n\t\tlength = d.number_of_rows\n\t\treturn length\n\tif MODE == 'fake':\n\t\timage_id = sample_id\n\t\td = valid_df[\n\t\t\t(valid_df['id']==image_id) & (valid_df['lead']=='II')\n\t\t].iloc[0]\n\t\tlength = d.number_of_rows\n\t\treturn length\n\n#valid_id = valid_id[:300]\nprint('valid_id:', len(valid_id))\nprint('\\t', valid_id[:3], '...')\nprint('setting ok!!!\\n')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:41.418401Z","iopub.execute_input":"2026-05-02T16:22:41.418673Z","iopub.status.idle":"2026-05-02T16:22:41.427108Z","shell.execute_reply.started":"2026-05-02T16:22:41.418651Z","shell.execute_reply":"2026-05-02T16:22:41.426434Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"print(\"MODE:\", MODE)\nprint(\"valid_df shape:\", valid_df.shape)\nprint(\"unique ids:\", valid_df['id'].nunique())\nprint(valid_df.head(20))\nprint(valid_df['id'].value_counts())","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:41.427862Z","iopub.execute_input":"2026-05-02T16:22:41.428058Z","iopub.status.idle":"2026-05-02T16:22:41.450355Z","shell.execute_reply.started":"2026-05-02T16:22:41.428043Z","shell.execute_reply":"2026-05-02T16:22:41.449632Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# stage0\nprint('*** STARTING STAGE0 ***')\n\nfrom stage0_model import Net as Stage0Net\nfrom stage0_common import *\n#image_to_batch,output_to_predict,normalise_by_homography,draw_results_stage0\n\nos.makedirs(f'{OUT_DIR}/normalised', exist_ok=True)\n\ndef run_stage0():\n   # Load the pre-trained Stage 0 model\n\tstage0_net = Stage0Net(pretrained=False)\n\tstage0_net = load_net(stage0_net, f'{WEIGHT_DIR}/stage0-last.checkpoint.pth')\n\tstage0_net.to(DEVICE)\n    # .to(DEVICE) moves all model parameters to GPU\n    # Essential: if model is on CPU but data is on GPU → error\n    # If model is on GPU but data is on CPU → error\n    # They must be on the same device\n\n\tstart_timer = timer()\n\tfor n, sample_id in enumerate(valid_id):\n\t\ttimestamp = time_to_str(timer() - start_timer, 'sec')\n\t\tprint(f'\\r\\t {n:4d} {sample_id}', timestamp, end='', flush=True)\n\n\t\timage = read_image(sample_id)\n        # image shape: (H, W, 3) — e.g., (3024, 4032, 3) for a phone photo\n        # dtype: uint8, values 0-255\n\t\tbatch = image_to_batch(image)\n        # Converts raw image to model-ready tensor:\n        # 1. Resize to fixed input size (e.g., 512×512 or 768×1024)\n        # 2. Normalize: (pixel/255 - mean) / std  (ImageNet mean and std stats cause we will use res/efficentNet that is trained on thatits what it excpects)\n        #eg:two channels are both 0.05 above their mean, that 0.05 should count as more important in a channel that usually varies very little than in a channel that naturally varies a lot. Dividing by std captures exactly that idea.\n        # 3. Rearrange: HWC → CHW (channels first, required by PyTorch)\n        # 4. Add batch dim: CHW → BCHW (B=1 for single image)\n        # 5. Move to DEVICE\n        # Result shape: (1, 3, H_model, W_model)\n\n\t\twith torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n        # autocast: automatically uses float16 for most operations\n        # → 2× less GPU memory, 1.5-2× faster\n        # → slight precision loss, negligible for inference\n\t\t\twith torch.no_grad():\n\t\t\t\toutput = stage0_net(batch)\n        # no_grad(): tells PyTorch not to build computation graph\n        # → saves memory (no gradient tensors allocated)\n        # → faster (skip gradient tracking)\n        # ALWAYS use no_grad() during inference\n\n\t\t\t\ttry:\n\t\t\t\t\trotated, keypoint = output_to_predict(image, batch, output)\n            # rotated: image after rotation correction (H×W×3, uint8)\n            # keypoint: detected corner keypoints of the ECG paper in the image (the four corners of the paper i think: no the 9 start of the middle leads actually thats whhy the homography will be 3x3 ig)\n\t\t\t\t\tnormalised, keypoint, homo = normalise_by_homography(rotated, keypoint)\n            # ithink the new keypoint are just the first ones and we apply to the the same transformation the homo and normalized to heve there position in the normlized by homo photo \n\t\t\t# normalised: image after full homography correction (from the 9 keypoints)\n            # homo: the 3×3 homography matrix H\n            # Shape of normalised: always (1700, 2200, 3) — fixed canonical size\n\t\t\t\t\tcv2.imwrite(f'{OUT_DIR}/normalised/{sample_id}.norm.png', cv2.cvtColor(normalised, cv2.COLOR_RGB2BGR))\n\t\t\t\t\tnp.save(f'{OUT_DIR}/normalised/{sample_id}.homo.npy', homo)\n            # Save to disk: Stage 1 reads these files later\n            # Saving homo.npy: needed if you want to reconstruct original coordinates\n\t\t\t\texcept:\n\t\t\t\t\tFAIL_ID.append(sample_id)\n\n\t\ttorch.cuda.empty_cache()\n    # empty_cache(): releases GPU memory that PyTorch cached\n    # Important after failures: a failed forward pass might leave\n    # partial tensors on GPU → memory leak → OOM on next iteration\n\t\t#if n<10: # optional: show results\n\t\t\t#overlay = draw_results_stage0(rotated, keypoint)\n            # So the overlay is based on the keypoint set after that normalization step has processed it and there still compatible so the new keypoints work with non normlized images to just rotated ?\n\t\t\t#print('')\n\t\t\t#print('demo results for stage0--------------')\n\t\t\t#print(sample_id)\n\t\t\t#plt.imshow(image);plt.show()\n\t\t\t#plt.imshow(overlay);plt.show()\n\t\t\t#plt.imshow(normalised);plt.show()\n\t\t\t\n\tprint('')\n\nrun_stage0()\nprint('FAIL_ID:', FAIL_ID)\nprint('run_stage0() ok!!!\\n')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:22:41.451298Z","iopub.execute_input":"2026-05-02T16:22:41.451494Z","iopub.status.idle":"2026-05-02T16:31:22.607297Z","shell.execute_reply.started":"2026-05-02T16:22:41.451479Z","shell.execute_reply":"2026-05-02T16:31:22.606478Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 4.3 Inside `image_to_batch` — what preprocessing actually happens\n\n`image_to_batch(image)` does **not** directly apply ImageNet normalization, and it does **not** resize the image to a fixed shape such as `(512, 512)`. [file:144][file:145]  \nInstead, it first reads the original image size `(H, W)`, computes a scale factor `scale = WIDTH / W` with `WIDTH = 1440`, and resizes the image while keeping the aspect ratio. [file:144]  \nThe resize uses `cv2.INTER_AREA`, which is a common choice for shrinking images. [file:144]  \n\nAfter resizing, the function pads the image with zeros so that the padded height and width are slightly larger than the resized image and aligned to multiples of 32. [file:144]  \nThis is useful because the Stage 0 network is a convolutional encoder-decoder, and such architectures usually work more cleanly when spatial dimensions are compatible with repeated downsampling and upsampling. [file:145][file:144]  \n\nNext, the function converts the image from HWC format to CHW format and adds a batch dimension, so one image becomes a tensor of shape `(1, 3, H, W)`. [file:144]  \nThen it creates 4 test-time augmentation versions: the original image, a vertical flip, a horizontal flip, and a flip in both directions. [file:144]  \nThese 4 versions are concatenated into one batch, and the function also stores metadata such as the original size, resized size, and scale factor. [file:144]  \n\nSo the real role of `image_to_batch` is: resize with aspect ratio preserved, pad, convert to tensor layout, and build the TTA batch. [file:144]  \nThe actual ImageNet normalization happens later inside the model's `forward()` function, where the tensor is converted with `x = image.float() / 255` and then `x = (x - mean) / std`. [file:145]  \n\n---\n\n### 4.4 Inside `output_to_predict` — decoding the real Stage 0 outputs\n\nThe Stage 0 network is a two-head model. [file:145]  \nOne head predicts `marker`, which is a segmentation-like per-pixel class map, and the other head predicts `orientation`, which is an image-level rotation class. [file:145]  \nSo Stage 0 does **not** predict the 4 paper corners directly. [file:145][file:144]  \n\nMore precisely, the model defines `self.marker = nn.Conv2d(..., 13 + 1, kernel_size=1)` and `self.orientation = nn.Linear(..., 8)`. [file:145]  \nDuring inference, both outputs are passed through `softmax`, so `marker` becomes per-pixel class probabilities and `orientation` becomes class probabilities over 8 rotation-related classes. [file:145]  \n\nInside `output_to_predict`, the code loops over the 4 TTA images. [file:144]  \nFor each flipped input, it undoes the flip on the predicted `marker` map and also reorders the `orientation` probabilities so they match the original image reference frame before averaging. [file:144]  \nAfter that, it averages the 4 predictions into one final `marker` map and one final `orientation` vector. [file:144]  \n\nThen the function calls `marker_to_keypoint(image, orientation, marker, scale)`. [file:144]  \nInside `marker_to_keypoint`, the first step is `k = orientation.argmax()`, which means the image rotation is chosen **directly from the orientation head**, not computed later from the keypoints. [file:144][file:145]  \nThe predicted marker map is then rotated with `np.rot90(...)` so that landmark extraction happens in the predicted upright orientation. [file:144]  \n\nNext, the code converts the soft marker probabilities into a hard label image with `thresh = marker.argmax(-1)`. [file:144]  \nFor each of the 9 labels `[2, 3, 4, 6, 7, 8, 10, 11, 12]`, it finds connected components with `cc3d.connected_components(thresh == label)`, computes the component statistics, sorts the components by size, and keeps the centroid of the largest one. [file:144]  \nThat centroid becomes the final landmark coordinate for that class. [file:144]  \n\nSo the model does **not** directly output 9 `(x, y)` coordinates either. [file:144][file:145]  \nInstead, it outputs a marker map, and the code converts that map into 9 keypoints by taking one centroid per landmark class. [file:144]  \n\nThe 9 landmarks correspond to the lead-start positions for `aVR`, `V1`, `V4`, `aVL`, `V2`, `V5`, `aVF`, `V3`, and `V6`. [file:144]  \nTherefore, the Stage 0 keypoints are internal ECG-layout landmarks, not the four corners of the paper. [file:144]  \n\nFinally, `output_to_predict` rotates the original image with the chosen rotation `k` and returns `(rotated, keypoint)`. [file:144]  \nAt this moment, each keypoint has the form `[x, y, label, leadname]`. [file:144]  \n\nA useful intuition is this: the marker head paints small class-specific blobs near important ECG landmarks, and `cc3d` turns each blob into one clean center point. [file:144][file:145]  \n\n---\n\n### 4.5 Inside `normalise_by_homography` — the real geometric normalization\n\n`normalise_by_homography(image, keypoint)` does **not** use 4 paper corners. [file:144]  \nInstead, it extracts the 9 detected landmark coordinates with `pt9 = [[k[0], k[1]] for k in keypoint]` and sends them to `normalise_image(image, pt9)`. [file:144]  \n\nInside `normalise_image`, the code computes a homography with `cv2.findHomography(pt9, ref_pt9, method=cv2.RANSAC)`. [file:144]  \nSo the source points are the 9 detected landmarks in the current rotated ECG image, and the destination points are 9 canonical reference points `REF_PT9`. [file:144]  \nThese reference points are built beforehand from `640106434-0001.gridpoint_xy.npy`, then scaled and shifted into the Stage 0 canonical frame. [file:144]  \n\nThis means the Stage 0 homography maps the detected lead-start landmarks in the current image to their standard template locations. [file:144]  \nIn other words, Stage 0 normalizes the ECG image by aligning its internal landmark layout to a canonical ECG layout. [file:144]  \n\nThe call uses `method=cv2.RANSAC`, which is important because not every detected keypoint is guaranteed to be perfect. [file:144]  \nRANSAC tries to find one homography that is supported by the most geometrically consistent correspondences, while rejecting inconsistent ones as outliers. [file:144]  \nThis is safer than forcing all 9 points to agree exactly when one or more detections may be slightly wrong. [file:144]  \n\nThe function then warps the image with `cv2.warpPerspective(image, homo, (WIDTH, HEIGHT))`, where `WIDTH = 1440` and `HEIGHT = 1152`. [file:144]  \nSo the Stage 0 normalized image has shape `(1152, 1440, 3)`, not `(1700, 2200, 3)`. [file:144][file:95]  \nThe larger canonical geometry such as `1700 x 2200` appears elsewhere in the pipeline context, but the actual Stage 0 normalization code in `stage0_common.py` outputs `(1152, 1440)`. [file:144][file:95]  \n\n`cv2.findHomography` also returns a `match` mask. [file:144]  \nThis mask tells which of the 9 correspondences were considered inliers to the final homography estimated by RANSAC. [file:144]  \n`normalise_by_homography` appends this flag to each keypoint, so each keypoint changes from `[x, y, label, leadname]` to `[x, y, label, leadname, match]`. [file:144]  \n\nThis is why `draw_results_stage0(rotated, keypoint)` still works correctly on the rotated image. [file:144]  \nThe keypoint coordinates are still in the rotated-image coordinate system, and the only new information added is whether each point was an inlier or outlier for the homography fit. [file:144]  \nIn the visualization, matched points are drawn as filled circles with the lead name written next to them, while unmatched points are drawn as outlined circles. [file:144]  \n\nSo the correct Stage 0 interpretation is: predict orientation, predict a marker map, extract 9 lead-start landmarks with connected components, fit a homography from those 9 detected points to 9 canonical template points using RANSAC, and warp the image into a standard Stage 0 normalized frame. [file:144][file:145]  \n\nWhy heatmaps instead of directly regressing (x,y) coordinates?\n-\n\nThis is a fundamental architectural choice that appears in many pose estimation problems (human body keypoints, face landmarks, etc.).\n\nDirect regression approach:\n\ntext\nCNN → flatten → dense layer → [x, y]  (2 numbers)\nProblem: the network must compress all spatial information into 2 numbers through a bottleneck. If the corner is in the top-right of the image, the network must somehow encode \"top-right\" as a single number. This is hard to learn because there's no spatial structure in the output.\n\nHeatmap approach:\n\nCNN → UNet decoder → heatmap (H_feat × W_feat)  (spatial output)\nThe network outputs a 2D map where the peak is at the corner location. The spatial structure of the output directly mirrors the spatial structure of the input. The network just needs to learn \"draw a peak where you see a corner mark\" — much more natural for a convolutional architecture.\n\nWhy does this work so well with CNNs?\nConvolutions are inherently spatial operations. Each spatial position in the feature map \"knows\" what's happening at the corresponding location in the input image (within its receptive field). Outputting a heatmap leverages this spatial awareness directly. Outputting a single (x,y) via a dense layer destroys this spatial structure.\n\n The heatmap approach preserves the topological structure of the problem — nearby pixels in the input correspond to nearby outputs. Direct regression via dense layers breaks this topology completely (all spatial information is mixed in the flattened vector).\n\n-\n\n ","metadata":{"_kg_hide-input":false}},{"cell_type":"markdown","source":"Why type 0001 images for training Stage 1?\n-\nType 0001 is the synthetic original — perfect image, no noise, no distortion, no color cast. The gridlines are perfectly sharp and perfectly positioned. This gives pixel-perfect ground truth. The network then generalizes to other types (3, 4, 5, etc.) because after Stage 0 normalisation, all types look somewhat similar (roughly aligned, similar scale). The remaining differences (noise, blur, color) don't prevent gridline detection — the network learned to find gridlines in clean images and generalizes well enough to noisy ones.\n\nUnderstanding the homography matrix H mathematically:\n\nA homography maps points between two projective planes. In homogeneous coordinates:\n\ntext\n[x']     [h00  h01  h02]   [x]\n[y']  =  [h10  h11  h12] × [y]\n[w']     [h20  h21  h22]   [1]\n\nThen: x_canonical = x'/w',  y_canonical = y'/w'\nThe division by w' is what makes it a projective (not just affine) transformation. When w'=1 everywhere, it reduces to an affine transformation. When w' varies spatially, you get perspective effects — lines remain straight but parallel lines can converge (like railroad tracks vanishing to a point).\n\nWhy exactly 4 points to determine H?\n\nH has 9 elements but only 8 degrees of freedom (the overall scale doesn't matter — you can multiply all elements by any constant). Each point correspondence gives 2 equations (x and y). So 4 points × 2 equations = 8 equations → exactly determines all 8 degrees of freedom. With fewer points the system is underdetermined. With more points you'd use RANSAC (as in Stage 1).\n\ncv2.warpPerspective — what it actually does pixel by pixel:\n\ntext\nFor each output pixel (x_out, y_out) in canonical space:\n    1. Compute source position: [x_src, y_src, w] = H_inv @ [x_out, y_out, 1]\n    2. Normalize: x_src /= w,  y_src /= w\n    3. Sample the rotated image at (x_src, y_src) using bilinear interpolation\n    4. Place that color at output position (x_out, y_out)\nThis is called inverse warping (you go from output to input to find what to sample). The alternative, forward warping (map each input pixel to its output position), creates holes (multiple input pixels might map to the same output, or some outputs get nothing). Inverse warping fills every output pixel exactly once — no holes, no double-filling.\n\nBilinear interpolation at step 3:\nThe source position is usually non-integer (e.g., 347.3, 892.7). Bilinear interpolation blends the 4 surrounding pixels using weights proportional to the fractional distances:\n\ntext\npixel(347.3, 892.7) =\n    0.7 × 0.3 × pixel(347, 892) +   (weight = (1-0.3)×(1-0.7))\n    0.7 × 0.7 × pixel(347, 893) +   (weight = (1-0.3)×0.7)\n    0.3 × 0.3 × pixel(348, 892) +   (weight = 0.3×(1-0.7))\n    0.3 × 0.7 × pixel(348, 893)     (weight = 0.3×0.7)\nThis gives smooth, anti-aliased results. Compare with nearest-neighbor interpolation (just round to nearest integer pixel) which produces blocky/jagged edges.","metadata":{"_kg_hide-input":false}},{"cell_type":"code","source":"# stage1\nprint('*** STARTING STAGE1 ***')\n\nfrom stage1_model import Net as Stage1Net\nfrom stage1_common import *\n# output_to_predict (different from stage0's version) , rectify_image,draw_mapping,draw_results_stage1\nos.makedirs(f'{OUT_DIR}/rectified', exist_ok=True)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:31:22.608281Z","iopub.execute_input":"2026-05-02T16:31:22.608674Z","iopub.status.idle":"2026-05-02T16:31:22.631466Z","shell.execute_reply.started":"2026-05-02T16:31:22.608648Z","shell.execute_reply":"2026-05-02T16:31:22.630929Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#gridpoint sparse mask improvment from solution 1 : THE GRIDPOINTS WEIGHTED AVERGE THAT HAS SUB PIXEL PRESISION ( after the model output we just dont use cc3d alone we averge around the confident detections ) THE SURFACE FITTING POLYNOMIALE( so the detetected points make sens together more)\n\ndef _poly_design(ii, jj):\n    ii = np.asarray(ii, np.float32)\n    jj = np.asarray(jj, np.float32)\n    return np.stack([\n        np.ones_like(ii),\n        ii,\n        jj,\n        ii * ii,\n        jj * jj,\n        ii * jj,\n    ], axis=1)\n# For each detected blob, it computes a weighted average of all pixel positions inside it, using the heatmap probability value at each pixel as the weight.\n\ndef weighted_cc_centroids(prob, thresh=0.5, min_area=3, eps=1e-8):\n    cc = cc3d.connected_components(prob > thresh)\n    pts = []\n    for k in range(1, cc.max() + 1):\n        ys, xs = np.where(cc == k)\n        if len(xs) < min_area:\n            continue\n\n        w = prob[ys, xs].astype(np.float64)\n        s = w.sum()\n        if s < eps:\n            continue\n\n        y = (ys * w).sum() / (s + eps)\n        x = (xs * w).sum() / (s + eps)\n        pts.append([y, x])\n\n    if len(pts) == 0:\n        return np.zeros((0, 2), np.float32)\n\n    return np.asarray(pts, np.float32)\n\ndef subpixel_point_to_label(y, x, hcc, vcc):\n    H, W = hcc.shape\n\n    y0 = int(np.clip(np.floor(y), 0, H - 1))\n    x0 = int(np.clip(np.floor(x), 0, W - 1))\n    y1 = min(y0 + 1, H - 1)\n    x1 = min(x0 + 1, W - 1)\n\n    candidates = [\n        (hcc[y0, x0], vcc[y0, x0]),\n        (hcc[y0, x1], vcc[y0, x1]),\n        (hcc[y1, x0], vcc[y1, x0]),\n        (hcc[y1, x1], vcc[y1, x1]),\n    ]\n    candidates = [(j, i) for j, i in candidates if j > 0 and i > 0]\n\n    if len(candidates) == 0:\n        return 0, 0\n\n    vals, counts = np.unique(np.asarray(candidates), axis=0, return_counts=True)\n    j, i = vals[counts.argmax()]\n    return int(j), int(i)\n\ndef refine_grid_surface(gridpoint_xy, sigma=3.0, fill_missing=True):\n    out = gridpoint_xy.copy()\n\n    valid = np.any(out != 0, axis=2)\n    jj, ii = np.where(valid)  # row, col on canonical grid\n\n    if len(ii) < 12:\n        return out\n\n    A = _poly_design(ii, jj).astype(np.float64)\n    tx = out[jj, ii, 0].astype(np.float64)\n    ty = out[jj, ii, 1].astype(np.float64)\n\n    cx, *_ = np.linalg.lstsq(A, tx, rcond=None)\n    cy, *_ = np.linalg.lstsq(A, ty, rcond=None)\n\n    px = A @ cx\n    py = A @ cy\n\n    resid = np.sqrt((tx - px) ** 2 + (ty - py) ** 2)\n    thr = resid.mean() + sigma * (resid.std() + 1e-6)\n    bad = resid > thr\n\n    if bad.any():\n        out[jj[bad], ii[bad], 0] = px[bad].astype(np.float32)\n        out[jj[bad], ii[bad], 1] = py[bad].astype(np.float32)\n\n    if fill_missing:\n        jg, ig = np.meshgrid(\n            np.arange(out.shape[0], dtype=np.float32),\n            np.arange(out.shape[1], dtype=np.float32),\n            indexing='ij'\n        )\n        Afull = _poly_design(ig.reshape(-1), jg.reshape(-1)).astype(np.float64)\n        fullx = (Afull @ cx).reshape(out.shape[:2]).astype(np.float32)\n        fully = (Afull @ cy).reshape(out.shape[:2]).astype(np.float32)\n\n        miss = ~valid\n        out[..., 0][miss] = fullx[miss]\n        out[..., 1][miss] = fully[miss]\n\n    return out\n\ndef output_to_predict_refined(image, batch, output, interpolate_mapping_fn,\n                              segment_to_endpoints_fitline_fn,\n                              canonical_x_order_fn,\n                              canonical_y_order_fn,\n                              compare_segment_fn):\n    marker = output['marker'][0]\n    gridpoint = output['gridpoint'][0, 0]\n    gridhline = output['gridhline'][0]\n    gridvline = output['gridvline'][0]\n\n    marker = marker.argmax(0).byte().data.cpu().numpy()\n    gridpoint = gridpoint.float().data.cpu().numpy()\n    gridhline = gridhline.float().data.cpu().numpy()\n    gridvline = gridvline.float().data.cpu().numpy()\n\n    gridvline = gridvline.argmax(0).astype(np.uint8)\n    gridhline = gridhline.argmax(0).astype(np.uint8)\n\n    # 1) sub-pixel point locations: weighted centroid instead of plain CC centroid\n    point_yx = weighted_cc_centroids(gridpoint, thresh=0.5, min_area=3)\n\n    # 2) same line filtering logic as baseline\n    gvfiltered = np.zeros_like(gridvline)\n    cc = cc3d.connected_components(gridvline != 0)\n    num_line = cc.max()\n    for l in range(1, num_line + 1):\n        t = (cc == l)\n        bincount = np.bincount(gridvline[t])\n        c = bincount.argmax()\n        gvfiltered[t] = c\n\n    gvreject = np.zeros_like(gridvline)\n    for l in range(1, num_line + 1):\n        cc2 = cc3d.connected_components(gvfiltered == l)\n        if cc2.max() > 1:\n            num = cc2.max() + 1\n            stats = cc3d.statistics(cc2)\n            area = stats['voxel_counts'][1:]\n            label = np.arange(1, num)\n\n            argsort = np.argsort(area)[::-1]\n            area = area[argsort]\n            label = label[argsort]\n\n            if area[0] < 7:\n                continue\n\n            main_segment = segment_to_endpoints_fitline_fn(cc2 == label[0])\n            main_segment = canonical_y_order_fn(*main_segment)\n\n            for j in range(1, len(label)):\n                if area[j] < 7:\n                    continue\n                segment = segment_to_endpoints_fitline_fn(cc2 == label[j])\n                segment = canonical_y_order_fn(*segment)\n                ang_dis, ori_dis, seg_dis = compare_segment_fn(main_segment, segment)\n\n                if ori_dis > 5:\n                    gvreject[cc2 == label[j]] = 255\n                else:\n                    gvfiltered[cc2 == label[j]] = l\n\n    vcc = gvfiltered.copy()\n\n    ghfiltered = np.zeros_like(gridhline)\n    cc = cc3d.connected_components(gridhline != 0)\n    num_line = cc.max()\n    for l in range(1, num_line + 1):\n        t = (cc == l)\n        bincount = np.bincount(gridhline[t])\n        c = bincount.argmax()\n        ghfiltered[t] = c\n\n    ghreject = np.zeros_like(gridhline)\n    for l in range(1, num_line + 1):\n        cc2 = cc3d.connected_components(ghfiltered == l)\n        if cc2.max() > 1:\n            num = cc2.max() + 1\n            stats = cc3d.statistics(cc2)\n            area = stats['voxel_counts'][1:]\n            label = np.arange(1, num)\n\n            argsort = np.argsort(area)[::-1]\n            area = area[argsort]\n            label = label[argsort]\n\n            if area[0] < 7:\n                continue\n\n            main_segment = segment_to_endpoints_fitline_fn(cc2 == label[0])\n            main_segment = canonical_x_order_fn(*main_segment)\n\n            for j in range(1, len(label)):\n                if area[j] < 7:\n                    continue\n                segment = segment_to_endpoints_fitline_fn(cc2 == label[j])\n                segment = canonical_x_order_fn(*segment)\n                ang_dis, ori_dis, seg_dis = compare_segment_fn(main_segment, segment)\n\n                if ori_dis > 5:\n                    ghreject[cc2 == label[j]] = 255\n                else:\n                    ghfiltered[cc2 == label[j]] = l\n\n    hcc = ghfiltered.copy()\n\n    # 3) build sparse (44,57,2) grid using sub-pixel points\n    gridpoint_xy = np.zeros((44, 57, 2), np.float32)\n\n    for y, x in point_yx:\n        j, i = subpixel_point_to_label(y, x, hcc, vcc)\n        if (j == 0) or (i == 0):\n            continue\n\n        jj = j - 1\n        ii = i - 1\n\n        if np.all(gridpoint_xy[jj, ii] == 0):\n            gridpoint_xy[jj, ii] = [x, y]\n        else:\n            oy, ox = int(round(gridpoint_xy[jj, ii, 1])), int(round(gridpoint_xy[jj, ii, 0]))\n            ny, nx = int(round(y)), int(round(x))\n\n            oy = np.clip(oy, 0, gridpoint.shape[0] - 1)\n            ox = np.clip(ox, 0, gridpoint.shape[1] - 1)\n            ny = np.clip(ny, 0, gridpoint.shape[0] - 1)\n            nx = np.clip(nx, 0, gridpoint.shape[1] - 1)\n\n            old_score = gridpoint[oy, ox]\n            new_score = gridpoint[ny, nx]\n            if new_score > old_score:\n                gridpoint_xy[jj, ii] = [x, y]\n\n    # 4) 1st-place-style surface fitting: replace outliers on the sparse lattice\n    gridpoint_xy = refine_grid_surface(gridpoint_xy, sigma=3.0, fill_missing=False)\n\n    # 5) keep baseline hole filling\n    gridpoint_xy = interpolate_mapping_fn(gridpoint_xy)\n\n    more = {\n        'ghfiltered': ghfiltered,\n        'gvfiltered': gvfiltered,\n    }\n    return gridpoint_xy, more","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:31:22.632156Z","iopub.execute_input":"2026-05-02T16:31:22.632319Z","iopub.status.idle":"2026-05-02T16:31:22.658014Z","shell.execute_reply.started":"2026-05-02T16:31:22.632306Z","shell.execute_reply":"2026-05-02T16:31:22.657404Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"| Step                    | Action                                                        | Why                                                 |\n| ----------------------- | ------------------------------------------------------------- | --------------------------------------------------- |\n| weighted_cc_centroids   | Get sub-pixel (y, x) for each detected grid point blob        | Reduces initial quantization error                  |\n| subpixel_point_to_label | Assign each point to its (row, col) in the 44×57 grid         | Maps detections to canonical grid positions         |\n| refine_grid_surface     | Fit polynomial, flag 3σ outliers, replace with surface values | Removes bad grid points that would warp incorrectly |\n| interpolate_mapping_fn  | Fill any remaining holes                                      | Inherited from the baseline                         |","metadata":{"_kg_hide-input":false}},{"cell_type":"code","source":"# \ndef run_stage1():\n    stage1_net = Stage1Net(pretrained=False)\n    stage1_net = load_net(stage1_net, f'{WEIGHT_DIR}/stage1-last.checkpoint.pth')\n    stage1_net.to(DEVICE)\n\n    start_timer = timer()\n    for n, sample_id in enumerate(valid_id):\n        timestamp = time_to_str(timer() - start_timer, 'sec')\n        print(f'\\r\\t {n:4d} {sample_id}', timestamp, end='', flush=True)\n\n        if sample_id in FAIL_ID:\n            continue\n\n        image = cv2.imread(f'{OUT_DIR}/normalised/{sample_id}.norm.png', cv2.IMREAD_COLOR_RGB)\n\n        batch = {\n            'image': torch.from_numpy(\n                np.ascontiguousarray(image.transpose(2, 0, 1))\n            ).unsqueeze(0),\n        }\n\n        num_tta = 1\n\n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                output = stage1_net(batch)\n\n        try:\n            gridpoint_xy, more = output_to_predict_refined(\n                image=image,\n                batch=batch,\n                output=output,\n                interpolate_mapping_fn=interpolate_mapping,\n                segment_to_endpoints_fitline_fn=segment_to_endpoints_fitline,\n                canonical_x_order_fn=canonical_x_order,\n                canonical_y_order_fn=canonical_y_order,\n                compare_segment_fn=compare_segment,\n            )\n           # gridpoint_xy, more = output_to_predict(image, batch, output)\n            rectified = rectify_image(image, gridpoint_xy)\n\n            cv2.imwrite(\n                f'{OUT_DIR}/rectified/{sample_id}.rect.png',\n                cv2.cvtColor(rectified, cv2.COLOR_RGB2BGR)\n            )\n            np.save(f'{OUT_DIR}/rectified/{sample_id}.gridpoint_xy.npy', gridpoint_xy)\n\n        except Exception as e:\n            print(f'\\nStage1 failed for {sample_id}: {e}')\n            FAIL_ID.append(sample_id)\n\n        torch.cuda.empty_cache()\n\n       # if n < 10:\n           # overlay = draw_mapping(image, gridpoint_xy)\n           # ghfiltered, gvfiltered = draw_results_stage1(more)\n\n          #  print('')\n           # print('demo results for stage1--------------')\n            #print(sample_id)\n           # plt.imshow(overlay)\n           # plt.show()\n           # plt.imshow(gvfiltered)\n           # plt.show()\n           # plt.imshow(ghfiltered)\n            #plt.show()\n           # plt.imshow(rectified)\n           # plt.show()\n\n    print('')\n\n\nrun_stage1()\nprint('FAIL_ID:', FAIL_ID)\nprint('run_stage1() ok!!!\\n')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T16:31:22.660460Z","iopub.execute_input":"2026-05-02T16:31:22.660670Z","iopub.status.idle":"2026-05-02T19:10:07.228234Z","shell.execute_reply.started":"2026-05-02T16:31:22.660654Z","shell.execute_reply":"2026-05-02T19:10:07.227601Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 5.3 The Stage 1 Neural Network Architecture\n\nStage 1 takes as input the **Stage 0 normalized image**, whose size is fixed at `(1152, 1440, 3)`. [file:144][file:146][file:147]  \nThe model then applies ImageNet-style normalization inside its `forward()` function with `x = image.float() / 255` followed by `(x - mean) / std`. [file:147]  \n\nStage 1 is not just a \"horizontal heatmap + vertical heatmap\" model. [file:147][file:146]  \nIt has **four output heads**: `gridpoint` with 1 channel, `gridhline` with `44+1=45` channels, `gridvline` with `57+1=58` channels, and `marker` with `13+1=14` channels. [file:147]  \nSo the network predicts candidate grid intersections, a per-pixel horizontal-line class, a per-pixel vertical-line class, and also a marker map. [file:147]  \n\nThe backbone is `resnet34.a3_in1k`, created with `num_classes=0` and `global_pool=''`, so the classifier head is removed and spatial feature maps are preserved. [file:147]  \nThe encoder feature dimensions are `[64, 128, 256, 512]`, and the decoder dimensions are `[256, 128, 64, 32]`. [file:147]  \nThe decoder is `MyUnetDecoder`, which is a standard UNet-style decoder with interpolation, skip connections, and convolution blocks; there is no coordinate-convolution module in this code. [file:147]  \n\nMore concretely, the decoder upsamples the deepest encoder feature map in four stages and combines it with skip features from earlier encoder levels. [file:147]  \nThe final decoder output has 32 channels, and the four prediction heads are applied with `1x1` convolutions on top of that shared feature map. [file:147]  \nThis means Stage 1 is a **multi-head dense prediction model**, not a model with a single 4-channel \"pixel\" output. [file:147]  \n\nAt inference time, `output_to_predict` converts `gridpoint` into a binary candidate-point map by thresholding `gridpoint > 0.5`, then extracts connected components and uses their centroids as candidate intersection locations. [file:146]  \nFor `gridhline` and `gridvline`, it does `argmax` over channels, so these are not simple heatmaps of \"line present or absent\"; they are pixelwise class maps that try to assign each pixel to one of the 44 horizontal lines or 57 vertical lines, plus background. [file:146][file:147]  \nThe code then cleans and links line segments using connected components, line fitting, and geometric consistency checks before keeping filtered horizontal and vertical line-label maps. [file:146]  \n\nAfter that, each candidate point is assigned a horizontal line index `j` and a vertical line index `i` by reading the filtered maps at that point location. [file:146]  \nThe result is stored in `gridpoint_xy` with shape `(44, 57, 2)`, where `gridpoint_xy[j-1, i-1] = [x, y]` gives the image-space coordinate of that canonical grid intersection. [file:146]  \nMissing entries are then filled by `interpolate_mapping(gridpoint_xy)` using cubic interpolation with `scipy.interpolate.griddata`. [file:146]  \n\nSo the main goal of Stage 1 is not merely to detect isolated intersections. [file:146][file:147]  \nIts real goal is to recover a **dense mapping** from the canonical ECG grid to the current image, using point candidates together with horizontal and vertical line identities. [file:146]  \nThat mapping is exactly what makes precise geometric rectification possible. [file:146]  \n\nFinally, `rectify_image(image, gridpoint_xy)` converts the sparse `(44, 57, 2)` mapping into a dense sampling map, upsamples it to a reference size `(1700, 2200)`, and uses `torch.nn.functional.grid_sample` to warp the image. [file:146]  \nSo the Stage 1 **output image** is rectified to `(1700, 2200)`, but the Stage 1 **input image** is still the Stage 0 normalized image of size `(1152, 1440)`. [file:146][file:147]  \nThis is why saying \"Stage 1 takes a 1700x2200 image\" is incorrect; 1700x2200 is the rectified output size, not the network input size. [file:146][file:147]\nthe training target used is the 001 scaled and added output to each cordinate to match the stage 0 output and stage 1 working space , when we do it for most photos the coordinantes of grid and all are veeery colse only rounding probs maybe , but for the mobile photos if its not perfectly centerd in the corners we may not have perfect training target but we would have moset good cause we did a homography rectification before in stage0 and the ones that are bad and not guussed even in the corners would be calculated by the interpolation to reach the number cause as long as we have a good number of grid pointes detecterd the rest can be deduced by cubic interpolating i guess . for the 001 we got every infos there iso. the transformation is 280 / 2200 and add (-6, +10) for each x,y in canonical . \nimprvment till now to add : per type preprocess more tta in stage 1 and subpixel targets in stage 1 \n\nWhy type 0001 images for training Stage 1?\n-\nType 0001 is the synthetic original — perfect image, no noise, no distortion, no color cast. The gridlines are perfectly sharp and perfectly positioned. This gives pixel-perfect ground truth. The network then generalizes to other types (3, 4, 5, etc.) because after Stage 0 normalisation, all types look somewhat similar (roughly aligned, similar scale). The remaining differences (noise, blur, color) don't prevent gridline detection — the network learned to find gridlines in clean images and generalizes well enough to noisy ones.\n\ncc3d role :\n-\n After sigmoid: gh is a (850, 1100) float map with values in [0,1]\n Threshold to binary\ngh_binary = (gh_raw > 0.5).astype(np.uint8)\n\n Label connected components\nlabels = cc3d.connected_components(gh_binary, connectivity=4)\n Each connected blob gets a unique integer label\n A genuine horizontal gridline → long, thin horizontal blob\n Noise / false positives → small, compact blobs\n\n Filter: keep only blobs with width >> height (horizontal lines)\nfor label_id in np.unique(labels):\n    blob = (labels == label_id)\n    blob_width  = blob.any(axis=0).sum()   # number of columns the blob spans\n    blob_height = blob.any(axis=1).sum()   # number of rows the blob spans\n    if blob_width < 50 or blob_height > 20:\n        gh_binary[blob] = 0   # remove non-line blobs\n Result: only genuine horizontal lines remain\n Why connected components instead of just thresholding?\nSimple thresholding removes weak detections. Connected component analysis removes geometrically wrong detections — even confident ones. A false positive with activation 0.8 but spanning only 10 columns (not a full-width line) gets removed. A genuine gridline spanning 800 columns gets kept. Geometry matters more than raw confidence.\n\nWhy cc3d instead of scipy.ndimage.label?\ncc3d is dramatically faster for large binary arrays (1700×2200). It uses a highly optimized C++ implementation with parallel processing. For the Stage 1 heatmaps at half resolution (850×1100), scipy.ndimage.label might take several seconds. cc3d does it in milliseconds. In a competition with runtime limits, this difference matters across thousands of test images.\n\nRMQ: The gridpointxy.npy File — Why It's Saved\nThis saved file is used by multiple downstream solutions:\n\nhengck23 Stage 2: Reads the file to know the zero_mv pixel positions (where the ECG baseline is in the rectified image) and the time span t0, t1 (which columns of the rectified image correspond to the start and end of the ECG signal).\n\n2nd place solution: Reads gridpointxy.npy directly to perform their own refined warping (they trust hengck23's grid detection but apply their own coordinate transformations on top, mapping from the competition's resampled signals to PTB-XL's 500Hz original).","metadata":{"_kg_hide-input":false}},{"cell_type":"code","source":"# stage2\nprint('*** STARTING STAGE2 ***')\n\nimport torch.nn as nn\nfrom stage2_smp_model import Net as WholeModel\nfrom stage2_lead_model import Net as LeadModel\nfrom stage2_model import prob_to_series_by_max   # keep this, still used\nfrom stage2_common import *\n\nos.makedirs(f'{OUT_DIR}/digitalised', exist_ok=True)\n#os.makedirs(f'{OUT_DIR}/debug', exist_ok=True)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T19:10:07.228970Z","iopub.execute_input":"2026-05-02T19:10:07.229152Z","iopub.status.idle":"2026-05-02T19:10:07.309131Z","shell.execute_reply.started":"2026-05-02T19:10:07.229138Z","shell.execute_reply":"2026-05-02T19:10:07.308519Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#global constants:\nfrom scipy import signal   # add to top imports if not present\n\nWINDOW_SIZE  = 240\nOFFSET       = 416\nIGNORE_EDGE  = 8\nx_scale      = 5000 / (2080 - 118)\nadd_x        = 1\ny_scale      = 1\nIMG_H, IMG_W = int(1700 * y_scale), int(2200 * x_scale) + add_x\n\ntta = [0]   # test-time augmentation: original + horizontal flip\n\nx0, x1 = 0, 5600\ny0, y1 = 0, 1696\nzero_mv = [703.5, 987.5, 1271.5, 1531.5]        # pixel row of 0mV per row\nzero_mv_trimed  = [pos - OFFSET for pos in zero_mv]\nzero_mv_croped  = [WINDOW_SIZE + 0.5 for _ in range(4)]\nmv_to_pixel     = 79.0\nt0, t1 = int(118 * x_scale) + add_x, int(2080 * x_scale) + add_x\n\n# Ensemble regions: for each lead row, which part of the trimmed image\n# and which part of the cropped lead image to fuse\nheight_after_trimed = y1 - OFFSET\nens_regions = []\nfor zmv in zero_mv_trimed:\n    trim_upper = int(zmv) - WINDOW_SIZE\n    trim_lower = int(zmv) + WINDOW_SIZE\n    lead_upper = IGNORE_EDGE\n    lead_lower = -IGNORE_EDGE\n    if trim_lower > height_after_trimed:\n        lead_lower = (trim_lower - height_after_trimed + IGNORE_EDGE) * -1\n        trim_lower = height_after_trimed\n    trim_upper += IGNORE_EDGE\n    trim_lower -= IGNORE_EDGE\n    ens_regions.append([trim_upper, trim_lower, lead_upper, lead_lower])\n#helper fonction : pixel_to_series_exp — this is the subpixel conversion. It replaces hengck23's pixel_to_series (which used argmax):\ndef pixel_to_series_exp(pixel, zero_mv, length):\n    \"\"\"Subpixel-precision pixel→series via weighted row expectation.\"\"\"\n    _, H, W = pixel.shape\n    eps = 1e-8\n    y_idx = np.arange(H, dtype=np.float32)[:, None]   # (H,1)\n    series = []\n    for j in [0, 1, 2, 3]:\n        p = pixel[j]\n        denom = p.sum(axis=0)           # (W,)\n        y_exp = (p * y_idx).sum(axis=0) / (denom + eps)   # soft expectation\n        series.append(y_exp)\n    series = np.stack(series).astype(np.float32)\n    if length is not None and length != W:\n        series = np.stack([\n            signal.resample(s, length).astype(np.float32) for s in series\n        ])\n    return series\n#read_images — loads the rectified image and builds the 4 cropped lead images for LeadModel:\ndef read_images(path):\n    image = cv2.imread(path, cv2.IMREAD_COLOR)\n    image = cv2.cvtColor(image, cv2.COLOR_BGR2RGB)\n    image = cv2.resize(image, (IMG_W, IMG_H), interpolation=cv2.INTER_LINEAR)\n    trim_image = image.copy()[OFFSET:y1, x0:x1]\n\n    image = image[y0:y1, x0:x1]\n    H, W, _ = image.shape\n    lead_images = []\n    for i, zmv in enumerate(zero_mv):\n        h0, h1 = int(zmv) - WINDOW_SIZE, int(zmv) + WINDOW_SIZE\n        src_h0, src_h1 = max(0, h0), min(H, h1)\n        dst_h0 = src_h0 - h0\n        dst_h1 = dst_h0 + (src_h1 - src_h0)\n        lead_img = np.zeros((WINDOW_SIZE * 2, W, 3))\n        lead_img[dst_h0:dst_h1, :, :] = image[src_h0:src_h1, :, :]\n        lead_images.append(lead_img)\n    lead_images = np.stack(lead_images)   # (4, H, W, 3)\n    return trim_image, lead_images    \n#loader fonctions    \ndef get_whole_model(encoder_name, weight_path, device):\n    model = WholeModel(encoder_name=encoder_name, encoder_weights=None,\n                       decoder_name=\"unet\", use_coord_conv=True, pretrained=False)\n    state_dict = torch.load(weight_path, map_location=lambda s, l: s)\n    print(model.load_state_dict(state_dict, strict=False))\n    model.to(device); model.eval(); model.output_type = ['infer']\n    return model\n\ndef get_lead_model(encoder_name, weight_path, fusion_type, device):\n    model = LeadModel(encoder_name=encoder_name, encoder_weights=None,\n                      fusion_type=fusion_type)\n    state_dict = torch.load(weight_path, map_location=lambda s, l: s)\n    print(model.load_state_dict(state_dict, strict=False))\n    model.to(device); model.eval(); model.output_type = ['infer']\n    return model\ndef get_attention_lead_model(weight_path, device):\n    model = LeadModel(\n        encoder_name='tu-efficientnet_b6',\n        encoder_weights=None,\n        fusion_type='cross_attn',\n        fusion_levels=[3, 4],\n    )\n    state = torch.load(weight_path, map_location='cpu')\n    msg = model.load_state_dict(state, strict=False)\n    print(msg)\n    model.to(device)\n    model.eval()\n    model.output_type = ['infer']\n    return model","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T19:10:07.309953Z","iopub.execute_input":"2026-05-02T19:10:07.310177Z","iopub.status.idle":"2026-05-02T19:10:07.323689Z","shell.execute_reply.started":"2026-05-02T19:10:07.310163Z","shell.execute_reply":"2026-05-02T19:10:07.322967Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def run_stage2(gpu_id=0, assigned_ids=None, prev_fail_ids=None, fail_id_file=None):\n    device = f'cuda:{gpu_id}'\n    if assigned_ids is None: assigned_ids = valid_id\n    if prev_fail_ids is None: prev_fail_ids = []\n    local_fail_id = []\n\n    # Load all models once\n    whole_models = [\n        get_whole_model(\"tu-timm/tf_efficientnet_b7.ns_jft_in1k\",\n            \"/kaggle/input/datasets/takashisomeya/physionet-final-submission-models/whole_b7_lb22.93.pth\", device),\n        get_whole_model(\"tu-timm/tf_efficientnetv2_l.in21k\",\n            \"/kaggle/input/datasets/takashisomeya/physionet-final-submission-models/whole_v2_l_lb22.60.pth\", device),\n    ]\n    lead_models = [\n        get_lead_model(\"tu-timm/tf_efficientnet_b6.ns_jft_in1k\",\n            \"/kaggle/input/datasets/takashisomeya/physionet-final-submission-models/series_b6_shared_conv2d_lb23.10.pth\", \"shared_conv2d\", device),\n        get_lead_model(\"tu-timm/tf_efficientnet_b6.ns_jft_in1k\",\n            \"/kaggle/input/datasets/takashisomeya/physionet-final-submission-models/series_b6_shared_conv2d_lb23.00.pth\", \"shared_conv2d\", device),\n        get_lead_model(\"tu-timm/tf_efficientnetv2_l.in21k\",\n            \"/kaggle/input/datasets/takashisomeya/physionet-final-submission-models/series_v2_l_conv3d_lb22.92.pth\", \"conv3d\", device),\n        get_lead_model(\"tu-timm/tf_efficientnetv2_l.in21k\",\n            \"/kaggle/input/datasets/takashisomeya/physionet-final-submission-models/series_v2_l_conv2d_lb22.85.pth\", \"conv2d\", device),\n        get_attention_lead_model(\"/kaggle/input/datasets/zahouaniyacine/attention-fusion-model-ecg/cross_attn_b6_final_1.pth\", device),\n    ]\n\n    start_timer = timer()\n    for n, sample_id in enumerate(assigned_ids):\n        timestamp = time_to_str(timer() - start_timer, 'sec')\n        print(f'\\r\\t {n:4d} {sample_id}', timestamp, end='', flush=True)\n        if sample_id in prev_fail_ids: continue\n\n        length = read_sampling_length(sample_id)\n        trim_image, lead_images = read_images(f'{OUT_DIR}/rectified/{sample_id}.rect.png')\n\n        pixel_ens = np.zeros((4, trim_image.shape[0], trim_image.shape[1]))\n\n        # --- Whole model inference (full trimmed image) ---\n        batch = {'image': torch.from_numpy(\n            np.ascontiguousarray(trim_image.transpose(2,0,1))).unsqueeze(0)}\n        batch_tta = {'image': torch.from_numpy(\n            np.ascontiguousarray(np.fliplr(trim_image).copy().transpose(2,0,1))).unsqueeze(0)}\n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                for model in whole_models:\n                    for flip in tta:\n                        if flip:\n                            pixel = model(batch_tta)['pixel'].float().cpu().numpy()[0]\n                            pixel = np.flip(pixel, axis=flip)\n                        else:\n                            pixel = model(batch)['pixel'].float().cpu().numpy()[0]\n                        pixel_ens += pixel\n\n        # --- Lead model inference (per-lead cropped windows) ---\n        lead_t = torch.from_numpy(lead_images.transpose(0,3,1,2)).contiguous()\n        batch = {'image': lead_t.unsqueeze(0)}\n        batch_tta = {'image': torch.flip(lead_t, dims=[3]).unsqueeze(0)}\n        with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n            with torch.no_grad():\n                for model in lead_models:\n                    for flip in tta:\n                        use_batch = batch_tta if flip else batch\n                \n                        if getattr(model, \"fusion_type\", \"\") == \"cross_attn\":\n                            with torch.no_grad():\n                                with torch.amp.autocast('cuda', enabled=False):\n                                    out = model({'image': use_batch['image'].float()})\n                        else:\n                            with torch.no_grad():\n                                with torch.amp.autocast('cuda', dtype=FLOAT_TYPE):\n                                    out = model(use_batch)\n                \n                        pixel = out['pixel'].float().cpu().numpy()[0].squeeze(1)\n                \n                        if not np.isfinite(pixel).all():\n                            print(\"NON-FINITE from\", getattr(model, \"fusion_type\", \"unknown\"))\n                            print(\"min/max:\", np.nanmin(pixel), np.nanmax(pixel))\n                            raise ValueError(\"NaN/Inf detected in lead model output\")\n                \n                        if flip:\n                            pixel = np.flip(pixel, axis=flip)\n                \n                        for i in range(4):\n                            tu, tl, lu, ll = ens_regions[i]\n                            pixel_ens[i][tu:tl] += pixel[i][lu:ll]\n\n        # --- Weighted average ---\n        ens_weight = np.ones((trim_image.shape[0], trim_image.shape[1])) * len(whole_models) * len(tta)\n        for i in range(4):\n            tu, tl, _, _ = ens_regions[i]\n            ens_weight[tu:tl] += len(lead_models) * len(tta)\n        pixel_ens /= ens_weight\n\n        try:\n            # Subpixel pixel→series conversion\n            series_in_pixel = pixel_to_series_exp(pixel_ens[..., t0:t1], zero_mv_trimed, length)\n            series = (np.array(zero_mv_trimed).reshape(4, 1) - series_in_pixel) / mv_to_pixel\n            np.save(f'{OUT_DIR}/digitalised/{sample_id}.series.npy', series)\n        except:\n            local_fail_id.append(sample_id)\n\n        torch.cuda.empty_cache()\n\n       # if n<10 and gpu_id==0: # optional: show results\n         #   print()\n          #  print(\"check max intensity : \", np.max(pixel_ens))\n         #   print()\n            \n          #  overlay = draw_lead_pixel(trim_image, pixel_ens)\n          #  plt.imshow(overlay); plt.show()\n    \n           # t = np.arange(len(series[0]))\n          #  fig, axes = plt.subplots(4, 1, figsize=(12, 10))\n           # for j in range(4):\n             #   snr=0\n            #    axes[j].plot(t, series[j], alpha=1.0, color='blue', linewidth=1, label='predict')\n              #  axes[j].set_title(f'snr {snr:8.3f}')\n              #  axes[j].legend()\n         #   plt.show()\n    print('')\n   \n\n    print(f'\\n[GPU{gpu_id}] Stage2 done. Failed: {len(local_fail_id)}')\n    if fail_id_file:\n        with open(fail_id_file, 'wb') as f:\n            pickle.dump(local_fail_id, f)\n    return local_fail_id\n\nrun_stage2()\nprint('FAIL_ID:', FAIL_ID)\nprint('run_stage2() ok!!!\\n')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T19:10:07.324493Z","iopub.execute_input":"2026-05-02T19:10:07.324742Z","iopub.status.idle":"2026-05-02T19:10:07.341442Z","shell.execute_reply.started":"2026-05-02T19:10:07.324727Z","shell.execute_reply":"2026-05-02T19:10:07.340683Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"def split_series_to_leads(series_4xL, lead_lengths):\n    \"\"\"\n    Reconstruct 12 leads from your saved Stage2 array exactly like make_submission().\n    series_4xL shape: (4, L)\n    \"\"\"\n    series_by_lead = {}\n\n    lead_groups = [\n        ['I',   'aVR', 'V1', 'V4'],\n        ['II',  'aVL', 'V2', 'V5'],\n        ['III', 'aVF', 'V3', 'V6'],\n    ]\n\n    for row_idx, group in enumerate(lead_groups):\n        lengths = [int(lead_lengths[g]) for g in group]\n\n        if group[0] == 'II':\n            lengths[0] = lengths[0] - sum(lengths[1:])\n\n        cut_idx = np.cumsum(lengths)[:-1]\n        pieces = np.split(series_4xL[row_idx], cut_idx)\n\n        for lead_name, arr in zip(group, pieces):\n            series_by_lead[lead_name] = np.asarray(arr, dtype=np.float32)\n\n    series_by_lead['II'] = np.asarray(series_4xL[3], dtype=np.float32)\n    return series_by_lead\n\n\ndef score_one_local_sample(sample_id, out_dir=OUT_DIR, kaggle_dir=KAGGLE_DIR, train_meta_df=None):\n    image_id = sample_id.split('-')[0] if '-' in sample_id else sample_id\n\n    pred_path = f\"{out_dir}/digitalised/{sample_id}.series.npy\"\n    truth_path = f\"{kaggle_dir}/train/{image_id}/{image_id}.csv\"\n\n    if not os.path.exists(pred_path):\n        raise FileNotFoundError(f\"Missing prediction file: {pred_path}\")\n    if not os.path.exists(truth_path):\n        raise FileNotFoundError(f\"Missing truth file: {truth_path}\")\n\n    pred_series = np.load(pred_path)  # expected shape (4, L)\n    truth_df = pd.read_csv(truth_path)\n\n    if train_meta_df is None:\n        train_meta_df = pd.read_csv(f\"{kaggle_dir}/train.csv\")\n        train_meta_df[\"id\"] = train_meta_df[\"id\"].astype(str)\n\n    meta_row = train_meta_df[train_meta_df[\"id\"] == image_id].iloc[0]\n    fs = int(meta_row[\"fs\"])\n\n    lead_lengths = {lead: int(truth_df[lead].count()) for lead in LEADS}\n    pred_by_lead = split_series_to_leads(pred_series, lead_lengths)\n\n    sum_signal = 0.0\n    sum_noise = 0.0\n    lead_scores = {}\n\n    for lead in LEADS:\n        label = truth_df[lead].dropna().to_numpy(np.float64)\n        pred = np.asarray(pred_by_lead[lead], dtype=np.float64)\n\n        if len(pred) != len(label):\n            pred = scipy.signal.resample(pred, len(label)).astype(np.float64)\n\n        aligned_pred = align_signals(label, pred, int(fs * MAX_TIME_SHIFT))\n        p_signal, p_noise = compute_power(label, aligned_pred)\n\n        sum_signal += p_signal\n        sum_noise += p_noise\n        lead_scores[lead] = compute_snr(p_signal, p_noise)\n\n    image_snr_linear = compute_snr(sum_signal, sum_noise)\n    image_score_db = max(float(10 * np.log10(image_snr_linear)), -PERFECT_SCORE)\n\n    return {\n        \"sample_id\": sample_id,\n        \"image_id\": image_id,\n        \"fs\": fs,\n        \"snr_linear\": image_snr_linear,\n        \"score_db\": image_score_db,\n        **{f\"{lead}_snr\": lead_scores[lead] for lead in LEADS}\n    }\n\n\ndef score_local_predictions(sample_ids=None, out_dir=OUT_DIR, kaggle_dir=KAGGLE_DIR):\n    if sample_ids is None:\n        sample_ids = valid_id\n\n    train_meta_df = pd.read_csv(f\"{kaggle_dir}/train.csv\")\n    train_meta_df[\"id\"] = train_meta_df[\"id\"].astype(str)\n\n    rows = []\n    failed = []\n\n    for sample_id in sample_ids:\n        try:\n            rows.append(\n                score_one_local_sample(\n                    sample_id,\n                    out_dir=out_dir,\n                    kaggle_dir=kaggle_dir,\n                    train_meta_df=train_meta_df\n                )\n            )\n        except Exception as e:\n            failed.append((sample_id, str(e)))\n\n    score_df = pd.DataFrame(rows)\n\n    if len(score_df) > 0:\n        final_score_db = max(float(10 * np.log10(score_df[\"snr_linear\"].mean())), -PERFECT_SCORE)\n    else:\n        final_score_db = np.nan\n\n    print(f\"Local validation score = {final_score_db:.4f}\")\n    print(f\"Scored samples        = {len(score_df)}\")\n    print(f\"Failed samples        = {len(failed)}\")\n\n    if failed:\n        print(\"\\nFailures:\")\n        for x in failed[:10]:\n            print(x)\n\n    return score_df, final_score_db, failed","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2026-05-02T19:10:07.342022Z","iopub.status.idle":"2026-05-02T19:10:07.342335Z","shell.execute_reply.started":"2026-05-02T19:10:07.342177Z","shell.execute_reply":"2026-05-02T19:10:07.342192Z"}}},{"cell_type":"markdown","source":"score_df, local_score, failed = score_local_predictions()\nscore_df.sort_values(\"score_db\").head(10)\n\nEXP_NAME = \"notta_noattn_surfacefitting_1whl_1lead_model\"\nscore_df, local_score, failed = score_local_predictions()\n\nlog_row = pd.DataFrame([{\n    \"exp_name\": EXP_NAME,\n    \"local_score_db\": local_score,\n    \"n_samples\": len(score_df),\n    \"n_failed\": len(failed),\n}])\n\nlog_path = \"/kaggle/working/ablation_log.csv\"\nif os.path.exists(log_path):\n    old = pd.read_csv(log_path)\n    log_row = pd.concat([old, log_row], ignore_index=True)\n\nlog_row.to_csv(log_path, index=False)\nprint(log_row.tail())\nprint(EXP_NAME)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2026-05-02T19:10:07.343354Z","iopub.status.idle":"2026-05-02T19:10:07.343676Z","shell.execute_reply.started":"2026-05-02T19:10:07.343501Z","shell.execute_reply":"2026-05-02T19:10:07.343518Z"}}},{"cell_type":"code","source":"#make sbmission csv\n#FAIL_ID=[1053922973, ]\ndef make_submission():\n\tprint('===========================================')\n\tprint('making submission csv ...')\n\n\tsubmit_df=[]\n\tgb = valid_df.groupby('id')\n\tfor i,(sample_id, df) in enumerate(gb):\n        \n\t\t#if sample_id in FAIL_ID:\n\t\t#\tseries_by_lead = {}\n\t\t#\tfor j,d in df.iterrows():\n\t\t#\t\tseries_by_lead[d.lead] = np.zeros(d.number_of_rows)\n        \n\t\ttry:\n\t\t\tseries = np.load(f'{OUT_DIR}/digitalised/{sample_id}.series.npy')\n\n\t\t\tseries_by_lead={}\n\t\t\tfor l in range(3):\n\t\t\t\tlead = [\n\t\t\t\t\t['I',   'aVR', 'V1', 'V4'],\n\t\t\t\t\t['II',  'aVL', 'V2', 'V5'],\n\t\t\t\t\t['III', 'aVF', 'V3', 'V6'],\n\t\t\t\t][l]\n\n\t\t\t\tlength=[\n\t\t\t\t\tdf[df['lead']==lead[j]].iloc[0].number_of_rows\n\t\t\t\t\tfor j in range(4)\n\t\t\t\t]\n\t\t\t\tif lead[0]=='II':\n\t\t\t\t\tlength[0] = length[0]-sum(length[1:])\n\n\t\t\t\tindex = np.cumsum(length)[:-1]\n\t\t\t\tsplit = np.split(series[l], index)\n\t\t\t\t#print(length)\n\t\t\t\tfor (k, s) in zip(lead, split):\n\t\t\t\t\tseries_by_lead[k] = s\n\t\t\t\t\t#print(k,len(s))\n\t\t\tseries_by_lead['II'] = series[3]\n\t\t\t#print(series_by_lead)\n    \n\t\texcept: \n\t\t\tseries_by_lead = {}\n\t\t\tfor j,d in df.iterrows():\n\t\t\t\tseries_by_lead[d.lead] = np.zeros(d.number_of_rows)\n\n\t\t#print('\\r\\t {sample_id}', end='', flush=True)\n\t\tfor j,d in df.iterrows():\n\t\t\t#print(d.lead, len(series_by_lead[d.lead]),d.number_of_rows)\n            \n\t\t\t#assert(len(series_by_lead[d.lead])==d.number_of_rows)\n            \n            \n            #probably error here ... ???\n\t\t\tseries_by_lead[d.lead] = np.concatenate([\n                series_by_lead[d.lead], np.zeros_like(series_by_lead[d.lead])\n            ])[:d.number_of_rows]\n\t\t\tassert(len(series_by_lead[d.lead])==d.number_of_rows) \n\t\t\tprint(f'\\r\\t {i} {sample_id} : {d.lead}', end='', flush=True)\n\n\t\t\trow_id = [\n\t\t\t\tf'{sample_id}_{i}_{d.lead}' for i in range(d.number_of_rows)\n\t\t\t]\n\t\t\tthis_df = pd.DataFrame({\n\t\t\t\t'id':row_id,\n\t\t\t\t'value': series_by_lead[d.lead].astype(np.float32),\n\t\t\t})\n\t\t\tsubmit_df.append(this_df)\n\n\tprint('')\n\tsubmit_df = pd.concat(submit_df, axis=0, ignore_index=True, sort=False, copy=False)\n\tprint(submit_df)\n\tsubmit_df.to_csv('submission.csv',index=False)\n\nif (MODE=='fake')|(MODE=='submit'):\n    make_submission()\n    print('make_submission() ok!!!\\n')\n    if MODE=='submit':\n        shutil.rmtree(OUT_DIR)\n    !ls\n    #!rm -rf {OUT_DIR}\n\n'''\nfake:\n[21618231 rows x 2 columns]\n'''","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-05-02T19:10:07.344949Z","iopub.status.idle":"2026-05-02T19:10:07.345232Z","shell.execute_reply.started":"2026-05-02T19:10:07.345071Z","shell.execute_reply":"2026-05-02T19:10:07.345081Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"sub = pd.read_csv('submission.csv')\nprint(sub.head())\nprint(sub.shape)\nprint('NaN count:', sub[\"value\"].isna().sum())\nprint('Finite:', np.isfinite(sub[\"value\"].to_numpy()).all())","metadata":{"execution":{"iopub.status.busy":"2026-05-02T13:51:27.982539Z","iopub.execute_input":"2026-05-02T13:51:27.982842Z","iopub.status.idle":"2026-05-02T13:51:28.093791Z","shell.execute_reply.started":"2026-05-02T13:51:27.982816Z","shell.execute_reply":"2026-05-02T13:51:28.092643Z"},"trusted":true},"outputs":[],"execution_count":null}]}