Skip to main content

fpm_rs/algorithms/
epry.rs

1use crate::{
2    Result,
3    algorithms::objective::LossType,
4    diagnostics::StepDiagnostics,
5    error::Error,
6    measurements::MeasurementRead,
7    reconstruction::{Batch, ReconstructionProblem, ReconstructionState},
8};
9
10use super::{
11    ReconstructionAlgorithm,
12    common::{ObjectDenominator, UpdateConfiguration, projection_update},
13};
14
15/// Embedded pupil-recovery reconstruction for Fourier ptychographic microscopy.
16///
17/// # Method
18///
19/// EPRY alternates detector-amplitude projection with two normalized updates:
20/// the exit-wave error is divided by pupil power to update the overlapping
21/// object-spectrum patch and by object-patch power to update the pupil. Jointly
22/// recovering these two complex functions lets the pupil absorb aberrations
23/// that would otherwise be imprinted on the reconstructed object. The pupil
24/// can be projected back onto the known aperture support after every update.
25///
26/// This implementation also optionally estimates a relative gain and an
27/// additive, spatially uniform background for each frame. Those calibration
28/// updates and incoherent multiplexing support are crate extensions to the
29/// reference EPRY method.
30///
31/// # Reference
32///
33/// X. Ou, G. Zheng, and C. Yang, “Embedded pupil function recovery for Fourier
34/// ptychographic microscopy,” *Optics Express* **22**(5), 4960–4972 (2014),
35/// [doi:10.1364/OE.22.004960](https://doi.org/10.1364/OE.22.004960).
36#[derive(Clone, Debug)]
37pub struct Epry {
38    /// Number of complete passes through the acquisition schedule.
39    pub iterations: usize,
40    /// Relaxation factor applied to each object-spectrum correction.
41    pub object_step: f64,
42    /// Relaxation factor applied to each pupil correction.
43    pub pupil_step: f64,
44    /// Number of measured frames supplied to each reconstruction step.
45    pub batch_size: usize,
46    /// Whether to update the complex pupil alongside the object.
47    pub recover_pupil: bool,
48    /// Whether to zero recovered pupil values outside the compiled aperture.
49    pub constrain_pupil_support: bool,
50    /// Whether to estimate one multiplicative intensity gain per frame.
51    pub recover_frame_gains: bool,
52    /// Fraction of each least-squares frame-gain estimate applied per update.
53    pub gain_step: f64,
54    /// Lower bound for recovered frame gains; must be positive.
55    pub minimum_gain: f64,
56    /// Upper bound for recovered frame gains.
57    pub maximum_gain: f64,
58    /// Whether to estimate one additive, spatially uniform background per frame.
59    pub recover_background: bool,
60    /// Fraction of the mean frame residual added to the background per update.
61    pub background_step: f64,
62    /// Lower bound for recovered background intensities.
63    pub minimum_background: f64,
64    /// Upper bound for recovered background intensities.
65    pub maximum_background: f64,
66    /// Positive numerical floor used in normalized updates.
67    pub epsilon: f64,
68    /// Loss used for diagnostics; the projection itself always enforces the
69    /// measured amplitude.
70    pub loss_type: LossType,
71}
72
73impl Default for Epry {
74    fn default() -> Self {
75        Self {
76            iterations: 100,
77            object_step: 0.8,
78            pupil_step: 0.1,
79            batch_size: 1,
80            recover_pupil: true,
81            constrain_pupil_support: true,
82            recover_frame_gains: false,
83            gain_step: 0.2,
84            minimum_gain: 1e-6,
85            maximum_gain: 1e6,
86            recover_background: false,
87            background_step: 0.2,
88            minimum_background: 0.0,
89            maximum_background: 1e12,
90            epsilon: 1e-10,
91            loss_type: LossType::AmplitudeMse,
92        }
93    }
94}
95
96impl Epry {
97    pub fn iterations(mut self, iterations: usize) -> Self {
98        self.iterations = iterations;
99        self
100    }
101
102    pub fn object_step(mut self, step: f64) -> Self {
103        self.object_step = step;
104        self
105    }
106
107    pub fn pupil_step(mut self, step: f64) -> Self {
108        self.pupil_step = step;
109        self
110    }
111
112    pub fn recover_pupil(mut self, recover: bool) -> Self {
113        self.recover_pupil = recover;
114        self
115    }
116
117    pub fn constrain_pupil_support(mut self, constrain: bool) -> Self {
118        self.constrain_pupil_support = constrain;
119        self
120    }
121
122    pub fn recover_frame_gains(mut self, recover: bool) -> Self {
123        self.recover_frame_gains = recover;
124        self
125    }
126
127    pub fn gain_step(mut self, step: f64) -> Self {
128        self.gain_step = step;
129        self
130    }
131
132    pub fn gain_bounds(mut self, minimum: f64, maximum: f64) -> Self {
133        self.minimum_gain = minimum;
134        self.maximum_gain = maximum;
135        self
136    }
137
138    pub fn recover_background(mut self, recover: bool) -> Self {
139        self.recover_background = recover;
140        self
141    }
142
143    pub fn background_step(mut self, step: f64) -> Self {
144        self.background_step = step;
145        self
146    }
147
148    pub fn background_bounds(mut self, minimum: f64, maximum: f64) -> Self {
149        self.minimum_background = minimum;
150        self.maximum_background = maximum;
151        self
152    }
153}
154
155impl ReconstructionAlgorithm for Epry {
156    fn validate(&self) -> Result<()> {
157        if !self.object_step.is_finite() || self.object_step <= 0.0 {
158            return Err(Error::InvalidParameter {
159                name: "object_step",
160                reason: "must be finite and positive".into(),
161            });
162        }
163        if !self.pupil_step.is_finite() || self.pupil_step < 0.0 {
164            return Err(Error::InvalidParameter {
165                name: "pupil_step",
166                reason: "must be finite and non-negative".into(),
167            });
168        }
169        if !self.epsilon.is_finite() || self.epsilon <= 0.0 || self.batch_size == 0 {
170            return Err(Error::InvalidParameter {
171                name: "epsilon/batch_size",
172                reason: "epsilon must be positive and batch size non-zero".into(),
173            });
174        }
175        if !self.gain_step.is_finite() || !(0.0..=1.0).contains(&self.gain_step) {
176            return Err(Error::InvalidParameter {
177                name: "gain_step",
178                reason: "must be finite and between zero and one".into(),
179            });
180        }
181        if !self.minimum_gain.is_finite()
182            || !self.maximum_gain.is_finite()
183            || self.minimum_gain <= 0.0
184            || self.maximum_gain <= self.minimum_gain
185        {
186            return Err(Error::InvalidParameter {
187                name: "gain_bounds",
188                reason: "must be finite, positive, and strictly increasing".into(),
189            });
190        }
191        if !self.background_step.is_finite() || !(0.0..=1.0).contains(&self.background_step) {
192            return Err(Error::InvalidParameter {
193                name: "background_step",
194                reason: "must be finite and between zero and one".into(),
195            });
196        }
197        if !self.minimum_background.is_finite()
198            || !self.maximum_background.is_finite()
199            || self.maximum_background <= self.minimum_background
200        {
201            return Err(Error::InvalidParameter {
202                name: "background_bounds",
203                reason: "must be finite and strictly increasing".into(),
204            });
205        }
206        Ok(())
207    }
208
209    fn step<M: MeasurementRead>(
210        &mut self,
211        problem: &ReconstructionProblem<M>,
212        state: &mut ReconstructionState,
213        batch: &Batch,
214        _iteration: usize,
215    ) -> Result<StepDiagnostics> {
216        projection_update(
217            problem,
218            state,
219            batch,
220            UpdateConfiguration {
221                object_step: self.object_step,
222                pupil_step: self.recover_pupil.then_some(self.pupil_step),
223                epsilon: self.epsilon,
224                loss_type: self.loss_type,
225                object_denominator: ObjectDenominator::Global,
226                constrain_pupil: self.constrain_pupil_support,
227                gain_update: self.recover_frame_gains.then_some(
228                    super::common::GainUpdateConfiguration {
229                        step: self.gain_step,
230                        minimum: self.minimum_gain,
231                        maximum: self.maximum_gain,
232                    },
233                ),
234                background_update: self.recover_background.then_some(
235                    super::common::BackgroundUpdateConfiguration {
236                        step: self.background_step,
237                        minimum: self.minimum_background,
238                        maximum: self.maximum_background,
239                    },
240                ),
241            },
242        )
243    }
244
245    fn iterations(&self) -> usize {
246        self.iterations
247    }
248
249    fn batch_size(&self) -> usize {
250        self.batch_size
251    }
252}