akima

Author: CamDavidsonPilon

core/pioreactor/calibrations/protocols/od_fusion_offset.py

@@ -102,25 +102,28 @@ def _load_estimator_for_worker(worker: str, estimator_name: str) -> structs.ODFu
     return json_decode(json_encode(payload), type=structs.ODFusionEstimator)
 
 
-def _affine_transform_spline_fit_data(
-    spline_data: structs.SplineFitData,
+def _affine_transform_cubic_fit_data(
+    curve_data: structs.AkimaFitData | structs.SplineFitData,
     scale_logc: float,
     offset_logc: float,
-) -> structs.SplineFitData:
+) -> structs.AkimaFitData | structs.SplineFitData:
     if scale_logc <= 0:
-        raise ValueError("Scale must be positive to transform spline data.")
-    return structs.SplineFitData(
-        knots=[float(scale_logc * knot + offset_logc) for knot in spline_data.knots],
-        coefficients=[
-            [
-                float(coeffs[0]),
-                float(coeffs[1] / scale_logc),
-                float(coeffs[2] / scale_logc**2),
-                float(coeffs[3] / scale_logc**3),
-            ]
-            for coeffs in spline_data.coefficients
-        ],
-    )
+        raise ValueError("Scale must be positive to transform curve data.")
+
+    knots = [float(scale_logc * knot + offset_logc) for knot in curve_data.knots]
+    coefficients = [
+        [
+            float(coeffs[0]),
+            float(coeffs[1] / scale_logc),
+            float(coeffs[2] / scale_logc**2),
+            float(coeffs[3] / scale_logc**3),
+        ]
+        for coeffs in curve_data.coefficients
+    ]
+
+    if isinstance(curve_data, structs.AkimaFitData):
+        return structs.AkimaFitData(knots=knots, coefficients=coefficients)
+    return structs.SplineFitData(knots=knots, coefficients=coefficients)
 
 
 def _apply_logc_affine_to_estimator(
@@ -135,11 +138,11 @@ def _apply_logc_affine_to_estimator(
     source_estimator_name: str,
 ) -> structs.ODFusionEstimator:
     mu_splines = {
-        angle: _affine_transform_spline_fit_data(estimator.mu_splines[angle], scale_logc, offset_logc)
+        angle: _affine_transform_cubic_fit_data(estimator.mu_splines[angle], scale_logc, offset_logc)
         for angle in estimator.angles
     }
     sigma_splines_log = {
-        angle: _affine_transform_spline_fit_data(estimator.sigma_splines_log[angle], scale_logc, offset_logc)
+        angle: _affine_transform_cubic_fit_data(estimator.sigma_splines_log[angle], scale_logc, offset_logc)
         for angle in estimator.angles
     }
 

core/pioreactor/calibrations/protocols/od_fusion_standards.py

@@ -425,8 +425,10 @@ def advance(self, ctx: SessionContext) -> SessionStep | None:
             recorded_data=fit.recorded_data,
             ir_led_intensity=float(config["od_reading.config"]["ir_led_intensity"]),
             angles=list(FUSION_ANGLES),
-            mu_splines=fit.mu_splines,
-            sigma_splines_log=fit.sigma_splines_log,
+            mu_splines=cast(dict[pt.PdAngle, structs.AkimaFitData | structs.SplineFitData], fit.mu_splines),
+            sigma_splines_log=cast(
+                dict[pt.PdAngle, structs.AkimaFitData | structs.SplineFitData], fit.sigma_splines_log
+            ),
             min_logc=fit.min_logc,
             max_logc=fit.max_logc,
             sigma_floor=fit.sigma_floor,

core/pioreactor/calibrations/protocols/od_standards.py

@@ -171,11 +171,9 @@ def _calculate_curve_data(
     od600_values: list[float],
     voltages: list[float],
 ) -> structs.CalibrationCurveData:
-    weights = [1.0] * len(voltages)
-    weights[0] = len(voltages) / 2
-    from pioreactor.utils.splines import spline_fit
+    from pioreactor.utils.akimas import akima_fit
 
-    return spline_fit(od600_values, voltages, knots="auto", weights=weights)
+    return akima_fit(od600_values, voltages)
 
 
 def _build_standards_chart_metadata(

core/pioreactor/calibrations/protocols/stirring_dc_based.py

@@ -113,15 +113,15 @@ def _build_stirring_calibration_from_measurements(
     logger = create_logger("stirring_calibration", experiment="$experiment")
     logger.debug(f"rpm = {alpha:.2f} * dc% + {beta:.2f}")
 
-    from pioreactor.utils.splines import spline_fit
+    from pioreactor.utils.akimas import akima_fit
 
     return SimpleStirringCalibration(
         pwm_hz=config.getfloat("stirring.config", "pwm_hz"),
         voltage=voltage,
         calibration_name=f"stirring-calibration-{current_utc_datetime().strftime('%Y-%m-%d_%H-%M')}",
         calibrated_on_pioreactor_unit=unit,
         created_at=current_utc_datetime(),
-        curve_data_=spline_fit(dcs, rpms, knots="auto"),
+        curve_data_=akima_fit(dcs, rpms),
         recorded_data={"x": dcs, "y": rpms},
     )
 

core/pioreactor/calibrations/utils.py

@@ -44,6 +44,9 @@ def curve_to_functional_form(curve_data: structs.CalibrationCurveData) -> str:
     elif curve_type == "spline":
         spline_data = cast(structs.SplineFitData, curve_data)
         return f"natural cubic spline with {len(spline_data.knots)} knots"
+    elif curve_type == "akima":
+        akima_data = cast(structs.AkimaFitData, curve_data)
+        return f"Akima interpolator with {len(akima_data.knots)} knots"
     else:
         raise NotImplementedError()
 
@@ -68,6 +71,16 @@ def curve_callable(x: float):
 
         return curve_callable
 
+    elif curve_type == "akima":
+        from pioreactor.utils.akimas import akima_eval
+
+        akima_data = cast(structs.AkimaFitData, curve_data)
+
+        def curve_callable(x: float):
+            return akima_eval(akima_data, x)
+
+        return curve_callable
+
     else:
         raise NotImplementedError()
 
@@ -144,6 +157,8 @@ def crunch_data_and_confirm_with_user(
             candidate_curve = structs.PolyFitCoefficients(coefficients=[])
         elif fit == "spline":
             candidate_curve = structs.SplineFitData(knots=[], coefficients=[])
+        elif fit == "akima":
+            candidate_curve = structs.AkimaFitData(knots=[], coefficients=[])
         else:
             raise ValueError()
 
@@ -169,8 +184,12 @@ def crunch_data_and_confirm_with_user(
                     initial_knots = knots
 
                 candidate_curve = spline_fit(x, y, knots=knots, weights=weights)
+            elif fit == "akima":
+                from pioreactor.utils.akimas import akima_fit
+
+                candidate_curve = akima_fit(x, y)
             else:
-                raise ValueError("only `poly` or `spline` supported")
+                raise ValueError("only `poly`, `spline`, or `akima` supported")
 
         curve_callable = curve_to_callable(candidate_curve)
         plot_data(
@@ -229,6 +248,8 @@ def _fit_prompt_choices(fit: str) -> list[str]:
         return ["y", "q", "d"]
     if fit == "spline":
         return ["y", "q", "k"]
+    if fit == "akima":
+        return ["y", "q"]
     return ["y", "q"]
 
 
@@ -239,6 +260,8 @@ def _fit_prompt_hint(fit: str, candidate_curve: structs.CalibrationCurveData) ->
     if fit == "spline":
         spline_data = cast(structs.SplineFitData, candidate_curve)
         return f"k: choose a new knot count (currently {len(spline_data.knots)})"
+    if fit == "akima":
+        return ""
     return ""
 
 
@@ -249,4 +272,7 @@ def _curve_data_is_empty(fit: str, curve_data: structs.CalibrationCurveData) ->
     if fit == "spline":
         spline_data = cast(structs.SplineFitData, curve_data)
         return len(spline_data.knots) == 0 or len(spline_data.coefficients) == 0
+    if fit == "akima":
+        akima_data = cast(structs.AkimaFitData, curve_data)
+        return len(akima_data.knots) == 0 or len(akima_data.coefficients) == 0
     return True

core/pioreactor/cli/calibrations.py

@@ -285,7 +285,7 @@ def delete_calibration(device: str, calibration_name: str) -> None:
     "--fit",
     "fit",
     default="poly",
-    type=click.Choice(["poly", "spline"]),
+    type=click.Choice(["poly", "spline", "akima"]),
     show_default=True,
     help="Curve fit type to use when analyzing.",
 )

core/pioreactor/structs.py

@@ -59,7 +59,16 @@ def type(self):
         return self.__struct_config__.tag
 
 
-type CalibrationCurveData = PolyFitCoefficients | SplineFitData
+class AkimaFitData(Struct, tag="akima"):
+    knots: list[float]
+    coefficients: list[list[float]]
+
+    @property
+    def type(self):
+        return self.__struct_config__.tag
+
+
+type CalibrationCurveData = PolyFitCoefficients | SplineFitData | AkimaFitData
 
 
 class AutomationSettings(JSONPrintedStruct):
@@ -272,6 +281,13 @@ def x_to_y(self, x: X) -> Y:
             if len(spline_data.knots) == 0 or len(spline_data.coefficients) == 0:
                 raise exc.NoSolutionsFoundError(f"calibration {self}'s curve_data_ is empty")
             return round(spline_eval(spline_data, x), 10)
+        if self.curve_data_.type == "akima":
+            from pioreactor.utils.akimas import akima_eval
+
+            akima_data = t.cast(AkimaFitData, self.curve_data_)
+            if len(akima_data.knots) == 0 or len(akima_data.coefficients) == 0:
+                raise exc.NoSolutionsFoundError(f"calibration {self}'s curve_data_ is empty")
+            return round(akima_eval(akima_data, x), 10)
 
         raise NotImplementedError(f"Unsupported curve_type: {self.curve_data_.type}")
 
@@ -295,6 +311,13 @@ def y_to_x(self, y: Y, enforce_bounds=False) -> X:
             if len(spline_data.knots) == 0 or len(spline_data.coefficients) == 0:
                 raise exc.NoSolutionsFoundError(f"calibration {self}'s curve_data_ is empty")
             plausible_sols_ = spline_solve(spline_data, y)
+        elif self.curve_data_.type == "akima":
+            from pioreactor.utils.akimas import akima_solve
+
+            akima_data = t.cast(AkimaFitData, self.curve_data_)
+            if len(akima_data.knots) == 0 or len(akima_data.coefficients) == 0:
+                raise exc.NoSolutionsFoundError(f"calibration {self}'s curve_data_ is empty")
+            plausible_sols_ = akima_solve(akima_data, y)
         else:
             raise NotImplementedError(f"Unsupported curve_type: {self.curve_data_.type}")
 
@@ -428,8 +451,8 @@ class OD600Calibration(ODCalibration, kw_only=True, tag="od600"):
 class ODFusionEstimator(EstimatorBase, kw_only=True, tag="od_fused_estimator"):
     ir_led_intensity: float
     angles: list[pt.PdAngle]
-    mu_splines: dict[pt.PdAngle, SplineFitData]
-    sigma_splines_log: dict[pt.PdAngle, SplineFitData]
+    mu_splines: dict[pt.PdAngle, AkimaFitData | SplineFitData]
+    sigma_splines_log: dict[pt.PdAngle, AkimaFitData | SplineFitData]
     min_logc: float
     max_logc: float
     sigma_floor: float