From fdac4e04c4e68b391382a3a990cf5f7fc760b40d Mon Sep 17 00:00:00 2001 From: "w.pomp" Date: Fri, 31 Jul 2026 15:14:44 +0200 Subject: [PATCH] - start with asgd, finish with l-bfgs --- src/register.rs | 22 ++++++++++++++---- src/transform.rs | 60 ++++++++++++++++++++++++------------------------ 2 files changed, 47 insertions(+), 35 deletions(-) diff --git a/src/register.rs b/src/register.rs index 53a5124..f714ba7 100644 --- a/src/register.rs +++ b/src/register.rs @@ -62,8 +62,12 @@ impl RegistrationStep { /// Default registration steps matching elastix `FixedSmoothingImagePyramid`. /// - /// Uses a multi-resolution approach with L-BFGS optimizer and all pixels on all levels - /// for deterministic, precise convergence. + /// Uses a hybrid multi-fidelity approach: + /// - 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 { if ndim == 1 { return vec![Self { @@ -89,8 +93,16 @@ impl RegistrationStep { .map(|(level, (&s, &d))| { let n_pixels = (n / (d * d)).max(4); let is_finest = level == n_levels - 1; - // Use all pixels on all levels for deterministic results - let (max_iter, tol) = if is_finest { (2048, 1e-8) } else { (512, 1e-6) }; + let is_coarse = level < 2; + // 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 { sigma: Sigma::Absolute(vec![s; ndim]), samples: SamplingArg::Random(n_pixels), @@ -100,7 +112,7 @@ impl RegistrationStep { max_iterations: max_iter, learning_rate: 1.0, downsample: d, - optimizer: Optimizer::LBFGS, + optimizer, } }) .collect() diff --git a/src/transform.rs b/src/transform.rs index d7b1702..e574d10 100644 --- a/src/transform.rs +++ b/src/transform.rs @@ -1025,22 +1025,22 @@ mod tests { } } 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( - std::env::home_dir() - .unwrap() - .join("tmp/register_real_images.tif"), - )?; - tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?; - tif.save( - t.interpolate_par::<1, _, _>(moving.view())? - .mapv(|i| i as u16), - 1, - 0, - 0, - )?; - tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?; + // let mut tif = IJTiffFile::new( + // std::env::home_dir() + // .unwrap() + // .join("tmp/register_real_images.tif"), + // )?; + // tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?; + // tif.save( + // t.interpolate_par::<3, _, _>(moving.view())? + // .mapv(|i| i as u16), + // 1, + // 0, + // 0, + // )?; + // tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?; assert!(mean_diff < 0.1); @@ -1079,22 +1079,22 @@ mod tests { } } 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( - std::env::home_dir() - .unwrap() - .join("tmp/register_real_images2.tif"), - )?; - tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?; - tif.save( - t.interpolate_par::<1, _, _>(moving.view())? - .mapv(|i| i as u16), - 1, - 0, - 0, - )?; - tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?; + // let mut tif = IJTiffFile::new( + // std::env::home_dir() + // .unwrap() + // .join("tmp/register_real_images2.tif"), + // )?; + // tif.save(fixed.mapv(|i| i as u16), 0, 0, 0)?; + // tif.save( + // t.interpolate_par::<3, _, _>(moving.view())? + // .mapv(|i| i as u16), + // 1, + // 0, + // 0, + // )?; + // tif.save(moving.mapv(|i| i as u16), 2, 0, 0)?; assert!(mean_diff < 0.1);