diff --git a/hera_cal/abscal.py b/hera_cal/abscal.py index af78f3d76..080eb25ad 100644 --- a/hera_cal/abscal.py +++ b/hera_cal/abscal.py @@ -992,6 +992,7 @@ def RFI_delay_slope_cal(reds, antpos, red_data, freqs, rfi_chans, rfi_headings, gains[ant] = np.exp(2j * np.pi * np.outer(dlys, freqs)) return gains + def cross_pol_phase_cal(model, data, model_bls, data_bls, wgts={}, return_gains=False, gain_ants=[], refpol="Jee"): """ Solve for the relative phase degeneracy between two polarizations of redundantly calibrated data. @@ -1009,7 +1010,7 @@ def cross_pol_phase_cal(model, data, model_bls, data_bls, wgts={}, return_gains= wgts : dict Dictionary mapping baseline tuples to weights. Default is no weights. return_gains : bool - If True, return gains instead of the relative phase difference between the two polarizations. + If True, return gains instead of the relative phase difference between the two polarizations. Default is False. gain_ants : list of tuples List of antenna-pol tuples to return gains for. Default is None. @@ -1021,7 +1022,7 @@ def cross_pol_phase_cal(model, data, model_bls, data_bls, wgts={}, return_gains= -------- If return_gains is False: delta : np.ndarray - Array of relative phase differences between the two polarizations in units of radians. + Array of relative phase differences between the two polarizations in units of radians. If return_gains is True: delta_gains : dict Dictionary mapping antenna keys like (0, 'Jee') to gains of the same shape of the data @@ -1050,26 +1051,27 @@ def cross_pol_phase_cal(model, data, model_bls, data_bls, wgts={}, return_gains= if refpol == data_pol1: summation += ( - data[data_bl] * np.conj(model[model_bl]) * wgts.get(data_bl, 1.0) + np.conj(data[data_bl]) * model[model_bl] * wgts.get(data_bl, 1.0) ) else: summation += ( - np.conj(data[data_bl]) * model[model_bl] * wgts.get(data_bl, 1.0) + data[data_bl] * np.conj(model[model_bl]) * wgts.get(data_bl, 1.0) ) if return_gains: delta_gains = {} for (ant, pol) in gain_ants: if pol == refpol: - delta_gains[(ant, pol)] = np.exp(1j * np.angle(summation)) - else: delta_gains[(ant, pol)] = np.ones_like(summation) + else: + delta_gains[(ant, pol)] = np.exp(1j * np.angle(summation)) return delta_gains - + else: return np.angle(summation) + def dft_phase_slope_solver(xs, ys, data, flags=None): '''Solve for spatial phase slopes across an array by looking for the peak in the DFT. This is analogous to the method in utils.fft_dly(), except its in 2D and does not diff --git a/hera_cal/lst_stack/calibration.py b/hera_cal/lst_stack/calibration.py index 270f8a415..8b73674c0 100644 --- a/hera_cal/lst_stack/calibration.py +++ b/hera_cal/lst_stack/calibration.py @@ -109,7 +109,7 @@ def _expand_degeneracies_to_ant_gains( # Map antenna-polarizations to visibility indices antpol_to_idx = { utils.split_pol(pol)[0]: polidx - for polidx, pol in enumerate(stack.pols) + for polidx, pol in enumerate(stack.pols) if utils.split_pol(pol)[0] == utils.split_pol(pol)[1] } @@ -157,8 +157,6 @@ def _expand_degeneracies_to_ant_gains( per_day_flags, 1.0 + 0.0j, smooth_ant_gain * rephasor.conj() ) - - return gains @@ -580,7 +578,7 @@ def lstbin_absolute_calibration( ) if run_cross_pol_phase_cal and not all_nights_flagged: - cross_pols_in_data = any(utils.split_pol(pol)[0] != utils.split_pol(pol)[1] for pol in stack.pols) + cross_pols_in_data = any(utils.split_pol(pol)[0] != utils.split_pol(pol)[1] for pol in stack.pols) if cross_pols_in_data: delta = _lstbin_cross_pol_phase_calibration(stack=stack, model=model) @@ -595,7 +593,7 @@ def lstbin_absolute_calibration( for ant in gain_ants: for pol in unique_pols: - if pol != refpol: + if pol == refpol: phase_gains[(ant, pol)] = phase_gains.get( (ant, pol), np.ones_like(delta) ) diff --git a/hera_cal/lst_stack/tests/test_calibration.py b/hera_cal/lst_stack/tests/test_calibration.py index d1879a764..d494d7541 100644 --- a/hera_cal/lst_stack/tests/test_calibration.py +++ b/hera_cal/lst_stack/tests/test_calibration.py @@ -306,8 +306,8 @@ def test_relative_phase_calibration(self): [ np.ones((20, 10, 100)), np.ones((20, 10, 100)), - np.exp(1j * delta), np.exp(-1j * delta), + np.exp(1j * delta), ] ) gains = np.transpose(gains, (1, 2, 3, 0)) diff --git a/hera_cal/tests/test_abscal.py b/hera_cal/tests/test_abscal.py index b127e4dd8..f55f8764f 100644 --- a/hera_cal/tests/test_abscal.py +++ b/hera_cal/tests/test_abscal.py @@ -755,7 +755,7 @@ def test_complex_phase_abscal(self): def test_cross_pol_phase_cal(self): rng = np.random.default_rng(42) antpos = hex_array(3, split_core=True, outriggers=0) - + model, data = {}, {} ntimes, nfreqs = 1, 2 @@ -775,21 +775,22 @@ def test_cross_pol_phase_cal(self): data[bls] = model[bls] * gain # Solve for the phase degeneracy - solved_delta = abscal.cross_pol_phase_cal(model, data, model_bls, data_bls) - + solved_delta = abscal.cross_pol_phase_cal(model, data, model_bls, data_bls, refpol='Jnn') + # Check that the phase degeneracy was solved for correctly np.testing.assert_array_almost_equal(solved_delta, delta, decimal=5) gain_ants = set() for bl in data_bls: gain_ants.update(set(utils.split_bl(bl))) gain_ants = list(gain_ants) - gains = abscal.cross_pol_phase_cal(model, data, model_bls, data_bls, return_gains=True, gain_ants=gain_ants) + gains = abscal.cross_pol_phase_cal(model, data, model_bls, data_bls, return_gains=True, gain_ants=gain_ants, refpol='Jnn') for k in gains: # Check that the gains are correct if k[-1] == 'Jnn': np.testing.assert_array_almost_equal(gains[k], 1.0 + 0.0j, decimal=5) + @pytest.mark.filterwarnings("ignore:The default for the `center` keyword has changed") @pytest.mark.filterwarnings("ignore:invalid value encountered in true_divide") @pytest.mark.filterwarnings("ignore:divide by zero encountered in true_divide")