CDA Calibration Notebook#
This notebook is designed to calibrate the parameters of the CDA model using a dataset of GPS tracks from cycling activities.
The goal is to find the optimal values for the parameters that minimize the difference between the predicted and actual speeds of the cyclist.
What this Notebook does:
uses the biwipy package interfaces ,
calibrates CdA per GPX case against a target average power,
uses hourly GRIB selection (pas=1) for better temporal precision,
supports an explicit list of GPX files with their corresponding start datetimes,
reports per-case results + mean CdA.
[2]:
import os
os.environ["OUTPUT_LANG"] = "en"
print("OUTPUT_LANG forced to en")
home = os.path.expanduser("~")
OUTPUT_LANG forced to en
[3]:
# required imports for the notebook, including core simulation components, weather data handling, analysis tools, and visualization functions.
from dataclasses import dataclass
from datetime import datetime
from statistics import mean
import sys
from typing import List, Optional
import os
from zoneinfo import ZoneInfo
from biwipy.analysis import RouteAnalyzer
from biwipy.core import CyclistBehavior, Simulator
from biwipy.core.bike_physics import calibrate_cda_from_power
from biwipy.weather import WeatherProvider
from biwipy.weather.grib_finder import build_grib_list
Parameters and settings (optional customization)#
Customize the paths if necessary and other parmeters if you want
[4]:
HDIR = home+"/data" # GRIB files stored in /home/data
GPX_DIR= home+"/data/gpx/" # GPX files stored in /home/data/gpx
ORIGIN_TZ = ZoneInfo("Europe/Paris")
# Physics parameters
MASS_KG = 98
CR = 0.0055
# CdA search settings
CDA_MIN = 0.30
CDA_MAX = 0.80
TOLERANCE_W = 1.0
MAX_ITER = 15
# Weather interpolation settings
GRIB_STEP_H = 1 # hourly precision
GRIB_INTERVAL_H = 9 # estimated ride duration for interpolation coverage
# Behavior used for calibration and replay consistency
BEHAVIOR = CyclistBehavior(uphill="realistic", downhill="realistic", corner="realistic")
Useful routines (no customization needed)#
Full preprocessing
Calibrate for one file
[5]:
@dataclass
class CalibrationCase:
gpx_file: str
target_power_w: float
def _build_segments(filegpx: str):
"""Build cleaned/merged segments using RouteAnalyzer workflow."""
analyzer = RouteAnalyzer(
smoothing_window=11,
smoothing_method="mediane",
merge_min_distance=50.0,
merge_max_bearing_diff=20.0,
merge_max_slope_diff=0.10,
merge_max_slope=0.15,
)
gpx_result = analyzer.process_gpx(
filegpx,
smooth=True,
detect_noise=True,
remove_noise=True,
merge_segments_flag=True,
filter_stops_flag=False,
max_dist_noise=5.0,
min_slope_threshold_noise=0.10,
normal_slope_threshold_noise=0.05,
verbose=False,
)
return gpx_result
def _run_single_case(case: CalibrationCase) -> dict:
filegpx = os.path.join(GPX_DIR, case.gpx_file)
#t_start = datetime.fromisoformat(case.start_datetime).replace(tzinfo=ORIGIN_TZ)
print("\n" + "=" * 72)
print(f" CALIBRATION CdA - {case.gpx_file}")
print(f" Target power: {case.target_power_w:.1f} W")
print("=" * 72)
print("📍 Preprocessing GPX via RouteAnalyzer...")
gpx_result = _build_segments(filegpx)
segments=gpx_result.segments
stats=gpx_result.stats
if gpx_result.has_timestamps:
t_start=gpx_result.t_start.astimezone(ORIGIN_TZ)
print(f"Starting date (locale) : {t_start} {t_start.tzname()}")
else:
print("⚠️ Gpx file does not contain timestamps!")
raise ValueError("Gpx file must contain timestamps for replay.")
sys.exit(1)
if "distance_km" in gpx_result.stats:
print(f" Route distance: {gpx_result.stats['distance_km']:.3f} km")
print("📡 Loading GRIB files (hourly interpolation)...")
gribs_list = build_grib_list(HDIR, t_start, GRIB_STEP_H, GRIB_INTERVAL_H)
weather = WeatherProvider(gribs_list, bcache=True)
mygrib = weather.grib
print("🎯 Running CdA calibration...\n")
print(
f" behavior: {BEHAVIOR.mode_uphill}/{BEHAVIOR.mode_downhill}/{BEHAVIOR.mode_corner}"
)
cda_opt, p0_opt, power_opt = calibrate_cda_from_power(
segments_in=gpx_result.segments,
grib=mygrib,
t_start=t_start,
target_power=case.target_power_w,
Cr=CR,
m=MASS_KG,
cda_min=CDA_MIN,
cda_max=CDA_MAX,
tolerance=TOLERANCE_W,
max_iterations=MAX_ITER,
behavior=BEHAVIOR,
v_max=22.0,
ratio_wind=0.25,
)
# Optional consistency replay check
sim_replay = Simulator(mygrib, behavior=BEHAVIOR, CdA=cda_opt, Cr=CR, m=MASS_KG)
replay_result = sim_replay.simulate_replay(segments, t_start, passes=2)
p0_replay = replay_result.power.P0_calibrated if replay_result.power else None
print("\n" + "-" * 72)
print(f"CdA optimal : {cda_opt:.4f} m²")
print(f"P0 optimal : {p0_opt:.1f} W")
print(f"Power obtained : {power_opt:.1f} W")
print(f"Power target : {case.target_power_w:.1f} W")
print(f"Delta : {power_opt - case.target_power_w:+.1f} W")
if p0_replay is not None:
print(f"Replay P0 check : {p0_replay:.1f} W")
print("-" * 72)
return {
"gpx_file": case.gpx_file,
"start_datetime": t_start,
"target_power_w": case.target_power_w,
"cda_opt": cda_opt,
"p0_opt": p0_opt,
"power_opt": power_opt,
"delta_w": power_opt - case.target_power_w,
"p0_replay": p0_replay,
}
Definition of the calibration case (mandatory customization)#
Customize this part with your gpx files (with timestamps) and their associated power (power meters measure, strava measurement , … )
NB: the 0.95 factor used is to into account power is lost in drivetrain transmission.
[6]:
# Explicit cases: GPX + corresponding start datetime + target power
CALIBRATION_CASES: List[CalibrationCase] = [
CalibrationCase(
gpx_file="S20260218-049K-0958.gpx",
target_power_w=156.0 *0.95,
),
CalibrationCase(
gpx_file="S20251213-045K-1506.gpx",
target_power_w=177.0*0.95,
),
]
Calibration (no customization needed)#
[7]:
# list the files in the GPX_DIR directory
cases = CALIBRATION_CASES
results = [_run_single_case(case) for case in cases]
cda_values = [r["cda_opt"] for r in results]
p0_values = [r["p0_opt"] for r in results]
delta_values = [r["delta_w"] for r in results]
print("\n" + "=" * 72)
print(" GLOBAL SUMMARY")
print("=" * 72)
print(f"Cases : {len(results)}")
print(f"CdA mean : {mean(cda_values):.4f} m²")
print(f"CdA min / max : {min(cda_values):.4f} / {max(cda_values):.4f}")
print(f"P0 mean : {mean(p0_values):.1f} W")
print(f"Delta power mean : {mean(delta_values):+.2f} W")
print("=" * 72)
========================================================================
CALIBRATION CdA - S20260218-049K-0958.gpx
Target power: 148.2 W
========================================================================
📍 Preprocessing GPX via RouteAnalyzer...
Starting date (locale) : 2026-02-18 09:58:14+01:00 CET
Route distance: 49.743 km
📡 Loading GRIB files (hourly interpolation)...
🎯 Running CdA calibration...
behavior: realistic/realistic/realistic
WARNING:root:✅ P0 calibrated: 148.4W (reproduces 148.5W observed)
WARNING:root:✅ P0 calibrated: 148.4W (reproduces 148.5W observed)
------------------------------------------------------------------------
CdA optimal : 0.5266 m²
P0 optimal : 148.4 W
Power obtained : 148.5 W
Power target : 148.2 W
Delta : +0.3 W
Replay P0 check : 148.4 W
------------------------------------------------------------------------
========================================================================
CALIBRATION CdA - S20251213-045K-1506.gpx
Target power: 168.2 W
========================================================================
📍 Preprocessing GPX via RouteAnalyzer...
Starting date (locale) : 2025-12-13 15:06:06+01:00 CET
Route distance: 45.075 km
📡 Loading GRIB files (hourly interpolation)...
🎯 Running CdA calibration...
behavior: realistic/realistic/realistic
WARNING:root:✅ P0 calibrated: 167.6W (reproduces 167.8W observed)
WARNING:root:✅ P0 calibrated: 167.6W (reproduces 167.8W observed)
------------------------------------------------------------------------
CdA optimal : 0.5188 m²
P0 optimal : 167.6 W
Power obtained : 167.8 W
Power target : 168.2 W
Delta : -0.4 W
Replay P0 check : 167.6 W
------------------------------------------------------------------------
========================================================================
GLOBAL SUMMARY
========================================================================
Cases : 2
CdA mean : 0.5227 m²
CdA min / max : 0.5188 / 0.5266
P0 mean : 158.0 W
Delta power mean : -0.07 W
========================================================================