src/forward_model/sensor_grid.rs

branch
dev
changeset 61
4f468d35fa29
parent 49
6b0db7251ebe
--- a/src/forward_model/sensor_grid.rs	Sun Apr 27 15:03:51 2025 -0500
+++ b/src/forward_model/sensor_grid.rs	Thu Feb 26 11:38:43 2026 -0500
@@ -2,13 +2,17 @@
 Sensor grid forward model
 */
 
-use nalgebra::base::{DMatrix, DVector};
-use numeric_literals::replace_float_literals;
-use std::iter::Zip;
-use std::ops::RangeFrom;
-
+use super::{BasicCurvatureBoundEstimates, ForwardModel};
+use crate::frank_wolfe::FindimQuadraticModel;
+use crate::kernels::{AutoConvolution, BoundedBy, Convolution};
+use crate::measures::{DiscreteMeasure, Radon, RNDM};
+use crate::preadjoint_helper::PreadjointHelper;
+use crate::seminorms::{ConvolutionOp, SimpleConvolutionKernel};
+use crate::types::*;
 use alg_tools::bisection_tree::*;
-use alg_tools::error::DynError;
+use alg_tools::bounds::Bounded;
+use alg_tools::error::{DynError, DynResult};
+use alg_tools::euclidean::Euclidean;
 use alg_tools::instance::Instance;
 use alg_tools::iter::{MapX, Mappable};
 use alg_tools::lingrid::*;
@@ -18,79 +22,74 @@
 use alg_tools::nalgebra_support::ToNalgebraRealField;
 use alg_tools::norms::{Linfinity, Norm, L1, L2};
 use alg_tools::tabledump::write_csv;
+use anyhow::anyhow;
+use nalgebra::base::{DMatrix, DVector};
+use numeric_literals::replace_float_literals;
+use std::iter::Zip;
+use std::ops::RangeFrom;
 
-use super::{AdjointProductBoundedBy, BoundedCurvature, ForwardModel};
-use crate::frank_wolfe::FindimQuadraticModel;
-use crate::kernels::{AutoConvolution, BoundedBy, Convolution};
-use crate::measures::{DiscreteMeasure, Radon};
-use crate::preadjoint_helper::PreadjointHelper;
-use crate::seminorms::{ConvolutionOp, SimpleConvolutionKernel};
-use crate::types::*;
-
-type RNDM<F, const N: usize> = DiscreteMeasure<Loc<F, N>, F>;
-
-pub type ShiftedSensor<F, S, P, const N: usize> = Shift<Convolution<S, P>, F, N>;
+pub type ShiftedSensor<F, S, P, const N: usize> = Shift<Convolution<S, P>, N, F>;
 
 /// Trait for physical convolution models. Has blanket implementation for all cases.
-pub trait Spread<F: Float, const N: usize>:
-    'static + Clone + Support<F, N> + RealMapping<F, N> + Bounded<F>
+pub trait Spread<const N: usize, F: Float = f64>:
+    'static + Clone + Support<N, F> + RealMapping<N, F> + Bounded<F>
 {
 }
 
-impl<F, T, const N: usize> Spread<F, N> for T
+impl<F, T, const N: usize> Spread<N, F> for T
 where
     F: Float,
-    T: 'static + Clone + Support<F, N> + Bounded<F> + RealMapping<F, N>,
+    T: 'static + Clone + Support<N, F> + Bounded<F> + RealMapping<N, F>,
 {
 }
 
 /// Trait for compactly supported sensors. Has blanket implementation for all cases.
-pub trait Sensor<F: Float, const N: usize>:
-    Spread<F, N> + Norm<F, L1> + Norm<F, Linfinity>
+pub trait Sensor<const N: usize, F: Float = f64>:
+    Spread<N, F> + Norm<L1, F> + Norm<Linfinity, F>
 {
 }
 
-impl<F, T, const N: usize> Sensor<F, N> for T
+impl<F, T, const N: usize> Sensor<N, F> for T
 where
     F: Float,
-    T: Spread<F, N> + Norm<F, L1> + Norm<F, Linfinity>,
+    T: Spread<N, F> + Norm<L1, F> + Norm<Linfinity, F>,
 {
 }
 
 pub trait SensorGridBT<F, S, P, const N: usize>:
-    Clone + BTImpl<F, N, Data = usize, Agg = Bounds<F>>
+    Clone + BTImpl<N, F, Data = usize, Agg = Bounds<F>>
 where
     F: Float,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
 {
 }
 
 impl<F, S, P, T, const N: usize> SensorGridBT<F, S, P, N> for T
 where
-    T: Clone + BTImpl<F, N, Data = usize, Agg = Bounds<F>>,
+    T: Clone + BTImpl<N, F, Data = usize, Agg = Bounds<F>>,
     F: Float,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
 {
 }
 
 // We need type alias bounds to access associated types
 #[allow(type_alias_bounds)]
 pub type SensorGridBTFN<F, S, P, BT: SensorGridBT<F, S, P, N>, const N: usize> =
-    BTFN<F, SensorGridSupportGenerator<F, S, P, N>, BT, N>;
+    BTFN<F, SensorGridSupportGenerator<S, F, P, N>, BT, N>;
 
 /// Sensor grid forward model
 #[derive(Clone)]
 pub struct SensorGrid<F, S, P, BT, const N: usize>
 where
     F: Float,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
     BT: SensorGridBT<F, S, P, N>,
 {
-    domain: Cube<F, N>,
+    domain: Cube<N, F>,
     sensor_count: [usize; N],
     sensor: S,
     spread: P,
@@ -102,16 +101,16 @@
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, BT::Agg, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F> + LocalAnalysis<F, BT::Agg, N>,
 {
     /// Create a new sensor grid.
     ///
     /// The parameter `depth` indicates the search depth of the created [`BT`]s
     /// for the adjoint values.
     pub fn new(
-        domain: Cube<F, N>,
+        domain: Cube<N, F>,
         sensor_count: [usize; N],
         sensor: S,
         spread: P,
@@ -141,12 +140,12 @@
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
     /// Return the grid of sensor locations.
-    pub fn grid(&self) -> LinGrid<F, N> {
+    pub fn grid(&self) -> LinGrid<N, F> {
         lingrid_centered(&self.domain, &self.sensor_count)
     }
 
@@ -157,7 +156,7 @@
 
     /// Constructs a sensor shifted by `x`.
     #[inline]
-    fn shifted_sensor(&self, x: Loc<F, N>) -> ShiftedSensor<F, S, P, N> {
+    fn shifted_sensor(&self, x: Loc<N, F>) -> ShiftedSensor<F, S, P, N> {
         self.base_sensor.clone().shift(x)
     }
 
@@ -180,44 +179,44 @@
     }
 }
 
-impl<F, S, P, BT, const N: usize> Mapping<RNDM<F, N>> for SensorGrid<F, S, P, BT, N>
+impl<F, S, P, BT, const N: usize> Mapping<RNDM<N, F>> for SensorGrid<F, S, P, BT, N>
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
     type Codomain = DVector<F>;
 
     #[inline]
-    fn apply<I: Instance<RNDM<F, N>>>(&self, μ: I) -> DVector<F> {
+    fn apply<I: Instance<RNDM<N, F>>>(&self, μ: I) -> DVector<F> {
         let mut y = self._zero_observable();
         self.apply_add(&mut y, μ);
         y
     }
 }
 
-impl<F, S, P, BT, const N: usize> Linear<RNDM<F, N>> for SensorGrid<F, S, P, BT, N>
+impl<F, S, P, BT, const N: usize> Linear<RNDM<N, F>> for SensorGrid<F, S, P, BT, N>
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
 }
 
 #[replace_float_literals(F::cast_from(literal))]
-impl<F, S, P, BT, const N: usize> GEMV<F, RNDM<F, N>, DVector<F>> for SensorGrid<F, S, P, BT, N>
+impl<F, S, P, BT, const N: usize> GEMV<F, RNDM<N, F>, DVector<F>> for SensorGrid<F, S, P, BT, N>
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
-    fn gemv<I: Instance<RNDM<F, N>>>(&self, y: &mut DVector<F>, α: F, μ: I, β: F) {
+    fn gemv<I: Instance<RNDM<N, F>>>(&self, y: &mut DVector<F>, α: F, μ: I, β: F) {
         let grid = self.grid();
         if β == 0.0 {
             y.fill(0.0)
@@ -227,38 +226,42 @@
         if α == 1.0 {
             self.apply_add(y, μ)
         } else {
-            for δ in μ.ref_instance() {
-                for &d in self.bt.iter_at(&δ.x) {
-                    let sensor = self.shifted_sensor(grid.entry_linear_unchecked(d));
-                    y[d] += sensor.apply(&δ.x) * (α * δ.α);
+            μ.eval_ref(|μr| {
+                for δ in μr {
+                    for &d in self.bt.iter_at(&δ.x) {
+                        let sensor = self.shifted_sensor(grid.entry_linear_unchecked(d));
+                        y[d] += sensor.apply(&δ.x) * (α * δ.α);
+                    }
                 }
-            }
+            })
         }
     }
 
-    fn apply_add<I: Instance<RNDM<F, N>>>(&self, y: &mut DVector<F>, μ: I) {
+    fn apply_add<I: Instance<RNDM<N, F>>>(&self, y: &mut DVector<F>, μ: I) {
         let grid = self.grid();
-        for δ in μ.ref_instance() {
-            for &d in self.bt.iter_at(&δ.x) {
-                let sensor = self.shifted_sensor(grid.entry_linear_unchecked(d));
-                y[d] += sensor.apply(&δ.x) * δ.α;
+        μ.eval_ref(|μr| {
+            for δ in μr {
+                for &d in self.bt.iter_at(&δ.x) {
+                    let sensor = self.shifted_sensor(grid.entry_linear_unchecked(d));
+                    y[d] += sensor.apply(&δ.x) * δ.α;
+                }
             }
-        }
+        })
     }
 }
 
-impl<F, S, P, BT, const N: usize> BoundedLinear<RNDM<F, N>, Radon, L2, F>
+impl<F, S, P, BT, const N: usize> BoundedLinear<RNDM<N, F>, Radon, L2, F>
     for SensorGrid<F, S, P, BT, N>
 where
     F: Float,
-    BT: SensorGridBT<F, S, P, N, Agg = Bounds<F>>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, BT::Agg, N>,
+    BT: SensorGridBT<F, S, P, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
     /// An estimate on the operator norm in $𝕃(ℳ(Ω); ℝ^n)$ with $ℳ(Ω)$ equipped
     /// with the Radon norm, and $ℝ^n$ with the Euclidean norm.
-    fn opnorm_bound(&self, _: Radon, _: L2) -> F {
+    fn opnorm_bound(&self, _: Radon, _: L2) -> DynResult<F> {
         // With {x_i}_{i=1}^n the grid centres and φ the kernel, we have
         // |Aμ|_2 = sup_{|z|_2 ≤ 1} ⟨z,Αμ⟩ = sup_{|z|_2 ≤ 1} ⟨A^*z|μ⟩
         // ≤ sup_{|z|_2 ≤ 1} |A^*z|_∞ |μ|_ℳ
@@ -271,22 +274,22 @@
         // = |φ|_∞ √N_ψ |μ|_ℳ.
         // Hence
         let n = self.max_overlapping();
-        self.base_sensor.bounds().uniform() * n.sqrt()
+        Ok(self.base_sensor.bounds().uniform() * n.sqrt())
     }
 }
 
-type SensorGridPreadjoint<'a, A, F, const N: usize> = PreadjointHelper<'a, A, RNDM<F, N>>;
+type SensorGridPreadjoint<'a, A, F, const N: usize> = PreadjointHelper<'a, A, RNDM<N, F>>;
 
-impl<F, S, P, BT, const N: usize> Preadjointable<RNDM<F, N>, DVector<F>>
+impl<F, S, P, BT, const N: usize> Preadjointable<RNDM<N, F>, DVector<F>>
     for SensorGrid<F, S, P, BT, N>
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, BT::Agg, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F> + LocalAnalysis<F, BT::Agg, N>,
 {
-    type PreadjointCodomain = BTFN<F, SensorGridSupportGenerator<F, S, P, N>, BT, N>;
+    type PreadjointCodomain = BTFN<F, SensorGridSupportGenerator<S, F, P, N>, BT, N>;
     type Preadjoint<'a>
         = SensorGridPreadjoint<'a, Self, F, N>
     where
@@ -303,10 +306,10 @@
 for SensorGridPreadjoint<'a, SensorGrid<F, S, P, BT, N>, F, N>
 where F : Float,
       BT : SensorGridBT<F, S, P, N>,
-      S : Sensor<F, N>,
-      P : Spread<F, N>,
-      Convolution<S, P> : Spread<F, N> + Lipschitz<L2, FloatType=F> + DifferentiableMapping<Loc<F,N>> + LocalAnalysis<F, BT::Agg, N>,
-      for<'b> <Convolution<S, P> as DifferentiableMapping<Loc<F,N>>>::Differential<'b> : Lipschitz<L2, FloatType=F>,
+      S : Sensor<N, F>,
+      P : Spread<N, F>,
+      Convolution<S, P> : Spread<N, F> + Lipschitz<L2, FloatType=F> + DifferentiableMapping<Loc<N, F>> + LocalAnalysis<F, BT::Agg, N>,
+      for<'b> <Convolution<S, P> as DifferentiableMapping<Loc<N, F>>>::Differential<'b> : Lipschitz<L2, FloatType=F>,
 {
 
     type FloatType = F;
@@ -332,55 +335,49 @@
 */
 
 #[replace_float_literals(F::cast_from(literal))]
-impl<'a, F, S, P, BT, const N: usize> BoundedCurvature for SensorGrid<F, S, P, BT, N>
+impl<'a, F, S, P, BT, const N: usize> BasicCurvatureBoundEstimates<F> for SensorGrid<F, S, P, BT, N>
 where
-    F: Float,
+    F: Float + ToNalgebraRealField,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    DVector<F>: Euclidean<F>,
+    Convolution<S, P>: Spread<N, F>
         + Lipschitz<L2, FloatType = F>
-        + DifferentiableMapping<Loc<F, N>>
+        + DifferentiableMapping<Loc<N, F>>
         + LocalAnalysis<F, BT::Agg, N>,
-    for<'b> <Convolution<S, P> as DifferentiableMapping<Loc<F, N>>>::Differential<'b>:
+    for<'b> <Convolution<S, P> as DifferentiableMapping<Loc<N, F>>>::Differential<'b>:
         Lipschitz<L2, FloatType = F>,
 {
-    type FloatType = F;
-
-    /// Returns factors $ℓ_F$ and $Θ²$ such that
-    /// $B_{F'(μ)} dγ ≤ ℓ_F c_2$ and $⟨F'(μ)+F'(μ+Δ)|Δ⟩ ≤ Θ²|γ|(c_2)‖γ‖$,
-    /// where $Δ=(π_♯^1-π_♯^0)γ$.
-    ///
-    /// See Lemma 3.8, Lemma 5.10, Remark 5.14, and Example 5.15.
-    fn curvature_bound_components(&self) -> (Option<Self::FloatType>, Option<Self::FloatType>) {
+    fn basic_curvature_bound_components(&self) -> (DynResult<F>, DynResult<F>) {
         let n_ψ = self.max_overlapping();
         let ψ_diff_lip = self.base_sensor.diff_ref().lipschitz_factor(L2);
         let ψ_lip = self.base_sensor.lipschitz_factor(L2);
-        let ℓ_F = ψ_diff_lip.map(|l| (2.0 * n_ψ).sqrt() * l);
-        let θ2 = ψ_lip.map(|l| 4.0 * n_ψ * l.powi(2));
+        let ℓ_F0 = ψ_diff_lip.map(|l| (2.0 * n_ψ).sqrt() * l);
+        let Θ2 = ψ_lip.map(|l| 4.0 * n_ψ * l.powi(2));
 
-        (ℓ_F, θ2)
+        (ℓ_F0, Θ2)
     }
 }
 
 #[derive(Clone, Debug)]
-pub struct SensorGridSupportGenerator<F, S, P, const N: usize>
+pub struct SensorGridSupportGenerator<S, F, P, const N: usize>
 where
     F: Float,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
 {
     base_sensor: Convolution<S, P>,
-    grid: LinGrid<F, N>,
+    grid: LinGrid<N, F>,
     weights: DVector<F>,
 }
 
-impl<F, S, P, const N: usize> SensorGridSupportGenerator<F, S, P, N>
+impl<F, S, P, const N: usize> SensorGridSupportGenerator<S, F, P, N>
 where
     F: Float,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
     #[inline]
     fn construct_sensor(&self, id: usize, w: F) -> Weighted<ShiftedSensor<F, S, P, N>, F> {
@@ -397,12 +394,12 @@
     }
 }
 
-impl<F, S, P, const N: usize> SupportGenerator<F, N> for SensorGridSupportGenerator<F, S, P, N>
+impl<F, S, P, const N: usize> SupportGenerator<N, F> for SensorGridSupportGenerator<S, F, P, N>
 where
     F: Float,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
 {
     type Id = usize;
     type SupportType = Weighted<ShiftedSensor<F, S, P, N>, F>;
@@ -434,14 +431,14 @@
     }
 }
 
-impl<F, S, P, BT, const N: usize> ForwardModel<DiscreteMeasure<Loc<F, N>, F>, F>
+impl<F, S, P, BT, const N: usize> ForwardModel<DiscreteMeasure<Loc<N, F>, F>, F>
     for SensorGrid<F, S, P, BT, N>
 where
     F: Float + ToNalgebraRealField<MixedType = F> + nalgebra::RealField,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, BT::Agg, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F> + LocalAnalysis<F, BT::Agg, N>,
 {
     type Observable = DVector<F>;
 
@@ -456,17 +453,17 @@
     }
 }
 
-impl<F, S, P, BT, const N: usize> FindimQuadraticModel<Loc<F, N>, F> for SensorGrid<F, S, P, BT, N>
+impl<F, S, P, BT, const N: usize> FindimQuadraticModel<Loc<N, F>, F> for SensorGrid<F, S, P, BT, N>
 where
     F: Float + ToNalgebraRealField<MixedType = F> + nalgebra::RealField,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, BT::Agg, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F> + LocalAnalysis<F, BT::Agg, N>,
 {
     fn findim_quadratic_model(
         &self,
-        μ: &DiscreteMeasure<Loc<F, N>, F>,
+        μ: &DiscreteMeasure<Loc<N, F>, F>,
         b: &Self::Observable,
     ) -> (DMatrix<F::MixedType>, DVector<F::MixedType>) {
         assert_eq!(b.len(), self.n_sensors());
@@ -483,29 +480,29 @@
     }
 }
 
-/// Implements the calculation a factor $L$ such that $A_*A ≤ L 𝒟$ for $A$ the forward model
+/// Implements the calculation a factor $√L$ such that $A_*A ≤ L 𝒟$ for $A$ the forward model
 /// and $𝒟$ a seminorm of suitable form.
 ///
 /// **This assumes (but does not check) that the sensors are not overlapping.**
 #[replace_float_literals(F::cast_from(literal))]
-impl<F, BT, S, P, K, const N: usize> AdjointProductBoundedBy<RNDM<F, N>, ConvolutionOp<F, K, BT, N>>
-    for SensorGrid<F, S, P, BT, N>
+impl<'a, F, BT, S, P, K, const N: usize>
+    BoundedLinear<RNDM<N, F>, &'a ConvolutionOp<F, K, BT, N>, L2, F> for SensorGrid<F, S, P, BT, N>
 where
-    F: Float + nalgebra::RealField + ToNalgebraRealField,
+    F: Float + ToNalgebraRealField,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N>,
-    K: SimpleConvolutionKernel<F, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F>,
+    K: SimpleConvolutionKernel<N, F>,
     AutoConvolution<P>: BoundedBy<F, K>,
+    Weighted<Shift<K, N, F>, F>: LocalAnalysis<F, BT::Agg, N>,
 {
-    type FloatType = F;
-
-    fn adjoint_product_bound(&self, seminorm: &ConvolutionOp<F, K, BT, N>) -> Option<F> {
+    fn opnorm_bound(&self, seminorm: &'a ConvolutionOp<F, K, BT, N>, _: L2) -> DynResult<F> {
         // Sensors should not take on negative values to allow
         // A_*A to be upper bounded by a simple convolution of `spread`.
+        // TODO: Do we really need this restriction?
         if self.sensor.bounds().lower() < 0.0 {
-            return None;
+            return Err(anyhow!("Sensor not bounded from below by zero"));
         }
 
         // Calculate the factor $L_1$ for betwee $ℱ[ψ * ψ] ≤ L_1 ℱ[ρ]$ for $ψ$ the base spread
@@ -517,33 +514,33 @@
         let l0 = self.sensor.norm(Linfinity) * self.sensor.norm(L1);
 
         // The final transition factor is:
-        Some(l0 * l1)
+        Ok((l0 * l1).sqrt())
     }
 }
 
 macro_rules! make_sensorgridsupportgenerator_scalarop_rhs {
     ($trait:ident, $fn:ident, $trait_assign:ident, $fn_assign:ident) => {
         impl<F, S, P, const N: usize> std::ops::$trait_assign<F>
-            for SensorGridSupportGenerator<F, S, P, N>
+            for SensorGridSupportGenerator<S, F, P, N>
         where
             F: Float,
-            S: Sensor<F, N>,
-            P: Spread<F, N>,
-            Convolution<S, P>: Spread<F, N>,
+            S: Sensor<N, F>,
+            P: Spread<N, F>,
+            Convolution<S, P>: Spread<N, F>,
         {
             fn $fn_assign(&mut self, t: F) {
                 self.weights.$fn_assign(t);
             }
         }
 
-        impl<F, S, P, const N: usize> std::ops::$trait<F> for SensorGridSupportGenerator<F, S, P, N>
+        impl<F, S, P, const N: usize> std::ops::$trait<F> for SensorGridSupportGenerator<S, F, P, N>
         where
             F: Float,
-            S: Sensor<F, N>,
-            P: Spread<F, N>,
-            Convolution<S, P>: Spread<F, N>,
+            S: Sensor<N, F>,
+            P: Spread<N, F>,
+            Convolution<S, P>: Spread<N, F>,
         {
-            type Output = SensorGridSupportGenerator<F, S, P, N>;
+            type Output = SensorGridSupportGenerator<S, F, P, N>;
             fn $fn(mut self, t: F) -> Self::Output {
                 std::ops::$trait_assign::$fn_assign(&mut self.weights, t);
                 self
@@ -551,14 +548,14 @@
         }
 
         impl<'a, F, S, P, const N: usize> std::ops::$trait<F>
-            for &'a SensorGridSupportGenerator<F, S, P, N>
+            for &'a SensorGridSupportGenerator<S, F, P, N>
         where
             F: Float,
-            S: Sensor<F, N>,
-            P: Spread<F, N>,
-            Convolution<S, P>: Spread<F, N>,
+            S: Sensor<N, F>,
+            P: Spread<N, F>,
+            Convolution<S, P>: Spread<N, F>,
         {
-            type Output = SensorGridSupportGenerator<F, S, P, N>;
+            type Output = SensorGridSupportGenerator<S, F, P, N>;
             fn $fn(self, t: F) -> Self::Output {
                 SensorGridSupportGenerator {
                     base_sensor: self.base_sensor.clone(),
@@ -575,14 +572,14 @@
 
 macro_rules! make_sensorgridsupportgenerator_unaryop {
     ($trait:ident, $fn:ident) => {
-        impl<F, S, P, const N: usize> std::ops::$trait for SensorGridSupportGenerator<F, S, P, N>
+        impl<F, S, P, const N: usize> std::ops::$trait for SensorGridSupportGenerator<S, F, P, N>
         where
             F: Float,
-            S: Sensor<F, N>,
-            P: Spread<F, N>,
-            Convolution<S, P>: Spread<F, N>,
+            S: Sensor<N, F>,
+            P: Spread<N, F>,
+            Convolution<S, P>: Spread<N, F>,
         {
-            type Output = SensorGridSupportGenerator<F, S, P, N>;
+            type Output = SensorGridSupportGenerator<S, F, P, N>;
             fn $fn(mut self) -> Self::Output {
                 self.weights = self.weights.$fn();
                 self
@@ -590,14 +587,14 @@
         }
 
         impl<'a, F, S, P, const N: usize> std::ops::$trait
-            for &'a SensorGridSupportGenerator<F, S, P, N>
+            for &'a SensorGridSupportGenerator<S, F, P, N>
         where
             F: Float,
-            S: Sensor<F, N>,
-            P: Spread<F, N>,
-            Convolution<S, P>: Spread<F, N>,
+            S: Sensor<N, F>,
+            P: Spread<N, F>,
+            Convolution<S, P>: Spread<N, F>,
         {
-            type Output = SensorGridSupportGenerator<F, S, P, N>;
+            type Output = SensorGridSupportGenerator<S, F, P, N>;
             fn $fn(self) -> Self::Output {
                 SensorGridSupportGenerator {
                     base_sensor: self.base_sensor.clone(),
@@ -612,13 +609,13 @@
 make_sensorgridsupportgenerator_unaryop!(Neg, neg);
 
 impl<'a, F, S, P, BT, const N: usize> Mapping<DVector<F>>
-    for PreadjointHelper<'a, SensorGrid<F, S, P, BT, N>, RNDM<F, N>>
+    for PreadjointHelper<'a, SensorGrid<F, S, P, BT, N>, RNDM<N, F>>
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, Bounds<F>, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F> + LocalAnalysis<F, Bounds<F>, N>,
 {
     type Codomain = SensorGridBTFN<F, S, P, BT, N>;
 
@@ -634,12 +631,12 @@
 }
 
 impl<'a, F, S, P, BT, const N: usize> Linear<DVector<F>>
-    for PreadjointHelper<'a, SensorGrid<F, S, P, BT, N>, RNDM<F, N>>
+    for PreadjointHelper<'a, SensorGrid<F, S, P, BT, N>, RNDM<N, F>>
 where
     F: Float,
     BT: SensorGridBT<F, S, P, N>,
-    S: Sensor<F, N>,
-    P: Spread<F, N>,
-    Convolution<S, P>: Spread<F, N> + LocalAnalysis<F, Bounds<F>, N>,
+    S: Sensor<N, F>,
+    P: Spread<N, F>,
+    Convolution<S, P>: Spread<N, F> + LocalAnalysis<F, Bounds<F>, N>,
 {
 }

mercurial