1use num_complex::Complex64;
4use serde::{Deserialize, Serialize};
5
6use crate::{
7 Array2, Result,
8 measurements::MeasurementRead,
9 metrics::{
10 complex_field::{ComplexFieldComparisonMetrics, compare_complex_fields_masked},
11 intensity::{IntensityComparisonMetrics, compare_intensity, compare_intensity_masked},
12 model::{
13 FrameGainComparisonMetrics, IlluminationPositionMetrics, PupilComparisonMetrics,
14 compare_frame_gains, compare_illumination_positions, compare_pupils,
15 },
16 },
17 model::{ForwardModel, FourierOffset, ImagePlaneModel},
18 reconstruction::{ReconstructionProblem, ReconstructionResult},
19};
20
21#[derive(Clone, Debug, Serialize, Deserialize)]
22pub struct FrameIntensityEvaluation {
23 pub per_frame: Vec<IntensityComparisonMetrics>,
24}
25
26#[derive(Clone, Debug, Serialize, Deserialize)]
27pub struct ReconstructionEvaluation {
28 pub object: ComplexFieldComparisonMetrics,
29 pub pupil: Option<PupilComparisonMetrics>,
30 pub illumination: Option<IlluminationPositionMetrics>,
31 pub frame_gains: Option<FrameGainComparisonMetrics>,
32 pub intensity: Option<FrameIntensityEvaluation>,
33}
34
35pub fn evaluate_reconstruction(
36 result: &ReconstructionResult,
37 reference_object: &Array2<Complex64>,
38 reference_model: Option<&ImagePlaneModel>,
39 valid_object_mask: Option<&Array2<u8>>,
40) -> Result<ReconstructionEvaluation> {
41 let object =
42 compare_complex_fields_masked(reference_object, &result.object, valid_object_mask)?;
43 let (pupil, illumination, frame_gains) = if let Some(reference_model) = reference_model {
44 let pupil = Some(compare_pupils(
45 &reference_model.pupil.values,
46 &result.recovered_pupil.values,
47 &reference_model.pupil.support,
48 )?);
49 let illumination = None;
50 let frame_gains = result
51 .recovered_frame_gains
52 .as_ref()
53 .map(|candidate| {
54 let reference = reference_model
55 .frame_gains
56 .as_deref()
57 .unwrap_or_else(|| &[]);
58 if reference.is_empty() {
59 compare_frame_gains(&vec![1.0; reference_model.frame_count()], candidate)
60 } else {
61 compare_frame_gains(reference, candidate)
62 }
63 })
64 .transpose()?;
65 (pupil, illumination, frame_gains)
66 } else {
67 (None, None, None)
68 };
69 Ok(ReconstructionEvaluation {
70 object,
71 pupil,
72 illumination,
73 frame_gains,
74 intensity: None,
75 })
76}
77
78pub fn evaluate_reconstruction_with_problem<M: MeasurementRead>(
79 result: &ReconstructionResult,
80 problem: &ReconstructionProblem<M>,
81 reference_object: &Array2<Complex64>,
82 reference_model: Option<&ImagePlaneModel>,
83 valid_object_mask: Option<&Array2<u8>>,
84) -> Result<ReconstructionEvaluation> {
85 let mut evaluation =
86 evaluate_reconstruction(result, reference_object, reference_model, valid_object_mask)?;
87 problem.validate()?;
88 let intensity = evaluate_frame_intensity(result, problem)?;
89 if let Some(reference_model) = reference_model {
90 evaluation.illumination = Some(compare_illumination_positions(
91 reference_model,
92 &problem.model,
93 result.calibrated_illumination.as_deref(),
94 )?);
95 }
96 evaluation.intensity = Some(intensity);
97 Ok(evaluation)
98}
99
100pub fn evaluate_frame_intensity<M: MeasurementRead>(
101 result: &ReconstructionResult,
102 problem: &ReconstructionProblem<M>,
103) -> Result<FrameIntensityEvaluation> {
104 problem.validate()?;
105 let model = model_with_result_calibration(problem, result)?;
106 let forward = ForwardModel::new(&model)?;
107 let mut workspace = forward.workspace();
108 let mut candidate = vec![0.0; model.image_shape.0 * model.image_shape.1];
109 let mut per_frame = Vec::with_capacity(model.frame_count());
110 for frame in 0..model.frame_count() {
111 if problem.measurements.frame_weight(frame)? == 0.0 {
112 per_frame.push(compare_intensity(&candidate, &candidate, None)?);
113 continue;
114 }
115 forward.forward_intensity_into(
116 &result.object_spectrum,
117 &result.recovered_pupil,
118 frame,
119 &mut workspace,
120 &mut candidate,
121 )?;
122 let reference = problem.measurements.frame(frame)?;
123 per_frame.push(compare_intensity_masked(
124 &reference,
125 &candidate,
126 problem.measurements.frame_mask(frame)?,
127 None,
128 )?);
129 }
130 Ok(FrameIntensityEvaluation { per_frame })
131}
132
133fn model_with_result_calibration<M: MeasurementRead>(
134 problem: &ReconstructionProblem<M>,
135 result: &ReconstructionResult,
136) -> Result<ImagePlaneModel> {
137 let mut model = problem.model.clone();
138 model.frame_gains = result.recovered_frame_gains.clone();
139 model.background = result.recovered_background.clone();
140 if let Some(corrections) = &result.calibrated_illumination {
141 model.subpixel_offsets = Some(
142 corrections
143 .iter()
144 .enumerate()
145 .map(|(source, &(row, column))| {
146 let base = problem.model.source_offset(source)?;
147 Ok(FourierOffset::new(base.row + row, base.column + column))
148 })
149 .collect::<Result<Vec<_>>>()?,
150 );
151 }
152 model.validate()?;
153 Ok(model)
154}