SWAPI: add L3A alpha and proton velocity vectors in GSE, GSM, and HAE-J2000 coordinates - #165
SWAPI: add L3A alpha and proton velocity vectors in GSE, GSM, and HAE-J2000 coordinates#165leowerneck wants to merge 8 commits into
Conversation
|
@jtniehof @hafarooki — Could you please review if this is the intended approach? I followed @hafarooki's suggestion and rotated the alpha and proton velocity vectors from the RTN frame. For now, I focused on adding these only to SWAPI. |
hafarooki
left a comment
There was a problem hiding this comment.
Looks good. Some minor things: a comment would help in the code for the einsum. And the test cases should perhaps flip around signs and axes instead of added 300 to all components…
|
@jtniehof requested for review of metadata aspects |
@hafarooki That's a good point. I updated the tests to add rotations. |
There was a problem hiding this comment.
Pull request overview
Adds SWAPI L3A proton and alpha velocity products in GSE, GSM, and HAE-J2000 coordinates.
Changes:
- Adds SPICE-based vector and covariance transformations.
- Extends CDF models, metadata, labels, and documentation.
- Updates processor and model tests for new products.
Reviewed changes
Copilot reviewed 11 out of 11 changed files in this pull request and generated 3 comments.
Show a summary per file
| File | Description |
|---|---|
imap_l3_processing/swapi/swapi_processor.py |
Generates target-frame products. |
imap_l3_processing/swapi/l3a/utils.py |
Implements coordinate transformations. |
imap_l3_processing/swapi/l3a/models.py |
Exposes new CDF variables and labels. |
imap_l3_processing/cdf/config/imap_swapi_l3a_variable_attrs.yaml |
Defines variable metadata. |
docs/swapi/proton-sw.md |
Documents proton outputs. |
docs/swapi/alpha-sw.md |
Documents alpha outputs. |
tests/test_utils.py |
Updates model fixtures. |
tests/swapi/test_swapi_processor.py |
Tests processor integration. |
tests/swapi/l3a/test_utils.py |
Tests transformation mathematics. |
tests/swapi/l3a/test_models.py |
Tests new model variables. |
tests/swapi/l3a/test_chunk_fits.py |
Isolates target-frame processing in existing tests. |
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
There was a problem hiding this comment.
Pull request overview
Copilot reviewed 11 out of 11 changed files in this pull request and generated 1 comment.
Suppressed comments (1)
imap_l3_processing/swapi/swapi_processor.py:68
- Each call computes TT-to-ET and the same RTN-to-target rotation matrix, and this function repeats that work for the Sun-relative vector and covariance at every epoch and frame. That causes three SPICE frame lookups per output frame per record. Compute the rotation once per epoch/frame (or batch a frame and fall back per epoch only on a SPICE coverage error), then apply it to both vectors and the covariance.
velocity = convert_velocity_rtn_to_frame(
epoch_slice,
result[f"{prefix}_rtn"][index : index + 1],
target_frame,
)[0]
There was a problem hiding this comment.
Pull request overview
Copilot reviewed 11 out of 11 changed files in this pull request and generated no new comments.
Suppressed comments (1)
imap_l3_processing/swapi/swapi_processor.py:65
- Each epoch/frame computes the same timestamp conversion and SPICE rotation three times: once for each velocity and once for the covariance. At the one-minute product cadence, this adds 9 SPICE rotations per record (12,960 per day/product) although only 3 unique rotations are needed. Compute one rotation matrix per epoch/frame and apply it to both vectors and the covariance, or pass that matrix into the conversion helpers.
velocity = convert_velocity_rtn_to_frame(
@leowerneck I'm not sure Copilot is right about this. Isn't that the behavior we had before (fill values with FIT_ERROR flag if any error occurs?) |
|
Something this PR could benefit from is integration test updates ensuring that the rotated vectors all make sense given our apriori understanding of these coordinate systems. |
Does fixing #171 fix this? |
Co-authored-by: Copilot Autofix powered by AI <175728472+Copilot@users.noreply.github.com>
…ame conversion failures.
|
@jtniehof this PR used the rotation matrix approach instead of the 6D transformation. Do you think thats good? |
Yes! Running: uv run --with pytest,pytest-xdist pytest -n auto tests/integration/test_swapi_processor_integration.pyresults in 4 failed tests. But if I rebase first: git rebase 171-172-174-fix-swapi-integration-testing
uv run --with pytest,pytest-xdist pytest -n auto tests/integration/test_swapi_processor_integration.pythen all 4 tests pass. |
jtniehof
left a comment
There was a problem hiding this comment.
Structurally this is pretty sound. I think we do need the full 6D transform, and to get GSM/GSE to be Earth-centered.
Given that this is adding so many more variables, I also think this needs a major version bump to 2, and a MODS global attribute describing the changes (once the set of changes is fully settled :) ). See #173 for more description.
| DataProductVariable(VELOCITY_GSE_LABEL_CDF_VAR_NAME, value=["X", "Y", "Z"]), | ||
| DataProductVariable(VELOCITY_GSM_LABEL_CDF_VAR_NAME, value=["X", "Y", "Z"]), | ||
| DataProductVariable(VELOCITY_HAE_LABEL_CDF_VAR_NAME, value=["X", "Y", "Z"]), |
There was a problem hiding this comment.
SPDF is concerned about having LABL_PTR_1 that was incomplete. They want the contents of LABL_PTR_1 to be understandable without having to look at LABLAXIS. See #140. It looks like this is the case with the RTN labels...I think we had left that one and just weren't using it.
But at any rate these (referring to following lines as well) should be something like Vp HAE X (quantity, coordinate, and dimension). Sun vs. S/C I'll get to in another comment :)
| def convert_velocity_rtn_to_frame( | ||
| epoch_tt2000_ns: ndarray, | ||
| velocity_rtn: ndarray, | ||
| target_frame: SpiceFrame, | ||
| ) -> ndarray: | ||
| """Express RTN velocity vectors in another coordinate frame.""" | ||
| ephemeris_times = ttj2000ns_to_et(epoch_tt2000_ns) | ||
| rotation_matrices = get_rotation_matrix( | ||
| ephemeris_times, SpiceFrame.IMAP_RTN, target_frame | ||
| ) | ||
| return np.einsum("nij,nj->ni", rotation_matrices, velocity_rtn) |
There was a problem hiding this comment.
There's two issues with this approach:
- GSM, GSE, and RTN are all rotating frames relative to the inertial frame; HAE is not. GSM and GSE are rotating at the same rate and RTN is close (being based at L1) but not exactly. So I do think you need
ckfxfm_cfor this and the full 6D transform. Fortunately the position in RTN is[0, 0, 0]. - GSM and GSE as full coordinate systems use Earth as the origin / standard of rest. I don't think there's any point to expressing spacecraft-relative or Sun-relative velocities in GSM or GSE, so we really only need the Earth-relative.
ckfxfm_cwill not correct for this (only for the coordinate rotation). So wherechunk_fitscallsget_spacecraft_velocity_rtn, there needs to be another call toimap_statethat gets the state relative to Earth, instead of relative to the Sun. I think everything after that should fall through. This means the variables do not need to have in the name / description that they're relative to Earth, there's only the one, although mentioning it in the VAR_NOTES might make sense. - HAE doesn't have additional velocity transforms to do; technically with HAE it should be Sun-referred and we shouldn't have a S/C version. But we do use S/C-referenced HAE with the mappers so I'm not opposed. I'd lean towards having the Sun-referenced one just be HAE, and the S/C-referenced one explicit, again using the VAR_NOTES to be more detailed on all of them, but I yield to @hafarooki 's opinion on that front.
As a minor quibble, I tend to find matmul easier to comprehend than einsum, but not enough to feel comfortable requesting a change unless somebody else says "oh yeah!"
There was a problem hiding this comment.
Actually I'm no longer convinced the appropriate RTN position is [0, 0, 0] because I keep forgetting that SPICE frames don't have an origin :) Probably best to use the full state vector out of imap_state relative to the Sun or Earth, as appropriate.
| _FRAME_ROTATIONS = { | ||
| "gse": np.array([[0, -1, 0], [1, 0, 0], [0, 0, 1]]), | ||
| "gsm": np.array([[0, 0, 1], [-1, 0, 0], [0, -1, 0]]), | ||
| "hae": np.array([[-1, 0, 0], [0, 0, 1], [0, 1, 0]]), | ||
| } |
There was a problem hiding this comment.
I take this to mean the unit tests are really faking the frames? In that case, I'd agree with an integration test for the real SPICE kernels, which I think is coming into the realm of possibility with #175
| proton_sw_velocity_gse_sun: | ||
| NAME: proton_sw_velocity_gse_sun | ||
| DATA_TYPE: CDF_REAL4 | ||
| FIELDNAM: Proton Solar Wind Velocity (GSE, Sun Frame) | ||
| CATDESC: Solar wind proton bulk velocity vector in GSE coordinates (Sun frame) | ||
| DEPEND_0: epoch | ||
| LABL_PTR_1: proton_sw_velocity_gse_sun_label | ||
| VAR_TYPE: data | ||
| RECORD_VARYING: RV | ||
| DISPLAY_TYPE: time_series | ||
| VARIABLE_PURPOSE: PRIMARY_VAR,SUMMARY | ||
| FORMAT: F19.3 | ||
| LABLAXIS: Vp GSE Sun | ||
| UNITS: km/s | ||
| VALIDMIN: -10000 | ||
| VALIDMAX: 10000 | ||
| FILLVAL: -1.000000e+31 | ||
| SCALETYP: linear | ||
| DELTA_PLUS_VAR: proton_sw_velocity_gse_uncert | ||
| DELTA_MINUS_VAR: proton_sw_velocity_gse_uncert |
There was a problem hiding this comment.
Only looked at this one in detail for metadata check :) -- actually, given that, it's maybe worth using the YAML "inheritance" for all these velocity variables. Units, format, min/max, scaletyp..
- Some of my notes regarding which standards of rest to use obviously have implications here
- I'm not sure we want PRIMARY_VAR,SUMMARY on all of these. @hafarooki , what's your feeling on the most significant 4 or 5 velocity variables? RTN s/c centered, GSE, HAE Sun centered? Maybe RTN Sun-centered?
- This applies to the existing ones as well, but I might want to have SCALEMIN/SCALEMAX of +/-800 (or maybe 1000). Feel free to bump that as out of scope for this PR
- Similarly I'd like F7.0 for the FORMAT on all of these. 19.3 is nine leading blanks that will never be used, and I don't think we need m/s precision if just dumping this out in a text table...again, applies to all, feel free to say out of scope
There was a problem hiding this comment.
@jtniehof I'm not quite sure. I think most people will use RTN and/or GSE. Both of them would be useful in Sun frame. For S/C frame, I'm not sure which one end users would find more helpful. I'm not quite sure what end users would do with it exactly... it would probably vary a lot case by case.
There was a problem hiding this comment.
Let's do both RTN, HAE Sun frame, and GSE (which should be in Earth frame). The rest should have no VARIABLE_PURPOSE then.
There was a problem hiding this comment.
I think we would want GSE in SC and sun frame too, perhaps in addition to Earth frame. WIND reports velocity in GSE but in the Sun frame for SWE, but in SC frame for 3DP. I’m not sure they do Earth frame at all
There was a problem hiding this comment.
Wind doing it is a good argument for what people might expect. Although I don't like it :)
The IMAP coordinate frames document defines GSE and GSM explicitly as geocentric. I think I'd better check in with Matina and Drew on this...we want to try and be consistent across the instruments as we add more, and I really don't want to be doing ten or so variables for every velocity. It's work and it's stuff for people to sort through when using the files.
Change Summary
Closes #21
Overview
Add SWAPI L3A proton and alpha velocity vectors, covariance matrices, and uncertainties in GSE, GSM, and HAE-J2000 coordinates.
File changes
Testing