fix: revision of the b-spline fitting code

Author: oesteban

sdcflows/interfaces/bspline.py

@@ -56,19 +56,19 @@ class _BSplineApproxInputSpec(BaseInterfaceInputSpec):
     in_data = File(exists=True, mandatory=True, desc="path to a fieldmap")
     in_mask = File(exists=True, desc="path to a brain mask")
     bs_spacing = InputMultiObject(
-        [DEFAULT_ZOOMS_MM],
+        [DEFAULT_HF_ZOOMS_MM],
         traits.Tuple(traits.Float, traits.Float, traits.Float),
         usedefault=True,
         desc="spacing between B-Spline control points",
     )
     ridge_alpha = traits.Float(
-        0.01, usedefault=True, desc="controls the regularization"
+        1e-4, usedefault=True, desc="controls the regularization"
     )
     recenter = traits.Enum(
+        False,
         "mode",
         "median",
         "mean",
-        False,
         usedefault=True,
         desc="strategy to recenter the distribution of the input fieldmap",
     )
@@ -80,7 +80,7 @@ class _BSplineApproxInputSpec(BaseInterfaceInputSpec):
     zooms_min = traits.Union(
         traits.Float,
         traits.Tuple(traits.Float, traits.Float, traits.Float),
-        default_value=4.0,
+        default_value=1.0,
         usedefault=True,
         desc="limit minimum image zooms, set 0.0 to use the original image",
     )
@@ -211,9 +211,7 @@ def _run_interface(self, runtime):
         )
 
         # Fit the model
-        model = lm.Ridge(
-            alpha=self.inputs.ridge_alpha, fit_intercept=False, solver="lsqr"
-        )
+        model = lm.Ridge(alpha=self.inputs.ridge_alpha, fit_intercept=True)
         for attempt in range(3):
             model.fit(colmat, data.reshape(-1))
             extreme = np.abs(model.coef_).max()
@@ -228,6 +226,8 @@ def _run_interface(self, runtime):
                 f"Extreme value {extreme:.2e} detected in spline coefficients."
             )
 
+        LOGGER.info(f"Model fit. Intercept = {model.intercept_}")
+
         # Store coefficients
         index = 0
         self._results["out_coeff"] = []

sdcflows/utils/wrangler.py

@@ -339,6 +339,7 @@ def find_estimators(
             "Dataset includes `B0FieldIdentifier` metadata."
             "Any data missing this metadata will be ignored."
         )
+
         for b0_id in b0_ids:
             # Found B0FieldIdentifier metadata entries
             b0_entities = base_entities.copy()

sdcflows/workflows/fit/fieldmap.py

@@ -86,7 +86,6 @@ def init_fmap_wf(omp_nthreads=1, sloppy=False, debug=False, mode="phasediff", na
     """
     from ...interfaces.bspline import (
         BSplineApprox,
-        DEFAULT_LF_ZOOMS_MM,
         DEFAULT_HF_ZOOMS_MM,
         DEFAULT_ZOOMS_MM,
     )
@@ -115,7 +114,7 @@ def _unzip(fmap_spec):
     bs_filter = pe.Node(BSplineApprox(), name="bs_filter")
     bs_filter.interface._always_run = debug
     bs_filter.inputs.bs_spacing = (
-        [DEFAULT_LF_ZOOMS_MM, DEFAULT_HF_ZOOMS_MM] if not sloppy else [DEFAULT_ZOOMS_MM]
+        [DEFAULT_HF_ZOOMS_MM] if not sloppy else [DEFAULT_ZOOMS_MM]
     )
     bs_filter.inputs.extrapolate = not debug