Bug Fixes for eht-imaging for 2021 project (1 minor and 1 major)

[ptiede]

Hi eht-imaging team,

I and others (primarily others) found a few bugs in eht-imaging that this pull-request attempts to fix.

The first is relatively minor and adds a check to read the MJD optionally from the uvfits DATE-OBS entry in the header. This is used for some CASA pipelines that modify the PZERO5 time offset so the mjd is scrambled.

The second bug fix is more critical and is an error in eht-imaging and could lead to the wrong errors getting assigned to the cross-hand polarization products. When the eht does baseline re-ordering so that everything is consistent, it currently does not also transpose the error matrix. This means that the error for RL and LR gets swapped. You can see this most easily in data that has a missing feed. For example, consider the baseline JCMT-GLT where JCMT has only the RCP product. Then the coherency and error look like

 coherency = [[a+ib c+id]
               [NaN NaN]]

error =  [ [e  f],
         [NaN NaN]

Now if we use obs.reorder_baselines() to swap the ordering to GLT-JCMT, the result looks like

 coherency = [[a-ib NaN]
            [c-ib NaN]]

error = [ [e   f],
         [NaN NaN]

I found three instances where this occurred in the code, but I may have missed more.

[achael]

Thanks for catching these!

For the mjd issue, can you change it so the current behavior is the default? I vaguely remember long ago there being issues/inconsistencies in some of the test data sets in the 'DATE-OBS' header parameter w/r/t to the data 'DATE' field. (or at least, I am worried such inconsistencies could possibly exist).

For the sigma issue, note that unless you load in data natively with polrep='circ' the RL sigma and LR sigma will always be equal anyway via the conversion to Stokes QU. Other places we should check are all consistent/correct are obsdata.switch_polrep, caltable.applycal, and various functions in obs_simulate & obs_helpers.

[ptiede]

Hi Andrew,

For the MJD stuff, the problem is I don't know of a way to automatically check whether the DATE field is "wrong". At least for the current CASA data, we end up in a situation where the DATE field is populated, but its value isn't compatible with what is assumed during the conversion to MJD. If I make the DATE field the default this won't really fix the issue unless I add an option to the load uvfits functions. If that is what you were thinking I can definitely do that.

For the sigma issue, note that unless you load in data natively with polrep='circ' the RL sigma and LR sigma will always be equal anyway via the conversion to Stokes QU. Other places we should check are all consistent/correct are obsdata.switch_polrep, caltable.applycal, and various functions in obs_simulate & obs_helpers.

I didn't realize that the errors for switch_polrep are slightly wrong. Given that obs.switch_polrep silently mishandles the errors, maybe that function should just be removed from the user API or at least flash a warning?

I also checked the other functions

1. obs.switch_polrep I don't believe I can fix this without a rewrite of how the stokes quantities are handled. Because Q, U are combinations of the same measured quantities they are correlated, but ehtim doesn't track this. If I want to go from the stokes errors back to the correlation coefficients I would need to know what this correlation is.
2. caltable.applycal: Looks fine to me but you may want to check it.
3. obs_simulate the only funky thing I saw here is in apply_jones_inverse where it appears the final error is identical for all correlation products. I am guessing this is because you don't track the correlations between the quantities?
4. obs_helpers: I don't think I saw anything but you should probably check it out.

[achael]

Hi Paul,

Sorry for the late follow up, the last month has been too busy (please send me emails if this is urgent).

Regarding the date, my concern is that the 'time' fields ehtim currently populates by subtracting off the integer MJD from the 'DATE' field will be inconsistent with the header date. I believe this is why we populate the MJD from the 'DATE' field instead of the 'DATE-OBS' field. I think in some of the BU datasets the 'DATE-OBS' was off from what would be calculated via the 'DATE' value. In any case, if the 'DATE' value is wrong, all of the time data is wrong too -- the fact that the CASA files have incorrect values in the 'DATE' field seems like a problem with the CASA pipeline that should be fixed before loading into ehtim, or with an optional flag to load.py for CASA data alone that doesn't break the default behavior for everything else.

Regarding neglecting the correlations in switching between IQUV -- RR,LL,LR,RL: is this behavior any different from what AIPS/Difmap/CASA/the uvfist standard themselves do? Since data in the uvfits standard (https://library.nrao.edu/public/memos/aips/memos/AIPSM_117.pdf) can be in any of the 3 polarization bases and the uvfits weight parameters do not carry around the full correlation matrix, I don't see how you can guarantee the full information hasn't been lost along the way if the data changed bases before being saved to uvfits. Are the differences between e.g. sigma_RR and sigma_LL or sigma_RL and sigma_LR ever significantly nonzero for EHT data? Certainly they are all zero in every synthetic data set we have ever generated.

If people in your team were interested i'd be happy to have a brief hackathon on improving data tables in ehtim to store the full polarization correlations. It would also be good to allow ehtim to read in linear feed & multi-frequency data. These projects have been on the backburner for me for a while, but if there were interest a dedicated group could probably accomplish them in a day or so.

[ptiede]

Hi Andrew,

Regarding the date, my concern is that the 'time' fields ehtim currently populates by subtracting off the integer MJD from the 'DATE' field will be inconsistent with the header date. I believe this is why we populate the MJD from the 'DATE' field instead of the 'DATE-OBS' field. I think in some of the BU datasets the 'DATE-OBS' was off from what would be calculated via the 'DATE' value. In any case, if the 'DATE' value is wrong, all of the time data is wrong too -- the fact that the CASA files have incorrect values in the 'DATE' field seems like a problem with the CASA pipeline that should be fixed before loading into ehtim, or with an optional flag to load.py for CASA data alone that doesn't break the default behavior for everything else.

I agree this is a good point. I'll raise this with the CASA team and ask them if they can fix how they output the date obs.

Regarding neglecting the correlations in switching between IQUV -- RR,LL,LR,RL: is this behavior any different from what AIPS/Difmap/CASA/the uvfist standard themselves do? Since data in the uvfits standard (https://library.nrao.edu/public/memos/aips/memos/AIPSM_117.pdf) can be in any of the 3 polarization bases and the uvfits weight parameters do not carry around the full correlation matrix, I don't see how you can guarantee the full information hasn't been lost along the way if the data changed bases before being saved to uvfits.

I agree, and I think this is a larger issue with uvfits, so I probably shouldn't have brought it up. We can leave it to a later PR.

Are the differences between e.g. sigma_RR and sigma_LL or sigma_RL and sigma_LR ever significantly nonzero for EHT data? Certainly they are all zero in every synthetic data set we have ever generated.

There are differences, but I agree they are small, and the issue only arises if the R and L feeds have significantly different SEFDs. The problem we had came up when looking at JCMT, which only measured a single polarization. In that case, both cross-hands got NaN'd because of the baseline ordering switch. This PR does fix that so I think we are good to go there, but the more general issue of how to deal with errors after switching polarization bases is more difficult and should probably be tackled in a later PR.

As a plan how about I roll back the date change until we get feedback from the CASA team, and keep the change I made to the error matrix? The other changes, e..g, fixing basis changes, should be left to a later PR.

[ptiede]

Hi Andrew,

I quickly checked how much the errors differ between the cross-hand ratios. Ignoring stations with missing feeds that can be quite a bit worse, I am seeing the following

So, while the differences are small, they aren't zero. At least the part of the PR that fixes this issue should be merged pretty quickly. Should I open another PR just focussing on that problem?

[achael]

Hi @ptiede -- yes, opening another pull request to focus just on this issue sounds good to me. If you have time tomorrow (Friday December 8) I'd be happy to just go through this with you for an hour -- I think it would be useful (at least for me) to talk it through to make sure I understand all of the issues (and hopefully we can agree on the fix in that time as well).

[ptiede]

Let's close this since the PR fixed the major bug and the time issue should be fixed in EHT CASA postproc