spedas / pyspedas

Python-based Space Physics Environment Data Analysis Software
https://pyspedas.readthedocs.io/
MIT License
143 stars 58 forks source link

ERG cotrans routine can introduce NaNs #920

Open jameswilburlewis opened 1 week ago

jameswilburlewis commented 1 week ago

It looks like something can go wrong converting from ERG dsi coordinates to sga coordinates (dsi->sgi is OK, but there seems to be an issue with sgi->sga). The input data has been cleaned of NaNs, but the output has NaNs, and values are mismatched if you try to do a round-trip comparison. See this excerpt from ERG cotrans_test.py:

import numpy as np
import pyspedas
from pyspedas.erg import mgf
from pyspedas.erg.satellite.erg.common.cotrans.erg_cotrans import erg_cotrans
from pytplot import store_data, get_data
[...]
        mgf_vars=pyspedas.erg.mgf()
        # Clean NaNs from input data
        orig_time, orig_data = get_data('erg_mgf_l2_mag_8sec_dsi')
        orig_meta = get_data('erg_mgf_l2_mag_8sec_dsi',metadata=True)
        good_locations = np.where(np.isfinite(orig_data))
        good_times_idx = np.unique(good_locations[0])
        good_times = orig_time[good_times_idx]
        good_vals = orig_data[good_times_idx,:]
        store_data('erg_mgf_l2_mag_8sec_clean_dsi',{'x':good_times,'y':good_vals}, attr_dict=orig_meta)

        valid_coords=['dsi','sgi','sga','j2000']
        erg_cotrans('erg_mgf_l2_mag_8sec_clean_dsi', out_coord='sgi')
        dsi_sgi = get_data('erg_mgf_l2_mag_8sec_clean_sgi')
        sgi_nan = np.where(np.isnan(dsi_sgi[1]))
        self.assertTrue(len(sgi_nan[0]) == 0)
        erg_cotrans('erg_mgf_l2_mag_8sec_clean_dsi', out_coord='sga')
        dsi_sga = get_data('erg_mgf_l2_mag_8sec_clean_sga')
        sga_nan = np.where(np.isnan(dsi_sga[1]))
        # This turns out to have NaNs in it
        # self.assertTrue(len(sga_nan[0]) == 0)

Maybe there are NaNs in the attitude data used for the coordinate transform?

jameswilburlewis commented 1 week ago

The problem seems to be in the attitude files themselves: erg_att_interpolate() gets RA and Dec variables from the attitude CDF for the axes SGI-Z, SGA-X, and SGA-Z. The differences between the SGI and SGA axes should be very small. SGA-Y is obtained by SGA-Z cross SGA-X. SGI-X comes from SGI-Z cross SGA-X. SGA-Y comes from SGA-Z cross SCA-X.

There is a period of several minutes on the test date, 2017-03-27, where SGI-Z and SGA-X are identical or nearly so, where for the rest of the day they are nearly orthogonal, as expected. The cross product of the parallel vectors gives 0-vectors as the SGI-X output, then when they're normalized to unit vectors, NaNs are generated and propagate through the cotrans calculation.

This doesn't seem recoverable without fixing the attitude files. Let's see what the ERG team says about it.