Merge branch 'fix-squared-channels-issue' into 'main'

Draft: Fix squared channel issue

See merge request !152

cpp/algorithms/multiscale_algorithm.cc

@@ -404,6 +404,15 @@ DeconvolutionResult MultiScaleAlgorithm::ExecuteMajorIteration(
         diverging = diverging || std::fabs(*peak_value) >
                                      initial_peak_value * DivergenceLimit();
       }
+      if (!peak_value) {
+        LogReceiver().Error << "Could not continue multi-scale clean, because "
+                               "the sub-minor loop failed to find\n"
+                               "components. This may be caused by combining "
+                               "multi-scale with squared-channel joining.\n"
+                               "It may help to turn off the sub-minor loop "
+                               "optimization with -no-fast-subminor.\n";
+        break;
+      }
 
       SetIterationNumber(subLoop.CurrentIteration());
       scale_infos_[scaleWithPeak].n_components_cleaned +=

cpp/radler.cc

@@ -165,9 +165,11 @@ void Radler::Perform(bool& reached_major_threshold,
   Image integrated(image_width_, image_height_);
   residual_set.GetLinearIntegrated(integrated);
   Logger::Debug << "Calculating standard deviation...\n";
-  double stddev = integrated.StdDevFromMAD();
+  const std::pair<float, float> median_and_stddev =
+      integrated.MedianAndStdDevFromMAD();
+  double stddev = median_and_stddev.second;
   Logger::Info << "Estimated standard deviation of background noise: "
-               << FluxDensity::ToNiceString(stddev) << '\n';
+               << FluxDensity::ToNiceString(median_and_stddev.second) << '\n';
   const bool auto_mask_is_enabled =
       settings_.auto_mask_sigma || settings_.absolute_auto_mask_threshold;
   if (auto_mask_is_enabled && auto_mask_is_finished_) {
@@ -206,15 +208,22 @@ void Radler::Perform(bool& reached_major_threshold,
       parallel_deconvolution_->SetRmsFactorImage(std::move(rms_image));
     }
   }
+
+  // When using squared joins, the stddev is not appropriate as a threshold,
+  // because the whole image will have a positive mean(/median). It is probably
+  // best to always add the mean, but at this point I didn't want to change the
+  // behaviour.
+  const double threshold_bias =
+      settings_.squared_joins ? median_and_stddev.first : 0.0;
   if (auto_mask_is_enabled && !auto_mask_is_finished_) {
-    const double combined_auto_mask_threshold =
-        std::max(stddev * settings_.auto_mask_sigma.value_or(0.0),
-                 settings_.absolute_auto_mask_threshold.value_or(0.0));
+    const double combined_auto_mask_threshold = std::max(
+        stddev * settings_.auto_mask_sigma.value_or(0.0) + threshold_bias,
+        settings_.absolute_auto_mask_threshold.value_or(0.0));
     parallel_deconvolution_->SetThreshold(
         std::max(combined_auto_mask_threshold, settings_.absolute_threshold));
   } else if (settings_.auto_threshold_sigma) {
     parallel_deconvolution_->SetThreshold(
-        std::max(stddev * (*settings_.auto_threshold_sigma),
+        std::max(stddev * (*settings_.auto_threshold_sigma) + threshold_bias,
                  settings_.absolute_threshold));
   }
   integrated.Reset();

aocommon @ 83959a31

-Subproject commit 2363df2242671d1835ca67952c16e3bde4a12148
+Subproject commit 83959a3141e9c040daae7cfe69e6162c34a29700