Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix display units of moments in velocity #2665

Merged
merged 6 commits into from
Jan 18, 2024

Conversation

rosteen
Copy link
Collaborator

@rosteen rosteen commented Jan 18, 2024

This now converts the spectral axis of the input cube to velocity prior to sending it through the moment calculation rather than attempting to convert the result afterward, thus fixing the display units of the result.

I believe that the reported unrealistic negative values for moment 2 are an artifact of doing continuum subtraction on pixels that are essentially noise around the continuum - the calculation is np.sum(flux * (dispersion - m1) ** order, axis=axis) / m0, so even though the dispersion term is getting squared, with negative flux values (especially far from m1) you can end up with a negative second moment in those cases. Everything looks ok when you do the moment calculation without continuum subtraction. Here's an example of the spectrum at a pixel that results in a negative second moment at that pixel:

Screenshot 2024-01-18 at 10 54 50 AM

@github-actions github-actions bot added cubeviz plugin Label for plugins common to multiple configurations labels Jan 18, 2024
@rosteen rosteen added this to the 3.8.2 milestone Jan 18, 2024
@rosteen rosteen added the bug Something isn't working label Jan 18, 2024
slab = Spectrum1D(slab.flux, slab_sa)

# Finally actually calculate the moment
self.moment = analysis.moment(slab, order=n_moment).T
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I feel like intermediate calculations shouldn't change the actual instance attribute, but it is just a personal preference, so no biggie.

@@ -165,7 +165,7 @@ def test_moment_velocity_calculation(cubeviz_helper, spectrum1d_cube, tmpdir):
label_mouseover = cubeviz_helper.app.session.application._tools['g-coords-info']
label_mouseover._viewer_mouse_event(uncert_viewer, {'event': 'mousemove',
'domain': {'x': 0, 'y': 0}})
assert label_mouseover.as_text() == ("Pixel x=00.0 y=00.0 Value -4.14382e+05 m / s",
assert label_mouseover.as_text() == ("Pixel x=00.0 y=00.0 Value -4.14668e+02 km / s",
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why did the decimal value also change?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think that I'm not terribly surprised that there's a small (<0.1%) difference between my manual conversion of the result vs converting the spectral axis of each spaxel beforehand. Likely due to my manual calculation before using the simpler optical equation, vs now using the full relativistic convention in the conversion through specutils.

assert label_mouseover.as_text() == ("Pixel x=00.0 y=00.0 Value -2.99792e+08 m / s",
"World 13h39m59.9731s +27d00m00.3600s (ICRS)",
"204.9998877673 27.0001000000 (deg)")
'domain': {'x': 1, 'y': 1}})
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why change it to x=1 y=1? Hard to compare the diff to see if results are the same or not before/after patching.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The new moment 2 values are indeed different, that calculation was incorrect before. And I switched to (1,1) because the new value at 0 is 0, which is boring 🙂 .

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That different, huh? I wonder how many scientists trusted the wrong values blindly. 😱

What about the moment 1 differences above?

Copy link

codecov bot commented Jan 18, 2024

Codecov Report

All modified and coverable lines are covered by tests ✅

Comparison is base (eeeeaf8) 91.53% compared to head (2c87000) 91.53%.

Additional details and impacted files
@@           Coverage Diff           @@
##             main    #2665   +/-   ##
=======================================
  Coverage   91.53%   91.53%           
=======================================
  Files         161      161           
  Lines       20107    20103    -4     
=======================================
- Hits        18404    18401    -3     
+ Misses       1703     1702    -1     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

CHANGES.rst Outdated Show resolved Hide resolved
@rosteen rosteen modified the milestones: 3.8.2, 3.9 Jan 18, 2024
Copy link
Contributor

@pllim pllim left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Code looks like code. Test updates are explained satisfactorily.

Do we need Patrick to do science verification and be the second approver?

@rosteen
Copy link
Collaborator Author

rosteen commented Jan 18, 2024

Code looks like code. Test updates are explained satisfactorily.

Do we need Patrick to do science verification and be the second approver?

I would like @PatrickOgle's confirmation that my explanation of the negative moment 2 values holds water, if possible.

@pllim pllim requested a review from PatrickOgle January 18, 2024 17:43
Copy link
Contributor

@PatrickOgle PatrickOgle left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The problems with units and the moment2 calculation are both fixed.

I checked the peak wavelengths of the ‘fastest’ pixels in the map, and they do match the velocities in the moment 0 map, to within ~20 km/s. (To do this, I selected single spaxels with the spaxel tool and used the slice tool to find the peak of line in these spaxel spectra, then computed velocity = c*(lambda_peak - lambda_ref)/lambda.) The ~20 km/s differences are related to the line shape and the accuracy to which I can center the slice marker.

Bottom line--I see no gross problems with the moment calculation for single spaxels. The science is ‘good’ in this case, and shows that stuff is moving faster around a black hole that is more massive than the ground measurements indicated. Makes sense, since JWST has better resolution.

@pllim pllim merged commit 9140858 into spacetelescope:main Jan 18, 2024
15 checks passed
gibsongreen pushed a commit to gibsongreen/jdaviz that referenced this pull request Feb 12, 2024
* Fixing units in result, trying to debug negative second moments

* Remove unneeded import

* Update test

* Add changelog

* Update jdaviz/configs/cubeviz/plugins/moment_maps/moment_maps.py

Co-authored-by: P. L. Lim <[email protected]>

* Fix changelog entry

---------

Co-authored-by: P. L. Lim <[email protected]>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug Something isn't working cubeviz plugin Label for plugins common to multiple configurations
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants