import numpy as np
from perda.core_data_structures.data_instance import (
DataInstance,
left_join_data_instances,
)
from perda.core_data_structures.single_run_data import SingleRunData
from perda.units import _from_seconds, _to_seconds
from suboptimumg.constants import RHO, G
from suboptimumg.vehicle.irl_vehicle import IrlCar
from .kinematics import GROUND_SPEED
from .utils import safe_gradient
CLA = "aero.cla"
COP = "aero.cop"
CDA = "aero.cda"
[docs]
def add_cla(
data: SingleRunData,
irl_car: IrlCar,
ref_time_range: tuple[float, float],
pot_fl: str = "ludwig.shockpot.frontLeft",
pot_fr: str = "ludwig.shockpot.frontRight",
pot_rl: str = "ludwig.shockpot.rearLeft",
pot_rr: str = "ludwig.shockpot.rearRight",
ax_col: str = "pcm.vnav.linearAccelBody.x",
low_speed_thresh: float = 5.0,
) -> SingleRunData:
"""Add aerodynamic lift coefficient (CLA) and center of pressure (CoP).
Estimates per-axle downforce from shock pot compression relative to a
zero-aero reference window, subtracting longitudinal weight transfer.
Stores results in ``aero.cla`` and ``aero.cop``, accesible via the constants
``CLA`` and ``COP``.
Parameters
----------
irl_car : IrlCar
ref_time_range : tuple[float, float]
``(t_start, t_end)`` in seconds of a zero-aero reference window
(e.g. car stationary or slow).
low_speed_thresh : float
Speed (m/s) below which CLA and CoP are zeroed / NaN.
"""
missing = [c for c in (GROUND_SPEED,) if c not in data]
if missing:
raise KeyError(
f"add_cla: missing variables(s) {missing}.\n"
f"Try running add_groundspeed first."
)
pfl = data[pot_fl]
pfl, pfr_aln, prl_aln, prr_aln, ax_aln, vel_aln = left_join_data_instances(
pfl,
[data[pot_fr], data[pot_rl], data[pot_rr], data[ax_col], data[GROUND_SPEED]],
)
ts = pfl.timestamp_np
front_avg = (pfl.value_np + pfr_aln.value_np) / 2.0
rear_avg = (prl_aln.value_np + prr_aln.value_np) / 2.0
t0_raw = _from_seconds(ref_time_range[0], data.timestamp_unit)
t1_raw = _from_seconds(ref_time_range[1], data.timestamp_unit)
ref_mask = (ts >= t0_raw) & (ts <= t1_raw)
front_delta = front_avg - np.mean(front_avg[ref_mask])
rear_delta = rear_avg - np.mean(rear_avg[ref_mask])
f_front, f_rear = irl_car.aero_force_from_pots(
front_delta, rear_delta, ax_aln.value_np / G
)
f_total = f_front + f_rear
v = vel_aln.value_np
q_A = 0.5 * RHO * v**2 * irl_car.params.aero.front_area
cla = np.zeros_like(v, dtype=np.float64)
cla_mask = np.abs(v) > low_speed_thresh
cla[cla_mask] = f_total[cla_mask] / q_A[cla_mask]
cop = np.full_like(v, np.nan, dtype=np.float64)
cop_mask = np.abs(f_total) > 1.0
cop[cop_mask] = f_front[cop_mask] / f_total[cop_mask]
data[CLA] = DataInstance(
timestamp_np=ts,
value_np=cla,
label="Aerodynamic CLA",
cpp_name=CLA,
)
data[COP] = DataInstance(
timestamp_np=ts,
value_np=cop,
label="Aerodynamic Center of Pressure (front fraction)",
cpp_name=COP,
)
return data
[docs]
def add_cda(
data: SingleRunData,
irl_car: IrlCar,
ax_col: str = "pcm.vnav.linearAccelBody.x",
low_speed_thresh: float = 5.0,
) -> SingleRunData:
"""Add aerodynamic drag coefficient (CDA) from coastdown deceleration.
Applies ``F_drag = -m·a_x - F_rolling``, then ``CDA = F_drag / (0.5·ρ·v²)``.
For meaningful results trim ``data`` to a coasting segment first (no throttle or brake).
Requires ``add_groundspeed``.
Stores results in ``aero.cda``, accessible via the constant ``CDA``.
Parameters
----------
data: SingleRunData
irl_car : IrlCar
ax_col : str
Column name of longitudinal acceleration in body frame.
low_speed_thresh : float
Speed (m/s) below which CDA is forced to zero.
"""
missing = [c for c in (GROUND_SPEED, ax_col) if c not in data]
if missing:
raise KeyError(
f"add_cda: missing variables(s) {missing}.\n"
f"Try running add_groundspeed first."
)
gs = data[GROUND_SPEED]
t_s = _to_seconds(gs.timestamp_np, data.timestamp_unit)
v = gs.value_np
a_x = safe_gradient(v, np.gradient(t_s))
m = irl_car.params.mass
f_drag = -m * a_x - m * G * irl_car.params.rolling_coeff
q = 0.5 * RHO * v**2
with np.errstate(divide="ignore", invalid="ignore"):
cda = np.where(np.abs(v) > low_speed_thresh, f_drag / q, 0.0)
data[CDA] = DataInstance(
timestamp_np=gs.timestamp_np,
value_np=cda,
label="Aerodynamic CDA",
cpp_name=CDA,
)
return data