diff --git a/indica/operators/extrapolate_impurity_density.py b/indica/operators/extrapolate_impurity_density.py
index ec85723e..772af426 100644
--- a/indica/operators/extrapolate_impurity_density.py
+++ b/indica/operators/extrapolate_impurity_density.py
@@ -837,7 +837,9 @@ def optimize_perturbation(
             orig_bolometry = orig_bolometry_data
 
         extrapolated_smooth_data_mean = np.mean(
-            extrapolated_smooth_data.loc[threshold_rho[0] :, :, :]
+            extrapolated_smooth_data.sel(
+                rho_poloidal=slice(threshold_rho.isel(t=0), None)
+            )
         )
 
         def objective_func(objective_array: Sequence, time: float):