diff --git a/include/schaapcommon/h5parm/jonesparameters.h b/include/schaapcommon/h5parm/jonesparameters.h index 54345038a6810d3b621f414a4b0a1ecdf1812acc..b354915e40b343bf6f701df8fdea33bdaae1b7c6 100644 --- a/include/schaapcommon/h5parm/jonesparameters.h +++ b/include/schaapcommon/h5parm/jonesparameters.h @@ -70,8 +70,6 @@ class JonesParameters { * (e.g. phase and amplitude). These are the parameters as they are * stored (e.g. TEC values). * \param invert (optional default=false) Invert the parameters - * \param sigma_mmse (optional default=0.) Minimum mean square error parameter - * (remnant from BBS times, leave at 0 unless you know what you're doing). */ JonesParameters(const std::vector<double>& freqs, const std::vector<double>& times, @@ -79,7 +77,7 @@ class JonesParameters { GainType gain_type, InterpolationType interpolation_type, hsize_t direction, std::vector<std::vector<std::vector<double>>>&& parm_values, - bool invert = false, float sigma_mmse = 0.); + bool invert = false); /** * Constructor for JonesParameters with given parm_values. To be used if @@ -91,14 +89,12 @@ class JonesParameters { * \param correct_type Correction type of the Jones matrix * \param solution Solution in format [n_ants * n_pols, n_chans] * \param invert (optional default=false) Invert the parameters - * \param sigma_mmse (optional default=0.) Minimum mean square error parameter - * (remnant from BBS times, leave at 0 unless you know what you're doing). */ JonesParameters( const std::vector<double>& freqs, const std::vector<double>& times, const std::vector<std::string>& antenna_names, GainType gain_type, const std::vector<std::vector<std::complex<double>>>& solution, - bool invert = false, float sigma_mmse = 0.); + bool invert = false); /** * Contructor for JonesParameters. JonesParameters will extract parameter @@ -113,8 +109,6 @@ class JonesParameters { * \param sol_tab2 (optional default=nullptr) soltab with parameters for * complex values. Shapes of sol_tab and sol_tab2 can differ * \param invert (optional default=false) Invert the parameters - * \param sigma_mmse (optional default=0.) Minimum mean square error parameter - * (remnant from BBS times, leave at 0 unless you know what you're doing) * \param parm_size (optional default=0) allows to override the vector size * for parm_values * \param missing_antenna_behavior (optional default=kError) what to do with @@ -126,8 +120,7 @@ class JonesParameters { GainType gain_type, InterpolationType interpolation_type, hsize_t direction, schaapcommon::h5parm::SolTab* sol_tab, schaapcommon::h5parm::SolTab* sol_tab2 = nullptr, - bool invert = false, float sigma_mmse = 0., - size_t parm_size = 0, + bool invert = false, size_t parm_size = 0, MissingAntennaBehavior missing_antenna_behavior = MissingAntennaBehavior::kError); @@ -214,12 +207,10 @@ class JonesParameters { * in MakeComplex. * \param parms Reference to the complex parameters that will be inverted * obtained via MakeComplex - * \param sigma_mmse (optional default=0.) Minimum mean square error parameter - * (remnant from BBS times, leave at 0 unless you know what you're doing). * \param correct_type Correction type of the Jones matrix */ static void Invert(casacore::Cube<std::complex<float>>& parms, - float sigma_mmse, GainType gain_type); + GainType gain_type); /// Stored Jones matrices, dimensions (nparms, nantenna, ntime*nfreq), /// frequency varies fastest. nparms is 2 for diagonal, 4 for full jones diff --git a/src/h5parm/jonesparameters.cc b/src/h5parm/jonesparameters.cc index 119c4c8a5077bd21f961771c97069eb31b202fff..eecd0ebe48a5b87d510e0a3d361dfc73eddd22d5 100644 --- a/src/h5parm/jonesparameters.cc +++ b/src/h5parm/jonesparameters.cc @@ -27,23 +27,22 @@ JonesParameters::JonesParameters( const std::vector<double>& freqs, const std::vector<double>& times, const std::vector<std::string>& antenna_names, GainType gain_type, InterpolationType interpolation_type, hsize_t direction, - std::vector<std::vector<std::vector<double>>>&& parm_values, bool invert, - float sigma_mmse) { + std::vector<std::vector<std::vector<double>>>&& parm_values, bool invert) { const size_t num_parms = GetNParms(gain_type); parms_.resize(num_parms, antenna_names.size(), freqs.size() * times.size()); for (size_t ant = 0; ant < antenna_names.size(); ++ant) { MakeComplex(parm_values, ant, freqs, gain_type, invert); } if (invert) { - Invert(parms_, sigma_mmse, gain_type); + Invert(parms_, gain_type); } } JonesParameters::JonesParameters( const std::vector<double>& freqs, const std::vector<double>& times, const std::vector<std::string>& antenna_names, GainType gain_type, - const std::vector<std::vector<std::complex<double>>>& solution, bool invert, - float sigma_mmse) { + const std::vector<std::vector<std::complex<double>>>& solution, + bool invert) { /* * Solution is in format [n_ants * n_pols, n_chans] * Expected format of the parm_values is [n_chans, n_ants, n_pols] @@ -62,7 +61,7 @@ JonesParameters::JonesParameters( } if (invert) { - Invert(parms_, sigma_mmse, gain_type); + Invert(parms_, gain_type); } } @@ -71,8 +70,8 @@ JonesParameters::JonesParameters( const std::vector<std::string>& antenna_names, GainType gain_type, InterpolationType interpolation_type, hsize_t direction, schaapcommon::h5parm::SolTab* sol_tab, - schaapcommon::h5parm::SolTab* sol_tab2, bool invert, float sigma_mmse, - size_t parm_size, MissingAntennaBehavior missing_antenna_behavior) { + schaapcommon::h5parm::SolTab* sol_tab2, bool invert, size_t parm_size, + MissingAntennaBehavior missing_antenna_behavior) { const size_t num_parms = GetNParms(gain_type); parms_.resize(num_parms, antenna_names.size(), freqs.size() * times.size()); @@ -124,7 +123,7 @@ JonesParameters::JonesParameters( } } if (invert) { - Invert(parms_, sigma_mmse, gain_type); + Invert(parms_, gain_type); } } @@ -315,7 +314,7 @@ void JonesParameters::FillParmValues( } void JonesParameters::Invert(casacore::Cube<casacore::Complex>& parms, - float sigma_mmse, GainType gain_type) { + GainType gain_type) { const size_t tf_size = parms.shape()[2]; const size_t ant_size = parms.shape()[1]; for (size_t tf = 0; tf < tf_size; ++tf) { @@ -326,12 +325,6 @@ void JonesParameters::Invert(casacore::Cube<casacore::Complex>& parms, } else if (gain_type == GainType::kFullJones || gain_type == GainType::kFullJonesRealImaginary) { aocommon::MC2x2F v(&parms(0, ant, tf)); - - // Add the variance of the nuisance term to the elements on the - // diagonal. - const float variance = sigma_mmse * sigma_mmse; - v[0] += variance; - v[3] += variance; v.Invert(); v.AssignTo(&parms(0, ant, tf)); } else if (gain_type == GainType::kRotationMeasure || diff --git a/src/h5parm/test/tjonesparameters.cc b/src/h5parm/test/tjonesparameters.cc index 19decd9f41f8716cd2f9dd281a6b7144ef648e3e..496f73251301eaa9bcaec83a8b2185f93c89bc6b 100644 --- a/src/h5parm/test/tjonesparameters.cc +++ b/src/h5parm/test/tjonesparameters.cc @@ -46,7 +46,6 @@ class SolTabMock : public SolTab { }; JonesParameters PrepareJonesParameters(GainType ct, bool invert = false, - float sigma_mmse = 0., size_t parm_size = 0) { std::vector<std::string> antNames; for (size_t i = 0; i < kNAnts; ++i) { @@ -59,8 +58,7 @@ JonesParameters PrepareJonesParameters(GainType ct, bool invert = false, SolTabMock mock2 = SolTabMock(); JonesParameters jones(kFreqs, kTimes, antNames, ct, kInterpolationType, - kDirection, &mock, &mock2, invert, sigma_mmse, - parm_size); + kDirection, &mock, &mock2, invert, parm_size); return jones; } @@ -91,17 +89,16 @@ BOOST_AUTO_TEST_CASE(make_scalar_gain) { } BOOST_AUTO_TEST_CASE(make_complex_fulljones) { - JonesParameters jones = - PrepareJonesParameters(GainType::kFullJones, true, 1.); + JonesParameters jones = PrepareJonesParameters(GainType::kFullJones, true); const auto parms = jones.GetParms(); BOOST_CHECK_EQUAL(parms.shape()[0], 4); - BOOST_CHECK_CLOSE(parms(0, 0, 0).real(), 0.5006091, 1e-3); - BOOST_CHECK_CLOSE(parms(0, 0, 0).imag(), 0.00108991, 1e-3); + BOOST_CHECK_CLOSE(parms(0, 0, 0).real(), 97.437, 1e-3); + BOOST_CHECK_CLOSE(parms(0, 0, 0).imag(), -174.659, 1e-3); } BOOST_AUTO_TEST_CASE(make_complex_tec) { - JonesParameters jones = PrepareJonesParameters(GainType::kTec, false, 0., 1U); + JonesParameters jones = PrepareJonesParameters(GainType::kTec, false, 1U); const auto parms = jones.GetParms(); BOOST_CHECK_EQUAL(parms.shape()[0], 2); @@ -110,7 +107,7 @@ BOOST_AUTO_TEST_CASE(make_complex_tec) { } BOOST_AUTO_TEST_CASE(make_complex_tec2) { - JonesParameters jones = PrepareJonesParameters(GainType::kTec, true, 0., 2U); + JonesParameters jones = PrepareJonesParameters(GainType::kTec, true, 2U); const auto parms = jones.GetParms(); BOOST_CHECK_CLOSE(parms(0, 0, 0).real(), -0.993747, 1e-3); @@ -146,7 +143,7 @@ BOOST_AUTO_TEST_CASE(make_complex_r_measure) { BOOST_AUTO_TEST_CASE(make_complex_r_measure_inverted) { JonesParameters jones = - PrepareJonesParameters(GainType::kRotationMeasure, true, 0.5, 4); + PrepareJonesParameters(GainType::kRotationMeasure, true, 4); const auto parms = jones.GetParms(); BOOST_CHECK_CLOSE(parms(0, 0, 0).real(), -0.185403, 1e-3); @@ -190,8 +187,7 @@ BOOST_AUTO_TEST_CASE(make_scalar_amplitude) { } BOOST_AUTO_TEST_CASE(make_complex_clock) { - JonesParameters jones = - PrepareJonesParameters(GainType::kClock, false, 0., 1U); + JonesParameters jones = PrepareJonesParameters(GainType::kClock, false, 1U); const auto parms = jones.GetParms(); BOOST_CHECK_EQUAL(parms.shape()[0], 2); @@ -200,8 +196,7 @@ BOOST_AUTO_TEST_CASE(make_complex_clock) { } BOOST_AUTO_TEST_CASE(make_complex_clock2) { - JonesParameters jones = - PrepareJonesParameters(GainType::kClock, true, 50, 2U); + JonesParameters jones = PrepareJonesParameters(GainType::kClock, true, 2U); const auto parms = jones.GetParms(); BOOST_CHECK_CLOSE(parms(0, 0, 0).real(), 1., 1e-3); @@ -272,7 +267,7 @@ BOOST_AUTO_TEST_CASE(missing_antenna_flag) { JonesParameters jones(kFreqs, kTimes, antNames, GainType::kDiagonalAmplitude, kInterpolationType, kDirection, &mock, nullptr, false, - 0, 0, JonesParameters::MissingAntennaBehavior::kFlag); + 0, JonesParameters::MissingAntennaBehavior::kFlag); const auto parms = jones.GetParms(); BOOST_CHECK_EQUAL(parms.shape()[1], kNAnts + 1); @@ -292,7 +287,7 @@ BOOST_AUTO_TEST_CASE(missing_antenna_unit_diag) { JonesParameters jones(kFreqs, kTimes, antNames, GainType::kDiagonalAmplitude, kInterpolationType, kDirection, &mock, nullptr, false, - 0, 0, JonesParameters::MissingAntennaBehavior::kUnit); + 0, JonesParameters::MissingAntennaBehavior::kUnit); const auto parms = jones.GetParms(); BOOST_CHECK_EQUAL(parms.shape()[0], 2); @@ -314,7 +309,7 @@ BOOST_AUTO_TEST_CASE(missing_antenna_unit_full) { JonesParameters jones(kFreqs, kTimes, antNames, GainType::kRotationAngle, kInterpolationType, kDirection, &mock, nullptr, false, - 0, 0, JonesParameters::MissingAntennaBehavior::kUnit); + 0, JonesParameters::MissingAntennaBehavior::kUnit); const auto parms = jones.GetParms(); BOOST_CHECK_EQUAL(parms.shape()[0], 4);