- start with asgd, finish with l-bfgs

This commit is contained in:
w.pomp
2026-07-31 15:14:44 +02:00
parent 36b24332cd
commit fdac4e04c4
2 changed files with 47 additions and 35 deletions
+17 -5
View File
@@ -62,8 +62,12 @@ impl RegistrationStep {
/// Default registration steps matching elastix `FixedSmoothingImagePyramid`. /// Default registration steps matching elastix `FixedSmoothingImagePyramid`.
/// ///
/// Uses a multi-resolution approach with L-BFGS optimizer and all pixels on all levels /// Uses a hybrid multi-fidelity approach:
/// for deterministic, precise convergence. /// - Coarse levels (sigma=4,2): ASGD optimizer with all pixels for speed
/// - Fine levels (sigma=1,0.5): L-BFGS optimizer with all pixels for accuracy
///
/// Note: Samples are cached per level (not randomized per iteration like SimpleElastix),
/// but ASGD's sigmoid momentum handles the stochastic gradients effectively.
pub fn default_steps(ndim: usize, n: usize) -> Vec<Self> { pub fn default_steps(ndim: usize, n: usize) -> Vec<Self> {
if ndim == 1 { if ndim == 1 {
return vec![Self { return vec![Self {
@@ -89,8 +93,16 @@ impl RegistrationStep {
.map(|(level, (&s, &d))| { .map(|(level, (&s, &d))| {
let n_pixels = (n / (d * d)).max(4); let n_pixels = (n / (d * d)).max(4);
let is_finest = level == n_levels - 1; let is_finest = level == n_levels - 1;
// Use all pixels on all levels for deterministic results let is_coarse = level < 2;
let (max_iter, tol) = if is_finest { (2048, 1e-8) } else { (512, 1e-6) }; // Hybrid: ASGD on coarse levels for speed, L-BFGS on fine levels for accuracy
// Use all pixels on all levels for consistent gradient estimates
let (optimizer, max_iter, tol) = if is_coarse {
(Optimizer::ASGD, 512, 1e-4)
} else if is_finest {
(Optimizer::LBFGS, 2048, 1e-8)
} else {
(Optimizer::LBFGS, 512, 1e-6)
};
Self { Self {
sigma: Sigma::Absolute(vec![s; ndim]), sigma: Sigma::Absolute(vec![s; ndim]),
samples: SamplingArg::Random(n_pixels), samples: SamplingArg::Random(n_pixels),
@@ -100,7 +112,7 @@ impl RegistrationStep {
max_iterations: max_iter, max_iterations: max_iter,
learning_rate: 1.0, learning_rate: 1.0,
downsample: d, downsample: d,
optimizer: Optimizer::LBFGS, optimizer,
} }
}) })
.collect() .collect()
+30 -30
View File
@@ -1025,22 +1025,22 @@ mod tests {
} }
} }
let mean_diff = sum_diff / count as f64; let mean_diff = sum_diff / count as f64;
println!("Our: {:?} mean_coord_diff: {:.6}", t, mean_diff); println!("Our: {:?}\nmean_coord_diff: {:.6}", t, mean_diff);
let mut tif = IJTiffFile::new( // let mut tif = IJTiffFile::new(
std::env::home_dir() // std::env::home_dir()
.unwrap() // .unwrap()
.join("tmp/register_real_images.tif"), // .join("tmp/register_real_images.tif"),
)?; // )?;
tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?; // tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?;
tif.save( // tif.save(
t.interpolate_par::<1, _, _>(moving.view())? // t.interpolate_par::<3, _, _>(moving.view())?
.mapv(|i| i as u16), // .mapv(|i| i as u16),
1, // 1,
0, // 0,
0, // 0,
)?; // )?;
tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?; // tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?;
assert!(mean_diff < 0.1); assert!(mean_diff < 0.1);
@@ -1079,22 +1079,22 @@ mod tests {
} }
} }
let mean_diff = sum_diff / count as f64; let mean_diff = sum_diff / count as f64;
println!("Our: {:?} mean_coord_diff: {:.6}", t, mean_diff); println!("Our: {:?}\nmean_coord_diff: {:.6}", t, mean_diff);
let mut tif = IJTiffFile::new( // let mut tif = IJTiffFile::new(
std::env::home_dir() // std::env::home_dir()
.unwrap() // .unwrap()
.join("tmp/register_real_images2.tif"), // .join("tmp/register_real_images2.tif"),
)?; // )?;
tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?; // tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?;
tif.save( // tif.save(
t.interpolate_par::<1, _, _>(moving.view())? // t.interpolate_par::<3, _, _>(moving.view())?
.mapv(|i| i as u16), // .mapv(|i| i as u16),
1, // 1,
0, // 0,
0, // 0,
)?; // )?;
tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?; // tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?;
assert!(mean_diff < 0.1); assert!(mean_diff < 0.1);