diff --git a/cpp/radler.cc b/cpp/radler.cc index 78ee214a3bbd0757721e9ca3935952f0894d1e86..7510eb3fc8409de556b55c957a8306d4c45676b0 100644 --- a/cpp/radler.cc +++ b/cpp/radler.cc @@ -154,7 +154,9 @@ void Radler::Perform(bool& reached_major_threshold, double stddev = integrated.StdDevFromMAD(); Logger::Info << "Estimated standard deviation of background noise: " << FluxDensity::ToNiceString(stddev) << '\n'; - if (settings_.auto_mask_sigma && auto_mask_is_finished_) { + const bool auto_mask_is_enabled = + settings_.auto_mask_sigma || settings_.absolute_auto_mask_threshold; + if (auto_mask_is_enabled && auto_mask_is_finished_) { // When we are in the second phase of automasking, don't use // the RMS background anymore parallel_deconvolution_->SetRmsFactorImage(Image()); @@ -208,12 +210,16 @@ void Radler::Perform(bool& reached_major_threshold, parallel_deconvolution_->SetRmsFactorImage(std::move(rms_image)); } } - if (settings_.auto_mask_sigma && !auto_mask_is_finished_) { + 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)); parallel_deconvolution_->SetThreshold( - std::max(stddev * (*settings_.auto_mask_sigma), settings_.threshold)); + 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), settings_.threshold)); + parallel_deconvolution_->SetThreshold( + std::max(stddev * (*settings_.auto_threshold_sigma), + settings_.absolute_threshold)); } integrated.Reset(); @@ -222,7 +228,7 @@ void Radler::Perform(bool& reached_major_threshold, residual_set.LoadAndAveragePsfs(); if (settings_.algorithm_type == AlgorithmType::kMultiscale) { - if (settings_.auto_mask_sigma) { + if (auto_mask_is_enabled) { if (auto_mask_is_finished_) { parallel_deconvolution_->SetAutoMaskMode(false, true); } else { @@ -230,7 +236,7 @@ void Radler::Perform(bool& reached_major_threshold, } } } else { - if (settings_.auto_mask_sigma && auto_mask_is_finished_) { + if (auto_mask_is_enabled && auto_mask_is_finished_) { if (auto_mask_.empty()) { // Generate the auto-mask from the model image(s) auto_mask_.assign(image_width_ * image_height_, false); @@ -251,7 +257,7 @@ void Radler::Perform(bool& reached_major_threshold, residual_set, model_set, psf_images, table_->PsfOffsets(), reached_major_threshold); - if (!reached_major_threshold && settings_.auto_mask_sigma && + if (!reached_major_threshold && auto_mask_is_enabled && !auto_mask_is_finished_) { Logger::Info << "Auto-masking threshold reached; continuing next major " "iteration with deeper threshold and mask.\n"; @@ -339,7 +345,7 @@ void Radler::InitializeDeconvolutionAlgorithm( } algorithm->SetMaxIterations(settings_.minor_iteration_count); - algorithm->SetThreshold(settings_.threshold); + algorithm->SetThreshold(settings_.absolute_threshold); algorithm->SetMinorLoopGain(settings_.minor_loop_gain); algorithm->SetMajorLoopGain(settings_.major_loop_gain); algorithm->SetCleanBorderRatio(settings_.border_ratio); diff --git a/cpp/settings.h b/cpp/settings.h index 7ddf08a0fbda4b8395a7931bed9dfad5acbb0716..1cc6e13643b6a24402ee7f2c218cec9a09f9495f 100644 --- a/cpp/settings.h +++ b/cpp/settings.h @@ -150,7 +150,7 @@ struct Settings { * which means that Radler will keep continuing until another criterion (e.g. * nr. of iterations) is reached. */ - double threshold = 0.0; + double absolute_threshold = 0.0; /** * Gain value for minor loop iterations. @@ -198,6 +198,12 @@ struct Settings { */ std::optional<double> auto_mask_sigma = std::nullopt; + /** + * @brief Like @ref auto_mask_sigma, but instead specifies an absolute + * level where to stop generation of the auto-mask. + */ + std::optional<double> absolute_auto_mask_threshold; + /** * If @c true, maintain a list of components while performing deconvolution. * This works with the @ref AlgorithmType::kGenericClean and @ref diff --git a/cpp/test/test_python_deconvolution.cc b/cpp/test/test_python_deconvolution.cc index 013a91b1d112440abacec191a3f4c54f96e774ef..b2e5ad686051e84c21d0352797145fcae5d878dc 100644 --- a/cpp/test/test_python_deconvolution.cc +++ b/cpp/test/test_python_deconvolution.cc @@ -32,7 +32,7 @@ class PythonFileFixture { settings.pixel_scale.x = kPixelScale; settings.pixel_scale.y = kPixelScale; settings.minor_iteration_count = 1000; - settings.threshold = 1.0e-8; + settings.absolute_threshold = 1.0e-8; settings.algorithm_type = AlgorithmType::kPython; settings.python.filename = kPythonFilename; } diff --git a/cpp/test/test_radler.cc b/cpp/test/test_radler.cc index 9ddd8319589da85c83b471d55193e932aa1c8607..e51bd8d57c85398519107bdd4869d69e5fea623a 100644 --- a/cpp/test/test_radler.cc +++ b/cpp/test/test_radler.cc @@ -86,7 +86,7 @@ struct SettingsFixture { settings.pixel_scale.x = kPixelScale; settings.pixel_scale.y = kPixelScale; settings.minor_iteration_count = 1000; - settings.threshold = 1.0e-8; + settings.absolute_threshold = 1.0e-8; } Settings settings; @@ -191,7 +191,7 @@ BOOST_AUTO_TEST_CASE(diffuse_source) { Settings settings; settings.algorithm_type = AlgorithmType::kMultiscale; - settings.threshold = 1.0e-8; + settings.absolute_threshold = 1.0e-8; settings.major_iteration_count = 30; settings.pixel_scale.x = imgReader.PixelSizeX(); settings.pixel_scale.y = imgReader.PixelSizeY(); @@ -199,7 +199,7 @@ BOOST_AUTO_TEST_CASE(diffuse_source) { settings.trimmed_image_height = imgReader.ImageHeight(); settings.minor_iteration_count = 300; settings.minor_loop_gain = 0.8; - settings.auto_mask_sigma = 4; + settings.auto_mask_sigma = 4.0; const double beamScale = imgReader.BeamMajorAxisRad(); @@ -209,18 +209,19 @@ BOOST_AUTO_TEST_CASE(diffuse_source) { const int major_iteration_count = 0; radler.Perform(reached_threshold, major_iteration_count); - BOOST_REQUIRE(radler.IterationNumber() <= settings.minor_iteration_count); + BOOST_CHECK_LE(radler.IterationNumber(), settings.minor_iteration_count); + BOOST_CHECK_GE(radler.IterationNumber(), 100); const double max_residual = residual_image.Max(); const double rms_residual = residual_image.RMS(); // Checks that RMS value in residual image went down at least by 25%. A // reduction is expected as we remove the peaks from the dirty image. - BOOST_REQUIRE(rms_residual < 0.75 * rms_dirty_image); + BOOST_CHECK_LT(rms_residual, 0.75 * rms_dirty_image); // Checks that the components in the dirty image are reduced in the first - // iteration. We expect thehighest peak to be reduced by 90% in this case. - BOOST_REQUIRE(max_residual < 0.1 * max_dirty_image); + // iteration. We expect the highest peak to be reduced by 90% in this case. + BOOST_CHECK_LT(max_residual, 0.1 * max_dirty_image); } BOOST_AUTO_TEST_SUITE_END() diff --git a/python/docstrings/settings_docstrings.h b/python/docstrings/settings_docstrings.h index 0881342db47503fa0f21043f04ddef4d0d78eab4..f76c0c8883934a5434a9bfc96f42f422e19cbf0f 100644 --- a/python/docstrings/settings_docstrings.h +++ b/python/docstrings/settings_docstrings.h @@ -205,6 +205,17 @@ static const char *__doc_radler_Settings_SpectralFitting_terms = R"doc(Number of spectral terms to constrain the channels to, or zero to disable.)doc"; +static const char *__doc_radler_Settings_absolute_auto_mask_threshold = +R"doc(Like auto_mask_sigma, but instead specifies an absolute level where to +stop generation of the auto-mask.)doc"; + +static const char *__doc_radler_Settings_absolute_threshold = +R"doc(Value in Jy that defines when to stop cleaning. Radler::Perform() will +stop its major iteration and set ``reached_major_threshold``=false +when the peak residual flux is below the given threshold. The default +value is 0.0, which means that Radler will keep continuing until +another criterion (e.g. nr. of iterations) is reached.)doc"; + static const char *__doc_radler_Settings_algorithm_type = R"doc(@{ These deconvolution settings are algorithm-specific. For each algorithm type, a single struct holds all algorithm-specific settings @@ -357,13 +368,6 @@ set ``reached_major_threshold``=false.)doc"; static const char *__doc_radler_Settings_thread_count = R"doc()doc"; -static const char *__doc_radler_Settings_threshold = -R"doc(Value in Jy that defines when to stop cleaning. Radler::Perform() will -stop its major iteration and set ``reached_major_threshold``=false -when the peak residual flux is below the given threshold. The default -value is 0.0, which means that Radler will keep continuing until -another criterion (e.g. nr. of iterations) is reached.)doc"; - static const char *__doc_radler_Settings_trimmed_image_height = R"doc(Trimmed image height)doc"; static const char *__doc_radler_Settings_trimmed_image_width = diff --git a/python/pysettings.cc b/python/pysettings.cc index 35b33be7f1dc6a4f8b65e9060ad094065ed0c65b..56854b9f9974aa30c4ed39b91b75f7a2971e2626 100644 --- a/python/pysettings.cc +++ b/python/pysettings.cc @@ -71,8 +71,9 @@ void init_settings(py::module& m) { DOC(radler_Settings_linked_polarizations)) .def_readwrite("parallel", &radler::Settings::parallel, DOC(radler_Settings_parallel)) - .def_readwrite("threshold", &radler::Settings::threshold, - DOC(radler_Settings_threshold)) + .def_readwrite("absolute_threshold", + &radler::Settings::absolute_threshold, + DOC(radler_Settings_absolute_threshold)) .def_readwrite("minor_loop_gain", &radler::Settings::minor_loop_gain, DOC(radler_Settings_minor_loop_gain)) .def_readwrite("major_loop_gain", &radler::Settings::major_loop_gain, @@ -82,6 +83,9 @@ void init_settings(py::module& m) { DOC(radler_Settings_auto_threshold_sigma)) .def_readwrite("auto_mask_sigma", &radler::Settings::auto_mask_sigma, DOC(radler_Settings_auto_mask_sigma)) + .def_readwrite("absolute_auto_mask_threshold", + &radler::Settings::absolute_auto_mask_threshold, + DOC(radler_Settings_absolute_auto_mask_threshold)) .def_readwrite("save_source_list", &radler::Settings::save_source_list, DOC(radler_Settings_save_source_list)) .def_readwrite("minor_iteration_count", diff --git a/python/test/test_radler.py b/python/test/test_radler.py index 414b95ed651c6174a418175dc83d8c38638ebe39..61d2743144649f4729a131c7fb1e97fb125f0e1a 100644 --- a/python/test/test_radler.py +++ b/python/test/test_radler.py @@ -45,7 +45,7 @@ def get_settings(): settings.pixel_scale.x = PIXEL_SCALE settings.pixel_scale.y = PIXEL_SCALE settings.minor_iteration_count = MINOR_ITERATION_COUNT - settings.threshold = 1e-8 + settings.absolute_threshold = 1e-8 return settings @@ -126,7 +126,7 @@ def test_input_dtype(settings): residual = residual.astype(np.float32) rd.Radler(settings, psf, residual, model, BEAM_SIZE) - model = model.astype(np.int) + model = model.astype(int) with pytest.raises(TypeError): rd.Radler(settings, psf, residual, model, BEAM_SIZE) model = model.astype(np.float32) diff --git a/python/test/test_settings.py b/python/test/test_settings.py index 3eb4e61bd9a52e2186eaea722d90e38a2ae28e63..3d54814eeadd8258b27780a6f4f8f679b23c7d48 100644 --- a/python/test/test_settings.py +++ b/python/test/test_settings.py @@ -25,7 +25,7 @@ def test_layout(): ) assert nested_structs == nested_structs_ref - n_properties_ref = 33 + n_properties_ref = 34 properties = set(filter(lambda x: re.match("^[a-z]+", x), dir(settings))) assert len(properties) == n_properties_ref @@ -45,7 +45,7 @@ def test_default(): assert settings.parallel.grid_width == 1 assert settings.parallel.grid_height == 1 assert settings.parallel.max_threads > 0 - assert settings.threshold == 0.0 + assert settings.absolute_threshold == 0.0 assert settings.minor_loop_gain == 0.1 assert settings.major_loop_gain == 1.0 assert settings.auto_threshold_sigma == None diff --git a/python/test/test_work_table.py b/python/test/test_work_table.py index 5a114dfee3d408e817722507ef30b73630f843e0..a8ca36fd0aaa4d0c2cd935dc3f1e4dc59c4793da 100644 --- a/python/test/test_work_table.py +++ b/python/test/test_work_table.py @@ -166,8 +166,8 @@ def test_add_entries_wrong_type(): entries = [rd.WorkTableEntry(), rd.WorkTableEntry(), rd.WorkTableEntry()] psf = np.ones((4, 4), dtype=np.float64) - residual = np.ones((4, 4), dtype=np.int) - model = np.ones((4, 4), np.complex) + residual = np.ones((4, 4), dtype=int) + model = np.ones((4, 4), complex) with pytest.raises(TypeError): entries[0].psfs.append(psf)