sim_lib_interference_core/
certificate.rs

1//! Sampling adequacy measurements carried with every field request.
2
3use crate::{Emitter, InterferenceError, InterferenceProblem, SamplingPlane};
4
5/// Explicit thresholds used to classify a [`SamplingCertificate`].
6///
7/// The default resolved carrier requirement is eight samples per wavelength,
8/// which gives four samples across the worst-case half-wavelength power
9/// fringe. Four carrier samples (two per power fringe) is the marginal Nyquist
10/// floor. Point-source `1/r` envelope change defaults to five percent for
11/// resolved and ten percent for marginal.
12#[derive(Clone, Copy, Debug, PartialEq)]
13pub struct SamplingThresholds {
14    /// Comfortable carrier samples per wavelength on both axes.
15    pub resolved_min_samples_per_wavelength: f64,
16    /// Marginal carrier samples per wavelength on both axes.
17    pub marginal_min_samples_per_wavelength: f64,
18    /// Comfortable maximum fractional `1/r` change across a cell.
19    pub resolved_max_envelope_fraction_per_cell: f64,
20    /// Marginal maximum fractional `1/r` change across a cell.
21    pub marginal_max_envelope_fraction_per_cell: f64,
22}
23
24impl SamplingThresholds {
25    /// Validates explicit resolved and marginal classification thresholds.
26    pub fn new(
27        resolved_min_samples_per_wavelength: f64,
28        marginal_min_samples_per_wavelength: f64,
29        resolved_max_envelope_fraction_per_cell: f64,
30        marginal_max_envelope_fraction_per_cell: f64,
31    ) -> Result<Self, InterferenceError> {
32        Self {
33            resolved_min_samples_per_wavelength,
34            marginal_min_samples_per_wavelength,
35            resolved_max_envelope_fraction_per_cell,
36            marginal_max_envelope_fraction_per_cell,
37        }
38        .validate()
39    }
40
41    /// Revalidates a threshold record, including one built with field syntax.
42    pub fn validate(self) -> Result<Self, InterferenceError> {
43        positive_threshold(
44            "resolved-min-samples-per-wavelength",
45            self.resolved_min_samples_per_wavelength,
46        )?;
47        positive_threshold(
48            "marginal-min-samples-per-wavelength",
49            self.marginal_min_samples_per_wavelength,
50        )?;
51        non_negative_threshold(
52            "resolved-max-envelope-fraction-per-cell",
53            self.resolved_max_envelope_fraction_per_cell,
54        )?;
55        non_negative_threshold(
56            "marginal-max-envelope-fraction-per-cell",
57            self.marginal_max_envelope_fraction_per_cell,
58        )?;
59        if self.resolved_min_samples_per_wavelength < self.marginal_min_samples_per_wavelength
60            || self.resolved_max_envelope_fraction_per_cell
61                > self.marginal_max_envelope_fraction_per_cell
62        {
63            return Err(InterferenceError::InconsistentSamplingThresholds {
64                resolved_min_samples_per_wavelength: self.resolved_min_samples_per_wavelength,
65                marginal_min_samples_per_wavelength: self.marginal_min_samples_per_wavelength,
66                resolved_max_envelope_fraction_per_cell: self
67                    .resolved_max_envelope_fraction_per_cell,
68                marginal_max_envelope_fraction_per_cell: self
69                    .marginal_max_envelope_fraction_per_cell,
70            });
71        }
72        Ok(self)
73    }
74}
75
76impl Default for SamplingThresholds {
77    fn default() -> Self {
78        Self {
79            resolved_min_samples_per_wavelength: 8.0,
80            marginal_min_samples_per_wavelength: 4.0,
81            resolved_max_envelope_fraction_per_cell: 0.05,
82            marginal_max_envelope_fraction_per_cell: 0.10,
83        }
84    }
85}
86
87/// Whether a non-resolved request is refused or explicitly annotated.
88#[derive(Clone, Copy, Debug, Default, PartialEq, Eq)]
89pub enum SamplingPolicy {
90    /// Require a fully resolved certificate before solving.
91    #[default]
92    Strict,
93    /// Admit the request while preserving its non-resolved certificate.
94    Annotate,
95}
96
97impl SamplingPolicy {
98    /// Applies this policy without allocating field storage.
99    pub fn admit(self, certificate: &SamplingCertificate) -> Result<(), InterferenceError> {
100        if self == Self::Strict && certificate.verdict != SamplingVerdict::Resolved {
101            Err(InterferenceError::SamplingRefused {
102                certificate: *certificate,
103            })
104        } else {
105            Ok(())
106        }
107    }
108}
109
110/// Physical sampling classification for both phase and point-source envelope.
111#[derive(Clone, Copy, Debug, PartialEq, Eq)]
112pub enum SamplingVerdict {
113    /// Meets the comfortable carrier/fringe and envelope requirements.
114    Resolved,
115    /// Meets the minimum requirements but not the comfortable requirements.
116    Marginal,
117    /// Misses carrier/fringe Nyquist or the marginal envelope requirement.
118    Aliased,
119}
120
121/// Measurements explaining whether a problem is resolved by a sampling plane.
122///
123/// Coherent squared magnitude can contain a worst-case fringe period of
124/// `wavelength / 2`, so its samples-per-period values are half the carrier
125/// samples-per-wavelength values. For point sources, the envelope bound covers
126/// one complete cell diagonal relative to the nearest point on the finite
127/// plane. A source touching that plane is rejected because no finite `1/r`
128/// envelope certificate exists.
129#[derive(Clone, Copy, Debug, PartialEq)]
130pub struct SamplingCertificate {
131    /// Explicit thresholds used for this classification.
132    pub thresholds: SamplingThresholds,
133    /// Unattenuated carrier wavelength.
134    pub wavelength_m: f64,
135    /// Carrier samples per wavelength along the plane's `u` axis.
136    pub samples_per_wavelength_u: f64,
137    /// Carrier samples per wavelength along the plane's `v` axis.
138    pub samples_per_wavelength_v: f64,
139    /// Samples per worst-case squared-magnitude fringe along `u`.
140    pub samples_per_power_fringe_u: f64,
141    /// Samples per worst-case squared-magnitude fringe along `v`.
142    pub samples_per_power_fringe_v: f64,
143    /// Distance from a point source to the nearest point on the finite plane.
144    ///
145    /// This is `None` when the problem has no point source.
146    pub nearest_point_source_distance_m: Option<f64>,
147    /// Conservative maximum fractional `1/r` change across one cell.
148    pub max_envelope_fraction_per_cell: f64,
149    /// Classification obtained from every preceding measurement.
150    pub verdict: SamplingVerdict,
151}
152
153impl SamplingCertificate {
154    /// Measures phase, worst-case power fringes, and point-source envelope.
155    ///
156    /// This performs only constant storage work and iterates over sources; it
157    /// does not allocate a grid.
158    pub fn measure(
159        problem: &InterferenceProblem,
160        plane: &SamplingPlane,
161    ) -> Result<Self, InterferenceError> {
162        Self::measure_with_thresholds(problem, plane, SamplingThresholds::default())
163    }
164
165    /// Measures using caller-supplied, validated classification thresholds.
166    pub fn measure_with_thresholds(
167        problem: &InterferenceProblem,
168        plane: &SamplingPlane,
169        thresholds: SamplingThresholds,
170    ) -> Result<Self, InterferenceError> {
171        let thresholds = thresholds.validate()?;
172        let wavelength_m = finite_metric("wavelength-m", problem.wavelength_metres())?;
173        let samples_per_wavelength_u = finite_metric(
174            "samples-per-wavelength-u",
175            wavelength_m / plane.cell_size_u_m(),
176        )?;
177        let samples_per_wavelength_v = finite_metric(
178            "samples-per-wavelength-v",
179            wavelength_m / plane.cell_size_v_m(),
180        )?;
181        let samples_per_power_fringe_u = samples_per_wavelength_u / 2.0;
182        let samples_per_power_fringe_v = samples_per_wavelength_v / 2.0;
183        let nearest_point_source_distance_m = nearest_point_source_distance(problem, plane)?;
184        let max_envelope_fraction_per_cell = match nearest_point_source_distance_m {
185            Some(distance) => conservative_envelope_change(distance, plane)?,
186            None => 0.0,
187        };
188
189        let minimum_carrier_samples = samples_per_wavelength_u.min(samples_per_wavelength_v);
190        let verdict = if minimum_carrier_samples >= thresholds.resolved_min_samples_per_wavelength
191            && max_envelope_fraction_per_cell <= thresholds.resolved_max_envelope_fraction_per_cell
192        {
193            SamplingVerdict::Resolved
194        } else if minimum_carrier_samples >= thresholds.marginal_min_samples_per_wavelength
195            && max_envelope_fraction_per_cell <= thresholds.marginal_max_envelope_fraction_per_cell
196        {
197            SamplingVerdict::Marginal
198        } else {
199            SamplingVerdict::Aliased
200        };
201
202        Ok(Self {
203            thresholds,
204            wavelength_m,
205            samples_per_wavelength_u,
206            samples_per_wavelength_v,
207            samples_per_power_fringe_u,
208            samples_per_power_fringe_v,
209            nearest_point_source_distance_m,
210            max_envelope_fraction_per_cell,
211            verdict,
212        })
213    }
214}
215
216fn positive_threshold(name: &'static str, value: f64) -> Result<(), InterferenceError> {
217    if value.is_finite() && value > 0.0 {
218        Ok(())
219    } else {
220        Err(InterferenceError::InvalidSamplingThreshold { name, value })
221    }
222}
223
224fn non_negative_threshold(name: &'static str, value: f64) -> Result<(), InterferenceError> {
225    if value.is_finite() && value >= 0.0 {
226        Ok(())
227    } else {
228        Err(InterferenceError::InvalidSamplingThreshold { name, value })
229    }
230}
231
232fn finite_metric(name: &'static str, value: f64) -> Result<f64, InterferenceError> {
233    value
234        .is_finite()
235        .then_some(value)
236        .ok_or(InterferenceError::NonFiniteSamplingMetric { name, value })
237}
238
239fn nearest_point_source_distance(
240    problem: &InterferenceProblem,
241    plane: &SamplingPlane,
242) -> Result<Option<f64>, InterferenceError> {
243    let mut nearest: Option<f64> = None;
244    for source in &problem.sources {
245        let Emitter::Point { position, .. } = source else {
246            continue;
247        };
248        let [origin_x, origin_y, origin_z] = plane.origin().coordinates_metres();
249        let [source_x, source_y, source_z] = position.coordinates_metres();
250        let displacement = [
251            source_x - origin_x,
252            source_y - origin_y,
253            source_z - origin_z,
254        ];
255        if displacement.iter().any(|component| !component.is_finite()) {
256            return Err(InterferenceError::NonFiniteSamplingMetric {
257                name: "point-source-plane-displacement-m",
258                value: displacement
259                    .into_iter()
260                    .find(|component| !component.is_finite())
261                    .unwrap_or(f64::NAN),
262            });
263        }
264        let [ux, uy, uz] = plane.u_axis().components();
265        let [vx, vy, vz] = plane.v_axis().components();
266        let [nx, ny, nz] = plane.normal().components();
267        let projected_u = displacement[0] * ux + displacement[1] * uy + displacement[2] * uz;
268        let projected_v = displacement[0] * vx + displacement[1] * vy + displacement[2] * vz;
269        let projected_normal = displacement[0] * nx + displacement[1] * ny + displacement[2] * nz;
270        let outside_u = outside_extent_distance(projected_u, plane.extent_u().get());
271        let outside_v = outside_extent_distance(projected_v, plane.extent_v().get());
272        let distance = finite_metric(
273            "nearest-point-source-distance-m",
274            outside_u.hypot(outside_v).hypot(projected_normal),
275        )?;
276        nearest = Some(nearest.map_or(distance, |current| current.min(distance)));
277    }
278    Ok(nearest)
279}
280
281fn outside_extent_distance(coordinate: f64, extent: f64) -> f64 {
282    if coordinate < 0.0 {
283        -coordinate
284    } else if coordinate > extent {
285        coordinate - extent
286    } else {
287        0.0
288    }
289}
290
291fn conservative_envelope_change(
292    distance: f64,
293    plane: &SamplingPlane,
294) -> Result<f64, InterferenceError> {
295    let cell_diagonal = finite_metric(
296        "sampling-cell-diagonal-m",
297        plane.cell_size_u_m().hypot(plane.cell_size_v_m()),
298    )?;
299    let value = cell_diagonal / distance;
300    if distance == 0.0 || !value.is_finite() {
301        Err(InterferenceError::UnboundedSamplingEnvelope {
302            nearest_point_source_distance_m: distance,
303            cell_diagonal_m: cell_diagonal,
304        })
305    } else {
306        Ok(value)
307    }
308}
sim_lib_interference_core/certificate.rs

sim_lib_interference_core/
certificate.rs
