Faster gamma calculations: a reproducible benchmark

How much time does the gamma implementation in PR #2066 save, and does it produce the same answer? This notebook compares two exact source revisions on the same computer, using the six synthetic workloads from the original DICOM coordinate demonstration. It measures complete, warmed-up gamma calls and checks every output value.

For physicists: start with the workloads and results below. Expand the code when you want to inspect or repeat the experiment. The companion DICOM coordinates and gamma, illustrated explains coordinate correctness; the ascending, regular grids here isolate performance from those geometry fixes.

The default run analyses the recorded measurements embedded in this notebook. It regenerates the figures without rerunning the long benchmark. Section 7 explains how to collect fresh measurements; configured runs must pass the source and numerical checks before their timings are plotted.

This is the earlier fixed-grid experiment. For the comprehensive four-way audit with calibrated scaling cases and specified volume examples, use the background runner and upload guide. Its continuous phantom differs from the rasterised phantom below; the two timing datasets must not be pooled. The recorded workstation evidence retains verified partial scaling measurements; the background runner now has a two-hour total limit, paired PyMedPhys/SciPy speed ratios and synthetic SABR/pelvic volume examples.

Recorded result: the three 3D workloads using the default interpolator saved 21.2–25.2% of median runtime on the measured Windows host, with exact elementwise agreement in all six workloads. This is a workload-specific saving, not a universal speed guarantee. Absolute runtimes, all individual measurements and paired-round variation are provided below.

Show code cell source

Hide code cell source

import hashlib
import json
import os
import platform
import re
import subprocess
import sys
import tempfile
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np
from html import escape
from IPython.display import HTML, Markdown, display

plt.rcParams.update({
    "figure.figsize": (8, 4.5), "figure.dpi": 120,
    "font.size": 11, "axes.spines.top": False,
    "axes.spines.right": False, "axes.axisbelow": True,
})
COLOURS = {"previous": "#B76121", "current": "#176B91"}
LABELS = {
    "3d": "3D · 3% / 3 mm",
    "3d-maxgamma": "3D · 3% / 3 mm · cap 2",
    "3d-local": "3D · 2% / 2 mm · local",
    "3d-scipy": "3D · 3% / 3 mm · SciPy",
    "2d": "2D · 3% / 3 mm",
    "2d-scipy": "2D · 3% / 3 mm · SciPy",
}


def table(headings, rows):
    """Render tables in both Jupyter and Sphinx's CommonMark output renderer."""
    def content(value):
        text = str(value)
        link = re.fullmatch(
            r"\[`([0-9a-f]{40})`\]\((https://github.com/pymedphys/pymedphys/tree/[0-9a-f]{40})\)",
            text,
        )
        if link:
            return f'<a href="{escape(link[2], quote=True)}"><code style="overflow-wrap:anywhere">{link[1]}</code></a>'
        return escape(text)

    heading = "".join(f'<th style="text-align:left;padding:0.4em">{escape(h)}</th>' for h in headings)
    body = "".join("<tr>" + "".join(f'<td style="padding:0.4em">{content(v)}</td>' for v in row) + "</tr>" for row in rows)
    display(HTML(f'<div style="overflow-x:auto"><table><thead><tr>{heading}</tr></thead><tbody>{body}</tbody></table></div>'))

Show code cell source

Hide code cell source

RECORDED = json.loads(r'''
{
  "times": {
    "3d": {
      "previous": [
        5.147435200000473,
        5.099064099998941,
        4.956633500001772,
        5.019999200001621,
        4.988596399998642,
        4.568466200002149
      ],
      "current": [
        4.030762700000196,
        3.986416600000666,
        3.7935049000006984,
        3.695042999999714,
        3.501579600000696,
        3.637928299998748
      ]
    },
    "3d-maxgamma": {
      "previous": [
        5.011578499997995,
        4.865375299999869,
        5.493412799998623,
        5.043869299999642,
        4.742424399999436,
        4.9588997999999265
      ],
      "current": [
        3.843475900001067,
        3.735443800000212,
        3.767609499998798,
        3.7530423999996856,
        3.7305885999994643,
        3.74873419999858
      ]
    },
    "3d-local": {
      "previous": [
        23.260586299998977,
        20.365133200000855,
        19.378479200000584,
        18.720706200001587,
        19.28715580000062,
        19.852144299999054
      ],
      "current": [
        15.385054399997898,
        22.84703680000166,
        15.764805399998295,
        15.54000200000155,
        14.848236699999688,
        14.8664858000011
      ]
    },
    "3d-scipy": {
      "previous": [
        15.590848099996947,
        18.868159800000285,
        14.566667700000835,
        15.356425999998464,
        16.16232559999844,
        15.786817000000156
      ],
      "current": [
        18.107391499997902,
        18.187399899998127,
        15.21309409999958,
        15.043153099999472,
        15.398140299999795,
        15.933863199999905
      ]
    },
    "2d": {
      "previous": [
        0.17002030000003288,
        0.16738600000098813,
        0.15945759999885922,
        0.16313910000098986,
        0.16855959999884362,
        0.16808160000073258
      ],
      "current": [
        0.14230329999918467,
        0.13905960000192863,
        0.1385162999977183,
        0.138943099998869,
        0.1437337999996089,
        0.14524040000105742
      ]
    },
    "2d-scipy": {
      "previous": [
        0.19588399999702233,
        0.20161559999905876,
        0.1873674000016763,
        0.18888019999940298,
        0.2020042000003741,
        0.2000906000030227
      ],
      "current": [
        0.2045572999995784,
        0.20080570000209264,
        0.18861909999759519,
        0.19107760000042617,
        0.19596799999999348,
        0.1928646000014851
      ]
    }
  },
  "records": [
    {
      "case": "3d",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 7.237295200000517,
      "times": [
        5.147435200000473,
        5.099064099998941
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 0,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 5.1677712999990035,
      "times": [
        4.030762700000196,
        3.986416600000666
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 0,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-maxgamma",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "max_gamma": 2
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 5.117816500001936,
      "times": [
        5.011578499997995,
        4.865375299999869
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 0,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-maxgamma",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "max_gamma": 2
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 3.7880229000002146,
      "times": [
        3.843475900001067,
        3.735443800000212
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 0,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-local",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 2,
        "distance_mm_threshold": 2,
        "lower_percent_dose_cutoff": 10,
        "local_gamma": true
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 46411,
      "max_gamma": 2.100060400755779,
      "warmup_seconds": 21.675853100001405,
      "times": [
        23.260586299998977,
        20.365133200000855
      ],
      "gamma_sum": 41758.367141341085,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 0,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "025b774c5f573195901404082633f44503ba71671696a701f3be36221d5dc0fe",
      "elementwise_equal": true
    },
    {
      "case": "3d-local",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 2,
        "distance_mm_threshold": 2,
        "lower_percent_dose_cutoff": 10,
        "local_gamma": true
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 46411,
      "max_gamma": 2.100060400755779,
      "warmup_seconds": 15.543868900000234,
      "times": [
        15.385054399997898,
        22.84703680000166
      ],
      "gamma_sum": 41758.367141341085,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 0,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "025b774c5f573195901404082633f44503ba71671696a701f3be36221d5dc0fe",
      "elementwise_equal": true
    },
    {
      "case": "3d-scipy",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 21.80124810000052,
      "times": [
        15.590848099996947,
        18.868159800000285
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 0,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "f553b5c0729c654aa53439f8a6fd655357d273e38e2eb0de46e9798e41d4ece9",
      "elementwise_equal": true
    },
    {
      "case": "3d-scipy",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 18.252869399999327,
      "times": [
        18.107391499997902,
        18.187399899998127
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 0,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "f553b5c0729c654aa53439f8a6fd655357d273e38e2eb0de46e9798e41d4ece9",
      "elementwise_equal": true
    },
    {
      "case": "2d",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 1.3336466000000655,
      "times": [
        0.17002030000003288,
        0.16738600000098813
      ],
      "gamma_sum": 14506.752795157714,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 0,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "e255c80d91fe82e96df01074f23a2eb8d4cee7593db9e12329269e50c6225ed9",
      "elementwise_equal": true
    },
    {
      "case": "2d",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 1.286129100000835,
      "times": [
        0.14230329999918467,
        0.13905960000192863
      ],
      "gamma_sum": 14506.752795157714,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 0,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "e255c80d91fe82e96df01074f23a2eb8d4cee7593db9e12329269e50c6225ed9",
      "elementwise_equal": true
    },
    {
      "case": "2d-scipy",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.37044710000191117,
      "times": [
        0.19588399999702233,
        0.20161559999905876
      ],
      "gamma_sum": 14506.752795157712,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 0,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "468498d5d66aa71e586c161add81a3d28e8b9395823b68228476e0d896867b64",
      "elementwise_equal": true
    },
    {
      "case": "2d-scipy",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.37463669999851845,
      "times": [
        0.2045572999995784,
        0.20080570000209264
      ],
      "gamma_sum": 14506.752795157712,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 0,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "468498d5d66aa71e586c161add81a3d28e8b9395823b68228476e0d896867b64",
      "elementwise_equal": true
    },
    {
      "case": "3d",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 4.10151229999974,
      "times": [
        3.7935049000006984,
        3.695042999999714
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 1,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 5.167867300002399,
      "times": [
        4.956633500001772,
        5.019999200001621
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 1,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-maxgamma",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "max_gamma": 2
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 3.9201573000027565,
      "times": [
        3.767609499998798,
        3.7530423999996856
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 1,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-maxgamma",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "max_gamma": 2
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 5.102866900000663,
      "times": [
        5.493412799998623,
        5.043869299999642
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 1,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-local",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 2,
        "distance_mm_threshold": 2,
        "lower_percent_dose_cutoff": 10,
        "local_gamma": true
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 46411,
      "max_gamma": 2.100060400755779,
      "warmup_seconds": 14.907378799998696,
      "times": [
        15.764805399998295,
        15.54000200000155
      ],
      "gamma_sum": 41758.367141341085,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 1,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "025b774c5f573195901404082633f44503ba71671696a701f3be36221d5dc0fe",
      "elementwise_equal": true
    },
    {
      "case": "3d-local",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 2,
        "distance_mm_threshold": 2,
        "lower_percent_dose_cutoff": 10,
        "local_gamma": true
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 46411,
      "max_gamma": 2.100060400755779,
      "warmup_seconds": 19.71661540000059,
      "times": [
        19.378479200000584,
        18.720706200001587
      ],
      "gamma_sum": 41758.367141341085,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 1,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "025b774c5f573195901404082633f44503ba71671696a701f3be36221d5dc0fe",
      "elementwise_equal": true
    },
    {
      "case": "3d-scipy",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 15.857877699996607,
      "times": [
        15.21309409999958,
        15.043153099999472
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 1,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "f553b5c0729c654aa53439f8a6fd655357d273e38e2eb0de46e9798e41d4ece9",
      "elementwise_equal": true
    },
    {
      "case": "3d-scipy",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 15.416765300000407,
      "times": [
        14.566667700000835,
        15.356425999998464
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 1,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "f553b5c0729c654aa53439f8a6fd655357d273e38e2eb0de46e9798e41d4ece9",
      "elementwise_equal": true
    },
    {
      "case": "2d",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.29505019999851356,
      "times": [
        0.1385162999977183,
        0.138943099998869
      ],
      "gamma_sum": 14506.752795157714,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 1,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "e255c80d91fe82e96df01074f23a2eb8d4cee7593db9e12329269e50c6225ed9",
      "elementwise_equal": true
    },
    {
      "case": "2d",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.3361216000012064,
      "times": [
        0.15945759999885922,
        0.16313910000098986
      ],
      "gamma_sum": 14506.752795157714,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 1,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "e255c80d91fe82e96df01074f23a2eb8d4cee7593db9e12329269e50c6225ed9",
      "elementwise_equal": true
    },
    {
      "case": "2d-scipy",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.35549689999970724,
      "times": [
        0.18861909999759519,
        0.19107760000042617
      ],
      "gamma_sum": 14506.752795157712,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 1,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "468498d5d66aa71e586c161add81a3d28e8b9395823b68228476e0d896867b64",
      "elementwise_equal": true
    },
    {
      "case": "2d-scipy",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.354619099998672,
      "times": [
        0.1873674000016763,
        0.18888019999940298
      ],
      "gamma_sum": 14506.752795157712,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 1,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "468498d5d66aa71e586c161add81a3d28e8b9395823b68228476e0d896867b64",
      "elementwise_equal": true
    },
    {
      "case": "3d",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 5.090694800001074,
      "times": [
        4.988596399998642,
        4.568466200002149
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 2,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 3.8940230999978667,
      "times": [
        3.501579600000696,
        3.637928299998748
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 2,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-maxgamma",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "max_gamma": 2
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 5.05139699999927,
      "times": [
        4.742424399999436,
        4.9588997999999265
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 2,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-maxgamma",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "max_gamma": 2
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 3.9362825999996858,
      "times": [
        3.7305885999994643,
        3.74873419999858
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 2,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "6654167b498f2d92581ed261599998958f37ec64c0d56094eeaf662f35767000",
      "elementwise_equal": true
    },
    {
      "case": "3d-local",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 2,
        "distance_mm_threshold": 2,
        "lower_percent_dose_cutoff": 10,
        "local_gamma": true
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 46411,
      "max_gamma": 2.100060400755779,
      "warmup_seconds": 20.23296140000093,
      "times": [
        19.28715580000062,
        19.852144299999054
      ],
      "gamma_sum": 41758.367141341085,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 2,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "025b774c5f573195901404082633f44503ba71671696a701f3be36221d5dc0fe",
      "elementwise_equal": true
    },
    {
      "case": "3d-local",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 2,
        "distance_mm_threshold": 2,
        "lower_percent_dose_cutoff": 10,
        "local_gamma": true
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 46411,
      "max_gamma": 2.100060400755779,
      "warmup_seconds": 15.265261000000464,
      "times": [
        14.848236699999688,
        14.8664858000011
      ],
      "gamma_sum": 41758.367141341085,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 2,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "025b774c5f573195901404082633f44503ba71671696a701f3be36221d5dc0fe",
      "elementwise_equal": true
    },
    {
      "case": "3d-scipy",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 16.490664699998888,
      "times": [
        16.16232559999844,
        15.786817000000156
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 2,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "f553b5c0729c654aa53439f8a6fd655357d273e38e2eb0de46e9798e41d4ece9",
      "elementwise_equal": true
    },
    {
      "case": "3d-scipy",
      "shape": [
        41,
        57,
        57
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "e49faed5c475154ee2b0a66a5bf685883210051617064855158655bfa58bef3f",
      "eligible_points": 63935,
      "finite_points": 63935,
      "pass_count": 62246,
      "max_gamma": 1.309989450836535,
      "warmup_seconds": 16.610222299997986,
      "times": [
        15.398140299999795,
        15.933863199999905
      ],
      "gamma_sum": 24884.948642116633,
      "nan_count": 69274,
      "points": 133209,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 2,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "f553b5c0729c654aa53439f8a6fd655357d273e38e2eb0de46e9798e41d4ece9",
      "elementwise_equal": true
    },
    {
      "case": "2d",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.33281129999886616,
      "times": [
        0.16855959999884362,
        0.16808160000073258
      ],
      "gamma_sum": 14506.752795157714,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 2,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "e255c80d91fe82e96df01074f23a2eb8d4cee7593db9e12329269e50c6225ed9",
      "elementwise_equal": true
    },
    {
      "case": "2d",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.30494479999833857,
      "times": [
        0.1437337999996089,
        0.14524040000105742
      ],
      "gamma_sum": 14506.752795157714,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 2,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "e255c80d91fe82e96df01074f23a2eb8d4cee7593db9e12329269e50c6225ed9",
      "elementwise_equal": true
    },
    {
      "case": "2d-scipy",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.3653806000002078,
      "times": [
        0.2020042000003741,
        0.2000906000030227
      ],
      "gamma_sum": 14506.752795157712,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "previous",
      "round": 2,
      "revision": "866f83edad8586a42a739094e488f45242b72c95",
      "gamma_sha256": "468498d5d66aa71e586c161add81a3d28e8b9395823b68228476e0d896867b64",
      "elementwise_equal": true
    },
    {
      "case": "2d-scipy",
      "shape": [
        401,
        401
      ],
      "options": {
        "dose_percent_threshold": 3,
        "distance_mm_threshold": 3,
        "lower_percent_dose_cutoff": 10,
        "interp_algo": "scipy"
      },
      "input_sha256": "ef1abac16a40cd5184d639e6adb5f32ce35f414089241fd188a255e6096d498e",
      "eligible_points": 52379,
      "finite_points": 52379,
      "pass_count": 52379,
      "max_gamma": 0.7631110418832897,
      "warmup_seconds": 0.3797915000031935,
      "times": [
        0.19596799999999348,
        0.1928646000014851
      ],
      "gamma_sum": 14506.752795157712,
      "nan_count": 108422,
      "points": 160801,
      "origins": {
        "pymedphys": "lib/pymedphys/__init__.py",
        "pymedphys._gamma.implementation.shell": "lib/pymedphys/_gamma/implementation/shell.py",
        "pymedphys._interp.interp": "lib/pymedphys/_interp/interp.py"
      },
      "python": "3.12.14 (main, Sep  1 2026, 14:17:39) [MSC v.1944 64 bit (AMD64)]",
      "platform": "Windows-11-10.0.26200-SP0",
      "processor": "Intel64 Family 6 Model 198 Stepping 2, GenuineIntel",
      "versions": {
        "numpy": "1.26.4",
        "scipy": "1.16.2",
        "numba": "0.61.2"
      },
      "numba_threads": 2,
      "version": "current",
      "round": 2,
      "revision": "d99893b46e3adb34e384b100d9932d87f5e38804",
      "gamma_sha256": "468498d5d66aa71e586c161add81a3d28e8b9395823b68228476e0d896867b64",
      "elementwise_equal": true
    }
  ],
  "elementwise_equal": true,
  "worker_sha256": "2d71d54740f11e85b82a1f392241bc6d0770b5930864f6320e86cfbcbb7b16db",
  "rounds": 3,
  "repeats": 2,
  "threads": 2,
  "measured_at_utc": "2026-09-26T14:55:14.446058+00:00",
  "hardware": "Intel Core Ultra 7 265HX; 20 cores / 20 logical processors; 2 Numba threads"
}
''')

1. What is held constant?

Gamma combines dose difference and spatial distance. With global 3% / 3 mm gamma, the dose denominator is 3% of the maximum reference dose: here 0.06 Gy because the reference maximum is 2 Gy. With local 2% / 2 mm gamma, it is 2% of the reference dose at each tested point, making the low-dose part of this workload more demanding. These cases change both normalisation and thresholds; their runtime difference does not isolate either factor.

Both versions use identical dose arrays, coordinates and settings. All reference points at or above 10% of the reference maximum are evaluated. There is no random subset and no early exit merely because a point passes. The interpolation fraction is 10; a 3 mm distance criterion therefore starts with a 0.3 mm radial search step. The RAM chunk budget is the unchanged default of 1.5 GiB, which is not a limit on total process memory.

The capped case uses max_gamma=2; the others leave it unset. This cap bounds the search and caps returned finite values at 2. It must be applied identically to both versions. We do not obtain the speed-up by loosening the gamma criterion or reducing sampling accuracy.

The synthetic dose model

Two rectangular fields and a smaller boost are sampled on the grid, then smoothed with a Gaussian of 6 mm standard deviation. The reference is normalised to a maximum of 2 Gy. The evaluation uses nominal box-centre offsets of (1.0, −0.7, 0.5) mm in 3D, or (1.0, −0.7) mm in 2D, and a maximum of 2.04 Gy. The boxes are rasterised before smoothing: those offsets change discrete membership at field edges; they are not an exact translation of a continuous dose distribution. This construction preserves the original benchmark workloads.

The 3D axes are (z, y, x), covering ±50, ±70 and ±70 mm at 2.5 mm spacing: 41 × 57 × 57 points. The 2D axes are (y, x), each covering ±100 mm at 0.5 mm spacing: 401 × 401 points. Both grids are ascending, evenly spaced and coincident. Synthetic Gy values provide scale; these are not patient plans or measurements.

2. Reproducible inputs and measurements

The following folded cells contain the complete worker and controller. Each timed worker imports one specified checkout in a fresh process and verifies the origins of PyMedPhys, its gamma implementation and its interpolator. The controller rejects missing paths, shortened or mismatched commit IDs, identical checkouts/revisions, and modified Python source or dependency files. It checks the source again after the run.

There are three rounds per case. Within each round, one process per revision performs an untimed full gamma warm-up followed by two timed calls. Revision order alternates between rounds: previous/current, current/previous, previous/current. Input generation, imports, compilation or cache loading, result comparison and file writing are excluded; preparation performed inside the gamma call is included. This measures repeated use after warm-up, not the wait for a first-ever call.

After every timed call, its complete array must equal the warm-up array. The controller also compares every element, shape and NaN position across revisions and rounds. It saves all individual timings and records the input and output hashes, source revisions, import locations, dependency versions and Numba thread count. Checksums make the records auditable; they do not replace the actual elementwise comparison.

Show code cell source

Hide code cell source

WORKER = r'''
"""Time pymedphys.gamma on fixed synthetic cases; run once per source tree.

Usage: python bench_worker.py CASE REPEATS CHECKOUT OUTPUT_NPY
Checks module origins, saves the full gamma array and prints timing/provenance JSON.
"""

import json
import hashlib
import importlib
import platform
from pathlib import Path
import numba
import scipy
import sys
import time
import warnings

import numpy as np
import scipy.ndimage

import pymedphys


def field(axes, centre_shift=(0.0, 0.0, 0.0), scale=1.0):
    """Smooth synthetic dose: two crossing fields and a boost."""
    grids = np.meshgrid(*axes, indexing="ij")
    ndim = len(axes)
    dose = np.zeros(grids[0].shape)

    def box(half_widths, centre, weight):
        inside = np.ones(grids[0].shape, dtype=bool)
        for g, h, c, s in zip(grids, half_widths, centre, centre_shift):
            inside &= np.abs(g - c - s) <= h
        return weight * inside

    centre = [0.0] * ndim
    dose += box([40, 60, 25][:ndim], centre, 1.0)
    dose += box([55, 30, 45][:ndim], centre, 0.8)
    dose += box([15, 15, 15][:ndim], [5, -8, 4][:ndim], 0.6)
    spacing = [a[1] - a[0] for a in axes]
    sigma = [6.0 / s for s in spacing]
    dose = scipy.ndimage.gaussian_filter(dose, sigma, mode="constant")
    return 2.0 * scale * dose / dose.max()


def case_inputs(case):
    if case.startswith("3d"):
        ref_axes = tuple(np.arange(-n, n + 1e-9, 2.5) for n in (50.0, 70.0, 70.0))
        eval_axes = ref_axes
    elif case.startswith("2d"):
        ref_axes = tuple(np.arange(-n, n + 1e-9, 0.5) for n in (100.0, 100.0))
        eval_axes = ref_axes
    else:
        raise ValueError(case)
    ndim = len(ref_axes)
    ref = field(ref_axes)
    ev = field(eval_axes, centre_shift=(1.0, -0.7, 0.5)[:ndim], scale=1.02)
    options = dict(
        dose_percent_threshold=3,
        distance_mm_threshold=3,
        lower_percent_dose_cutoff=10,
    )
    if "scipy" in case:
        options["interp_algo"] = "scipy"
    if "local" in case:
        options.update(dose_percent_threshold=2, distance_mm_threshold=2, local_gamma=True)
    if "maxgamma" in case:
        options["max_gamma"] = 2
    return ref_axes, ref, eval_axes, ev, options


def main():
    case, repeats = sys.argv[1], int(sys.argv[2])
    ref_axes, ref, eval_axes, ev, options = case_inputs(case)
    root = Path(sys.argv[3]).resolve()
    origins = {}
    for name in ("pymedphys", "pymedphys._gamma.implementation.shell", "pymedphys._interp.interp"):
        module = importlib.import_module(name)
        origin = Path(module.__file__).resolve()
        if not origin.is_relative_to(root / "lib" / "pymedphys"):
            raise RuntimeError(f"Wrong imported module: {name} from {origin}; expected {root}")
        origins[name] = str(origin)
    warnings.simplefilter("error")
    # Warm-up: loads cached Numba kernels and SciPy code paths.
    warmup_start = time.perf_counter()
    gamma = pymedphys.gamma(ref_axes, ref, eval_axes, ev, **options)
    warmup_seconds = time.perf_counter() - warmup_start
    baseline = gamma.copy()
    times = []
    for _ in range(repeats):
        start = time.perf_counter()
        gamma = pymedphys.gamma(ref_axes, ref, eval_axes, ev, **options)
        times.append(time.perf_counter() - start)
        np.testing.assert_array_equal(gamma, baseline)
    np.save(sys.argv[4], gamma, allow_pickle=False)
    print(
        json.dumps(
            {
                "case": case,
                "shape": list(ref.shape),
                "options": options,
                "input_sha256": hashlib.sha256(b"".join(
                    array.tobytes() for array in (*ref_axes, ref, *eval_axes, ev)
                )).hexdigest(),
                "eligible_points": int(np.count_nonzero(ref >= 0.1 * ref.max())),
                "finite_points": int(np.isfinite(gamma).sum()),
                "pass_count": int(np.count_nonzero(gamma <= 1)),
                "max_gamma": float(np.nanmax(gamma)),
                "warmup_seconds": warmup_seconds,
                "times": times,
                "gamma_sum": float(np.nansum(gamma)),
                "nan_count": int(np.isnan(gamma).sum()),
                "points": int(np.size(gamma)),
                "origins": origins,
                "python": sys.version,
                "platform": platform.platform(),
                "processor": platform.processor(),
                "versions": {"numpy": np.__version__, "scipy": scipy.__version__, "numba": numba.__version__},
                "numba_threads": numba.get_num_threads(),
            }
        )
    )


if __name__ == "__main__":
    main()
'''
worker_functions = {'__name__': 'benchmark_worker'}
exec(WORKER, worker_functions)

Show code cell source

Hide code cell source

CASES = ["3d", "3d-maxgamma", "3d-local", "3d-scipy", "2d", "2d-scipy"]


def checked_checkout(path, revision):
    root = Path(path).resolve(strict=True)
    if not (root / "lib" / "pymedphys" / "__init__.py").is_file():
        raise ValueError(f"No PyMedPhys source tree at {root}")
    if not re.fullmatch(r"[0-9a-f]{40}", revision or ""):
        raise ValueError("Specify the exact 40-character commit ID, not a branch name")

    def git(*args):
        return subprocess.check_output(
            ["git", "-C", str(root), *args], text=True
        ).strip()

    if git("rev-parse", "HEAD") != revision:
        raise ValueError(f"{root} is not at the requested revision {revision}")
    # Notebook outputs may be dirty; imported source and dependencies may not.
    changes = git(
        "status",
        "--porcelain",
        "--",
        ":(glob)lib/pymedphys/**/*.py",
        "pyproject.toml",
        "uv.lock",
    )
    if changes:
        raise ValueError(f"Benchmark source/dependency files are modified: {changes}")
    return root


def run_benchmark(
    previous_path,
    previous_revision,
    current_path,
    current_revision,
    cases=CASES,
    rounds=3,
    repeats=2,
    threads=2,
    progress_path=None,
):
    roots = {
        "previous": checked_checkout(previous_path, previous_revision),
        "current": checked_checkout(current_path, current_revision),
    }
    if roots["previous"] == roots["current"] or previous_revision == current_revision:
        raise ValueError("The benchmark requires two distinct checkouts and revisions")
    if not cases or any(case not in CASES for case in cases):
        raise ValueError("Select at least one supported benchmark case")
    if min(rounds, repeats, threads) < 1:
        raise ValueError("Rounds, repeats and threads must be positive")
    records, arrays, timings = (
        [],
        {},
        {case: {version: [] for version in roots} for case in cases},
    )
    with tempfile.TemporaryDirectory() as directory:
        directory = Path(directory)
        worker_path = directory / "worker.py"
        worker_path.write_text(WORKER, encoding="utf-8")
        for round_index in range(rounds):
            for case in cases:
                order = (
                    ["previous", "current"]
                    if round_index % 2 == 0
                    else ["current", "previous"]
                )
                for version in order:
                    array_path = directory / f"{case}-{version}-{round_index}.npy"
                    root = roots[version]
                    print(f"Round {round_index + 1}/{rounds}: {case}, {version}", flush=True)
                    result = subprocess.run(
                        [
                            sys.executable,
                            str(worker_path),
                            case,
                            str(repeats),
                            str(root),
                            str(array_path),
                        ],
                        env=dict(
                            os.environ,
                            PYTHONPATH=str(root / "lib"),
                            NUMBA_NUM_THREADS=str(threads),
                            PYTHONDONTWRITEBYTECODE="1",
                        ),
                        cwd=directory,
                        capture_output=True,
                        text=True,
                        check=True,
                        timeout=600,
                    )
                    record = json.loads(result.stdout)
                    print(f"  {record['times']} s", flush=True)
                    for origin in record["origins"].values():
                        if (
                            not Path(origin)
                            .resolve()
                            .is_relative_to(root / "lib" / "pymedphys")
                        ):
                            raise RuntimeError(
                                f"Unexpected import provenance: {origin}"
                            )
                    gamma_values = np.load(array_path, allow_pickle=False)
                    if case in arrays:
                        # Compare every element and NaN position, across versions and rounds.
                        np.testing.assert_array_equal(gamma_values, arrays[case])
                    else:
                        arrays[case] = gamma_values
                    record.update(
                        version=version,
                        round=round_index,
                        revision=previous_revision
                        if version == "previous"
                        else current_revision,
                    )
                    record["origins"] = {
                        name: Path(origin).relative_to(root).as_posix()
                        for name, origin in record["origins"].items()
                    }
                    record["gamma_sha256"] = hashlib.sha256(gamma_values.tobytes()).hexdigest()
                    record["elementwise_equal"] = True
                    records.append(record)
                    if progress_path is not None:
                        Path(progress_path).write_text(json.dumps(records, indent=2), encoding="utf-8")
                    timings[case][version].extend(record["times"])
    for version, root in roots.items():
        checked_checkout(root, previous_revision if version == "previous" else current_revision)
    return {
        "times": timings,
        "records": records,
        "elementwise_equal": True,
        "worker_sha256": hashlib.sha256(WORKER.encode()).hexdigest(),
        "rounds": rounds,
        "repeats": repeats,
        "threads": threads,
    }

Show code cell source

Hide code cell source

previous_path = os.environ.get("PYMEDPHYS_BENCH_PREVIOUS")
current_path = os.environ.get("PYMEDPHYS_BENCH_CURRENT")
if previous_path is not None or current_path is not None:
    if not previous_path or not current_path:
        raise ValueError("Configure both benchmark checkout paths")
    result = run_benchmark(
        previous_path, os.environ.get("PYMEDPHYS_BENCH_PREVIOUS_SHA"),
        current_path, os.environ.get("PYMEDPHYS_BENCH_CURRENT_SHA"),
        threads=int(os.environ.get("PYMEDPHYS_BENCH_THREADS", "2")),
    )
    source_label = "Fresh measurements collected by this notebook run"
else:
    result = RECORDED
    source_label = "Recorded measurements embedded in this notebook"

assert result["elementwise_equal"]
assert result["worker_sha256"] == hashlib.sha256(WORKER.encode()).hexdigest()
records = result["records"]
assert len(records) == len(CASES) * 2 * result["rounds"]
for case in CASES:
    selected = [r for r in records if r["case"] == case]
    assert len({r["input_sha256"] for r in selected}) == 1
    assert len({r["gamma_sha256"] for r in selected}) == 1
    assert all(r["elementwise_equal"] for r in selected)
    assert all(r["finite_points"] == r["eligible_points"] for r in selected)
    for version in ("previous", "current"):
        by_version = [r for r in selected if r["version"] == version]
        assert sorted(r["round"] for r in by_version) == list(range(result["rounds"]))
        values = [value for r in by_version for value in r["times"]]
        np.testing.assert_array_equal(values, result["times"][case][version])
        assert len(values) == result["rounds"] * result["repeats"]
        assert np.isfinite(values).all() and np.all(np.array(values) > 0)
assert len({json.dumps(r["versions"], sort_keys=True) for r in records}) == 1
assert all(r["numba_threads"] == result["threads"] for r in records)
assert all(len({r["revision"] for r in records if r["version"] == v}) == 1
           for v in ("previous", "current"))
assert records[0]["revision"] != next(r["revision"] for r in records if r["version"] == "current")
display(Markdown(f"**Data used below: {source_label}.**"))

medians = {case: {v: float(np.median(result["times"][case][v]))
                  for v in ("previous", "current")} for case in CASES}

Data used below: Recorded measurements embedded in this notebook.

Show code cell source

Hide code cell source

first = records[0]
revisions = {v: next(r["revision"] for r in records if r["version"] == v)
             for v in ("previous", "current")}
table(["Provenance", "Value"], [
    ["Previous source", f"[`{revisions['previous']}`](https://github.com/pymedphys/pymedphys/tree/{revisions['previous']})"],
    ["Current source", f"[`{revisions['current']}`](https://github.com/pymedphys/pymedphys/tree/{revisions['current']})"],
    ["Recorded run completed (UTC)", result.get("measured_at_utc", "Fresh run; see this execution")],
    ["Processor", result.get("hardware", first["processor"])],
    ["Platform", first["platform"]],
    ["Python", first["python"].split()[0]],
    ["NumPy / SciPy / Numba", " / ".join(first["versions"][v] for v in ("numpy", "scipy", "numba"))],
    ["Numba threads", result["threads"]],
    ["Measurements", f"{result['rounds']} rounds × {result['repeats']} timed calls per revision per case"],
])

Provenance	Value
Previous source	866f83edad8586a42a739094e488f45242b72c95
Current source	d99893b46e3adb34e384b100d9932d87f5e38804
Recorded run completed (UTC)	2026-09-26T14:55:14.446058+00:00
Processor	Intel Core Ultra 7 265HX; 20 cores / 20 logical processors; 2 Numba threads
Platform	Windows-11-10.0.26200-SP0
Python	3.12.14
NumPy / SciPy / Numba	1.26.4 / 1.16.2 / 0.61.2
Numba threads	2
Measurements	3 rounds × 2 timed calls per revision per case

The recorded comparison is previous main 866f83e versus PR revision d99893b. Both run in the same Python environment, so this tests the source changes rather than a simultaneous dependency upgrade. Numba’s thread setting is not a limit on every thread in the process. The recorded host is a Windows laptop with a hybrid-core processor; core placement, power management and other activity can affect timings.

The next picture regenerates the workload from the embedded worker, outside any timed region. It shows the central z slice of the 3D reference and the evaluation-minus-reference difference. Colours in the difference panel have a common zero and symmetric limits.

Show code cell source

Hide code cell source

axes, reference, _, evaluation, _ = worker_functions["case_inputs"]("3d")
z, y, x = axes
plane = len(z) // 2
extent = [x[0] - 1.25, x[-1] + 1.25, y[0] - 1.25, y[-1] + 1.25]
difference = evaluation[plane] - reference[plane]
limit = np.max(np.abs(difference))
fig, axs = plt.subplots(1, 2, figsize=(8.4, 4.1), layout="constrained")
image = axs[0].imshow(reference[plane], extent=extent, origin="lower", vmin=0, vmax=2, cmap="viridis")
fig.colorbar(image, ax=axs[0], label="Reference dose (Gy)", shrink=0.8)
image = axs[1].imshow(difference, extent=extent, origin="lower", vmin=-limit, vmax=limit, cmap="RdBu_r")
fig.colorbar(image, ax=axs[1], label="Evaluation − reference (Gy)", shrink=0.8)
for ax, title in zip(axs, ["Reference, z = 0 mm", "Difference, z = 0 mm"]):
    ax.set(title=title, xlabel="x (mm)", ylabel="y (mm)", aspect="equal")
plt.show()

3. Results: complete gamma-call runtime

The table reports the median of all six timed calls per revision. Time saved is 100 × (1 − current / previous); speed-up is previous / current. For example, reducing a 10 s calculation to 6 s saves 40% of the time and is a 1.67× speed-up. A negative time saving means the current version took longer.

Show code cell source

Hide code cell source

rows = []
for case in CASES:
    old, new = (medians[case][v] for v in ("previous", "current"))
    rows.append([LABELS[case], f"{old:.3f}", f"{new:.3f}",
                 f"{100 * (1 - new / old):+.1f}%", f"{old / new:.2f}×"])
table(["Workload", "Previous (s)", "Current (s)", "Time saved", "Speed-up"], rows)
default_cases = ["3d", "3d-maxgamma", "3d-local"]
savings = [100 * (1 - medians[c]["current"] / medians[c]["previous"]) for c in default_cases]
speedups = [medians[c]["previous"] / medians[c]["current"] for c in default_cases]
display(Markdown(
    f"The three **3D cases using the default interpolator** save **{min(savings):.1f}–{max(savings):.1f}%** "
    f"of median runtime on this host (**{min(speedups):.2f}–{max(speedups):.2f}×** speed-up). "
    "Every compared gamma array is elementwise equal, including NaN positions. "
    "The other cases show how the benefit depends on the workload and interpolator."
))

Workload	Previous (s)	Current (s)	Time saved	Speed-up
3D · 3% / 3 mm	5.004	3.744	+25.2%	1.34×
3D · 3% / 3 mm · cap 2	4.985	3.751	+24.8%	1.33×
3D · 2% / 2 mm · local	19.615	15.463	+21.2%	1.27×
3D · 3% / 3 mm · SciPy	15.689	15.666	+0.1%	1.00×
2D · 3% / 3 mm	0.168	0.141	+16.1%	1.19×
2D · 3% / 3 mm · SciPy	0.198	0.194	+1.8%	1.02×

The three 3D cases using the default interpolator save 21.2–25.2% of median runtime on this host (1.27–1.34× speed-up). Every compared gamma array is elementwise equal, including NaN positions. The other cases show how the benefit depends on the workload and interpolator.

Show code cell source

Hide code cell source

def runtime_plot(cases, title):
    fig, ax = plt.subplots(figsize=(8, 0.9 * len(cases) + 1.6), layout="constrained")
    positions = np.arange(len(cases))
    for version, offset, hatch in [("previous", -0.18, "//"), ("current", 0.18, None)]:
        values = [medians[c][version] for c in cases]
        ax.barh(positions + offset, values, height=0.32,
                color=COLOURS[version], hatch=hatch, label=version.capitalize(), alpha=0.9)
        for pos, case, value in zip(positions + offset, cases, values):
            samples = result["times"][case][version]
            ax.plot([min(samples), max(samples)], [pos, pos], color="#222222", lw=1.4)
            ax.text(max(samples) + max(values) * 0.025, pos, f"{value:.3f}", va="center", fontsize=10)
    longest = max(max(result["times"][c][v]) for c in cases for v in ("previous", "current"))
    ax.set(yticks=positions, yticklabels=[LABELS[c] for c in cases],
           xlabel="Warmed gamma-call runtime (s); lower is faster", title=title, xlim=(0, longest * 1.25))
    ax.invert_yaxis()
    ax.grid(axis="x", alpha=0.2)
    ax.legend(loc="upper center", bbox_to_anchor=(0.5, -0.2), ncol=2, frameon=False)
    plt.show()


runtime_plot(CASES[:4], "3D workloads: median runtime and observed range")

Bars start at zero. The thin black lines span the minimum to maximum individual call; they are observed ranges, not confidence intervals. The 2D results use a separate seconds axis so their shorter runtimes remain visible. Do not compare bar lengths across the two figures without reading the axes.

Show code cell source

Hide code cell source

runtime_plot(CASES[4:], "2D workloads: median runtime and observed range")

4. Does the improvement recur across rounds?

Each mark below is the ratio of the two-call median for the current revision to the two-call median for the previous revision in the same round. Ratios below 1 favour the current version; 1 means equal time. The large hollow diamond shows the ratio of the six-call medians used in the table; it need not equal the median of the three round ratios.

These are repeated observations on one machine, not independent patient datasets or independent machines. Alternating order reduces systematic order bias, but three rounds cannot characterise all scheduling, thermal or background-load effects. Small changes near 1 should not be promoted as a reliable benefit or regression from this experiment alone.

Show code cell source

Hide code cell source

fig, ax = plt.subplots(figsize=(8, 5), layout="constrained")
markers = ["o", "s", "^"]
for round_index in range(result["rounds"]):
    ratios = []
    for case in CASES:
        pair = {v: next(r for r in records if r["case"] == case and r["version"] == v and r["round"] == round_index)
                for v in ("previous", "current")}
        ratios.append(np.median(pair["current"]["times"]) / np.median(pair["previous"]["times"]))
    ax.scatter(ratios, np.arange(len(CASES)) + (round_index - 1) * 0.12,
               marker=markers[round_index % len(markers)], s=40,
               label=f"Round {round_index + 1}", zorder=3)
overall = [medians[c]["current"] / medians[c]["previous"] for c in CASES]
ax.scatter(overall, np.arange(len(CASES)), marker="D", s=100, facecolors="none",
           edgecolors="black", label="Ratio of six-call medians", zorder=4)
ax.axvline(1, color="#555555", ls="--", lw=1.2)
ax.set(yticks=np.arange(len(CASES)), yticklabels=[LABELS[c] for c in CASES],
       xlabel="Current / previous runtime; lower is faster", title="Paired rounds and overall runtime ratio")
ax.invert_yaxis()
ax.grid(axis="x", alpha=0.2)
ax.legend(loc="best", fontsize=9)
plt.show()

5. Is the numerical answer unchanged?

Yes, for every call in these six workloads: the comparison uses exact elementwise equality, without an absolute or relative tolerance. It also checks the NaN pattern, rather than comparing only a pass rate or a gamma sum. A NaN is not a passing or failing gamma value. Here all eligible reference points have finite gamma, and all excluded points are below the dose cutoff.

The pass fraction below is the number of finite values with gamma ≤ 1 divided by the number of eligible reference points. Equal pass fractions alone would be a much weaker check than comparing the complete arrays. Results are compared between revisions within each case; the default and SciPy algorithms are not required to be bitwise identical to each other.

Show code cell source

Hide code cell source

rows = []
for case in CASES:
    record = next(r for r in records if r["case"] == case)
    assert record["nan_count"] + record["finite_points"] == record["points"]
    rows.append([LABELS[case], f"{record['points']:,}", f"{record['eligible_points']:,}",
                 f"{record['nan_count']:,}", f"{100 * record['pass_count'] / record['eligible_points']:.2f}%", "Exact"])
table(["Workload", "Grid points", "Eligible / finite", "Excluded (NaN)", "Pass fraction", "Equality"], rows)

Workload	Grid points	Eligible / finite	Excluded (NaN)	Pass fraction	Equality
3D · 3% / 3 mm	133,209	63,935	69,274	97.36%	Exact
3D · 3% / 3 mm · cap 2	133,209	63,935	69,274	97.36%	Exact
3D · 2% / 2 mm · local	133,209	63,935	69,274	72.59%	Exact
3D · 3% / 3 mm · SciPy	133,209	63,935	69,274	97.36%	Exact
2D · 3% / 3 mm	160,801	52,379	108,422	100.00%	Exact
2D · 3% / 3 mm · SciPy	160,801	52,379	108,422	100.00%	Exact

This establishes numerical preservation on these synthetic inputs. It does not establish accuracy for every possible gamma problem or clinical suitability of a chosen criterion. The coordinate and interpolation regression tests cover different questions; see the validation note. In particular, descending or uneven evaluation axes and non-overlapping grids must not be mixed into this speed comparison: the former implementation could fail or fail to terminate on those inputs.

6. Why can the same calculation take less time?

Gamma searches progressively larger shells around reference points and interpolates the evaluation dose at many candidate positions. In 3D each shell is a sampled surface, whereas in 2D it is a sampled circumference. The number of candidate positions, the points still searching, and memory traffic can dominate runtime; the number of voxels alone is a poor predictor. This helps explain why the denser 2D grid can be quicker than the smaller 3D grid.

The changed default path reshapes the already-contiguous query coordinates as a view of the same memory. Previously it copied the coordinates into another array on each shell and RAM chunk. The current path prepares and validates the evaluation grid once, then uses a private interpolation entry point. The old gamma call already used skip_checks=True; it would be misleading to attribute the observed gain solely to removing old axis checks. The SciPy path now constructs its fixed-grid interpolator once per gamma calculation. How much any of these changes matters depends on the work inside interpolation and the shell search.

The example below demonstrates the avoided query-coordinate copy. It checks that the old and new layouts hold exactly the same values and that only the new one shares memory. The displayed bytes refer only to this illustrative coordinate buffer; they are not a measurement of total gamma memory use. This experiment compares the combined revisions, so it cannot assign a percentage of the overall saving to an individual optimisation.

Show code cell source

Hide code cell source

query_points = np.arange(5000 * 50 * 3, dtype=np.float64).reshape(5000, 50, 3)
copied = np.column_stack([query_points[..., i].ravel() for i in range(3)])
view = query_points.reshape(-1, 3)
np.testing.assert_array_equal(copied, view)
assert not np.shares_memory(copied, query_points)
assert np.shares_memory(view, query_points)
table(["Coordinate layout", "Values", "Shares source memory", "Additional coordinate-buffer allocation"], [
    ["Previous: copy", "Equal", "No", f"{copied.nbytes / 1024**2:.2f} MiB"],
    ["Current: view", "Equal", "Yes", "0 MiB (array metadata still exists)"],
])

Coordinate layout	Values	Shares source memory	Additional coordinate-buffer allocation
Previous: copy	Equal	No	5.72 MiB
Current: view	Equal	Yes	0 MiB (array metadata still exists)

7. Repeat the benchmark on your computer

One command, with a chart ready to share

From the root of a checkout containing PR #2066, activate your PyMedPhys Python environment and run:

python examples/gamma_performance.py

The runner creates temporary detached Git worktrees for previous main 866f83e and the current checkout’s committed HEAD, then removes those temporary worktrees when finished. It preserves your working files. Git and the gamma dependencies plus Matplotlib must already be available; no patient data is needed. Both commit objects must be present locally. With a shallow clone, fetch the required revision first. Uncommitted library edits are not benchmarked.

A new, timestamped directory contains comparison.png and comparison.svg, plus results.json, timings.csv and summary.csv. Use the saved PNG directly instead of a screenshot: it includes the revisions, processor, software versions, thread count, timing method and full-array equality result. Keep the JSON with it so someone else can inspect the underlying measurements. The script reports progress by case and round.

To choose another thread count, use --threads 4. Use --previous-ref and --current-ref to select different committed revisions. --help lists options. A short installation check is --cases 2d --rounds 1 --repeats 1; this reduced run is not enough to support the full performance claim. --plot-only path/to/results.json regenerates the figure from saved evidence without timing gamma again. All outputs must go to a new directory to avoid replacing a previous run.

Run or customise this notebook

To regenerate the figures, run all cells with the benchmark path variables unset. Python, NumPy, SciPy, Numba, Matplotlib, IPython and an importable PyMedPhys are needed. The worker and recorded data are embedded, so no benchmark data file or patient download is required.

For fresh timing measurements, use two separate Git checkouts and the same Python environment for both. Copy the full revision IDs from the provenance table above, check out each revision, then set these environment variables before starting the notebook kernel:

Variable	Meaning
PYMEDPHYS_BENCH_PREVIOUS	Existing repository root for the previous revision, not its lib directory
PYMEDPHYS_BENCH_PREVIOUS_SHA	Exact 40-character previous commit ID
PYMEDPHYS_BENCH_CURRENT	Existing repository root for the current revision
PYMEDPHYS_BENCH_CURRENT_SHA	Exact 40-character current commit ID
PYMEDPHYS_BENCH_THREADS	Positive Numba thread count, default 2, used for both revisions

Leave both checkout paths unset to use the embedded record. A partially configured or invalid comparison raises an error; it does not silently fall back to an installed package or relabel the recorded timings as fresh. The controller gives each worker a 600-second timeout, so a much slower host may require an explicit increase. Each case has three untimed warm-ups and six timed calls per revision. Allow time for initial Numba compilation as well.

Use a stable power mode, avoid other substantial computation, record the hardware and repeat the whole experiment. Re-run all cells after changing settings. Inspect the runtime ranges and paired rounds before quoting a result. To retain a new run, save result as JSON as well as the executed notebook; the individual times and provenance are needed to audit the summary. The embedded RECORDED record remains the original experiment until explicitly replaced.

What these results support

This comparison quantifies a useful reduction in warmed 3D gamma time on the measured workloads without changing their answers. It supports a separate performance release-note entry. It does not promise the same percentage on another computer, a different thread count, a cold first call, clinical plans, altered dose thresholds, tighter sampling or different RAM chunking. No peak-memory or isolated-optimisation benchmark was performed. For every case, read the absolute savings and observed variation alongside the median ratio.