fpm_rs/metrics/intensity/
comparison.rs1use serde::{Deserialize, Serialize};
4
5#[derive(Clone, Debug, Serialize, Deserialize)]
6pub struct IntensityComparisonMetrics {
7 #[serde(rename = "measured_sum", alias = "reference_sum")]
9 pub reference_sum: f64,
10 #[serde(rename = "predicted_sum", alias = "candidate_sum")]
12 pub candidate_sum: f64,
13 pub residual_l1: f64,
14 pub residual_l2: f64,
15 pub residual_mean: f64,
16 pub residual_std: f64,
17 pub residual_max_abs: f64,
18 pub normalized_l2: f64,
19 pub saturated_pixels: Option<usize>,
20}
21
22pub fn compare_intensity(
23 reference: &[f64],
24 candidate: &[f64],
25 saturation_value: Option<f64>,
26) -> crate::Result<IntensityComparisonMetrics> {
27 compare_intensity_masked(reference, candidate, None, saturation_value)
28}
29
30pub fn compare_intensity_masked(
31 reference: &[f64],
32 candidate: &[f64],
33 mask: Option<&[u8]>,
34 saturation_value: Option<f64>,
35) -> crate::Result<IntensityComparisonMetrics> {
36 if reference.len() != candidate.len() || reference.is_empty() {
37 return Err(crate::Error::InvalidShape(format!(
38 "intensity comparison inputs have lengths {} and {}",
39 reference.len(),
40 candidate.len()
41 )));
42 }
43 if mask.is_some_and(|mask| mask.len() != reference.len()) {
44 return Err(crate::Error::LengthMismatch {
45 actual: mask.map_or(0, <[u8]>::len),
46 expected: reference.len(),
47 shape: (1, reference.len()),
48 });
49 }
50 let mut reference_sum = 0.0;
51 let mut candidate_sum = 0.0;
52 let mut residual_l1 = 0.0;
53 let mut residual_squared = 0.0;
54 let mut residual_sum = 0.0;
55 let mut residual_max_abs: f64 = 0.0;
56 let mut reference_squared = 0.0;
57 let mut saturated_pixels = 0;
58 let mut count = 0usize;
59 for (index, (&reference, &candidate)) in reference.iter().zip(candidate).enumerate() {
60 if mask.is_some_and(|mask| mask[index] == 0) {
61 continue;
62 }
63 if !reference.is_finite() || !candidate.is_finite() {
64 return Err(crate::Error::Numerical(
65 "intensity comparison contains a non-finite value".into(),
66 ));
67 }
68 let residual = candidate - reference;
69 reference_sum += reference;
70 candidate_sum += candidate;
71 residual_l1 += residual.abs();
72 residual_squared += residual * residual;
73 residual_sum += residual;
74 residual_max_abs = residual_max_abs.max(residual.abs());
75 reference_squared += reference * reference;
76 saturated_pixels += usize::from(saturation_value.is_some_and(|limit| reference >= limit));
77 count += 1;
78 }
79 let count_f64 = count as f64;
80 let residual_mean = (count > 0)
81 .then_some(residual_sum / count_f64)
82 .unwrap_or(0.0);
83 Ok(IntensityComparisonMetrics {
84 reference_sum,
85 candidate_sum,
86 residual_l1,
87 residual_l2: residual_squared.sqrt(),
88 residual_mean,
89 residual_std: if count == 0 {
90 0.0
91 } else {
92 (residual_squared / count_f64 - residual_mean * residual_mean)
93 .max(0.0)
94 .sqrt()
95 },
96 residual_max_abs,
97 normalized_l2: if count == 0 {
98 0.0
99 } else {
100 residual_squared.sqrt() / (reference_squared.sqrt() + f64::EPSILON)
101 },
102 saturated_pixels: saturation_value.map(|_| saturated_pixels),
103 })
104}