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

flux_conversion and indirect_conversion tests to use flux_conversion_general #3398

Draft
wants to merge 1 commit into
base: main
Choose a base branch
from
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
108 changes: 105 additions & 3 deletions jdaviz/core/tests/test_unit_conversion_utils.py
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
import astropy.units as u
from itertools import combinations
import numpy as np
from numpy.testing import assert_allclose
import pytest

from jdaviz.core.custom_units_and_equivs import PIX2, SPEC_PHOTON_FLUX_DENSITY_UNITS
Expand All @@ -9,6 +9,7 @@
combine_flux_and_angle_units,
flux_conversion_general,
handle_squared_flux_unit_conversions)
from jdaviz.utils import _eqv_pixar_sr


@pytest.mark.parametrize("unit, is_solid_angle", [
Expand Down Expand Up @@ -87,11 +88,112 @@ def test_general_flux_conversion():
equivalencies = all_flux_unit_conversion_equivs(pixar_sr=4, cube_wave=1*u.nm)
for orig, targ, truth in units_and_expected:
converted_value = flux_conversion_general([1], orig, targ, equivalencies)
np.testing.assert_allclose(converted_value[0].value, truth)
assert_allclose(converted_value[0].value, truth)
assert converted_value.unit == targ

# as a bonus, also test the function that converts squared flux units
# (relevant in aperture photometry)
sq = handle_squared_flux_unit_conversions(1, orig**2, targ**2, equivalencies)
np.testing.assert_allclose(sq.value, truth**2, rtol=1e-06)
assert_allclose(sq.value, truth**2, rtol=1e-06)
assert sq.unit == targ ** 2


def test_general_flux_conversion_continuted():

"""
This is effectivley the same test as 'test_general_flux_conversion',
and is included to retain test coverage of a previous flux conversion
function that tested these specific cases, as to not remove any test coverage
since these tests are quite fast.
"""

# the values to be translated in every test case
values = [10, 20, 30]

# unit conversion equivalencies for u.spectral, and flux<>SB for all flux
# conversions supported. Using PIXAR_SR = 0.1 for flux to SB per steradian,
# and a spectral axis for u.spectral equivalency. depending on the conversion,
# all of these equivs might not be required but it is harmless to pass them anyway.
equivs = all_flux_unit_conversion_equivs(pixar_sr=0.1,
cube_wave=[1, 2, 3] * u.um)

test_combos = [(u.Jy / u.sr, u.Jy, [1, 2, 3]),
(u.Jy, u.Jy / u.sr, [100, 200, 300]),
(u.Jy / PIX2, u.Jy, [10, 20, 30]),
(u.Jy, u.Jy / PIX2, [10, 20, 30]),
(u.Jy / u.sr, u.erg / (u.Angstrom * u.s * u.cm**2 * u.sr),
[2.99792458e-12, 1.49896229e-12, 9.99308193e-13]),
(u.Jy / PIX2, u.erg / (u.Angstrom * u.s * u.cm**2 * PIX2),
[2.99792458e-12, 1.49896229e-12, 9.99308193e-13]),
(u.erg / (u.Angstrom * u.s * u.cm**2 * u.sr), u.Jy / u.sr,
[3.33564095e+13, 2.66851276e+14, 9.00623057e+14]),
(u.erg / (u.Angstrom * u.s * u.cm**2 * PIX2), u.Jy / PIX2,
[3.33564095e+13, 2.66851276e+14, 9.00623057e+14]),
(u.erg / (u.Angstrom * u.s * u.cm**2 * u.sr), u.Jy / u.sr,
[3.33564095e+13, 2.66851276e+14, 9.00623057e+14]),
(u.erg / (u.Angstrom * u.s * u.cm**2 * PIX2), u.Jy / PIX2,
[3.33564095e+13, 2.66851276e+14, 9.00623057e+14]),
(u.MJy, u.ph / (u.s * u.cm**2 * u.Hz * u.sr),
[0.00050341, 0.00201365, 0.0045307]),
(u.MJy, u.ph / (u.s * u.cm**2 * u.Hz),
[5.03411657e-05, 2.01364663e-04, 4.53070491e-04]),
(u.ph / (u.s * u.cm**2 * u.Hz * u.sr), u.MJy,
[198644.58571489, 198644.58571489, 198644.58571489]),
(u.ph / (u.s * u.cm**2 * u.Hz * PIX2), u.MJy,
[1986445.85714893, 1986445.85714893, 1986445.85714893]),
(u.Jy / u.sr, u.MJy, [1.e-06, 2.e-06, 3.e-06]),
(u.MJy / PIX2, u.MJy, [10, 20, 30])
]

for orig, targ, expected in test_combos:
converted = flux_conversion_general(values, orig, targ, equivs)
assert_allclose(converted, expected * targ, rtol=1e-05)

# some other misc test cases from previous test
converted = flux_conversion_general(1., u.MJy / u.sr,
u.erg / (u.s * u.cm**2 * u.Hz * u.sr),
equivs)
assert_allclose(converted, 1.e-17 * u.erg / (u.s * u.cm**2 * u.Hz * u.sr))

# another old test case
converted = flux_conversion_general(10., u.MJy / u.sr, u.Jy / u.sr)
assert_allclose(converted, 10000000. * u.Jy / u.sr)

# quantity scale factor
eqv = _eqv_pixar_sr(0.1 * (u.sr / u.pix))
assert_allclose(flux_conversion_general(values, u.Jy / u.sr, u.Jy, eqv),
[1, 2, 3] * u.Jy)
assert_allclose(flux_conversion_general(values, u.Jy, u.Jy / u.sr, eqv),
[100, 200, 300] * u.Jy / u.sr)

# values == 2
assert_allclose(flux_conversion_general([10, 20], u.Jy / u.sr, u.Jy, eqv),
[1, 2] * u.Jy)
assert_allclose(flux_conversion_general([10, 20], u.Jy, u.Jy / u.sr, eqv),
[100, 200] * u.Jy / u.sr)

# array scale factor, same length arrays
res = flux_conversion_general(values, u.Jy / u.sr, u.Jy,
_eqv_pixar_sr([0.1, 0.2, 0.3]))
assert_allclose(res, [1., 4., 9.] * u.Jy)

res = flux_conversion_general(values, u.Jy, u.Jy / u.sr,
_eqv_pixar_sr([0.1, 0.2, 0.3]))
assert_allclose(res, [100, 100, 100] * u.Jy / u.sr)

# array + quantity scale factor, same length arrays
eqv = _eqv_pixar_sr([0.1, 0.2, 0.3] * (u.sr / u.pix))
assert_allclose(flux_conversion_general(values, u.Jy / u.sr, u.Jy, eqv),
[1., 4., 9.] * u.Jy)
assert_allclose(flux_conversion_general(values, u.Jy, u.Jy / u.sr, eqv),
[100, 100, 100] * u.Jy / u.sr)

# values != 2 but _pixel_scale_factor size mismatch
with pytest.raises(ValueError, match="operands could not be broadcast together"):
eqv = _eqv_pixar_sr([0.1, 0.2, 0.3, 0.4])
flux_conversion_general(values, u.Jy / u.sr, u.Jy, eqv)

# # Other kind of flux conversion unrelated to _pixel_scale_factor.
# # The answer was obtained from synphot unit conversion.
targ = [2.99792458e-12, 1.49896229e-12, 9.99308193e-13] * (u.erg / (u.AA * u.cm * u.cm * u.s)) # FLAM # noqa: E501
assert_allclose(flux_conversion_general(values, u.Jy, targ.unit, equivs), targ)
144 changes: 2 additions & 142 deletions jdaviz/tests/test_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,153 +4,13 @@
import photutils
import pytest
from asdf.exceptions import AsdfWarning
from astropy import units as u
from astropy.utils import minversion
from astropy.wcs import FITSFixedWarning
from numpy.testing import assert_allclose
from specutils import Spectrum1D

from jdaviz.core.custom_units_and_equivs import PIX2
from jdaviz.utils import (alpha_index, download_uri_to_path, flux_conversion,
_indirect_conversion, _eqv_pixar_sr)

PHOTUTILS_LT_1_12_1 = not minversion(photutils, "1.12.1.dev")
from jdaviz.utils import (alpha_index, download_uri_to_path)


def test_spec_sb_flux_conversion():
# Actual spectrum content does not matter, just the meta is used here.
spec = Spectrum1D(flux=[1, 1, 1] * u.Jy, spectral_axis=[1, 2, 3] * u.um)

# values != 2
values = [10, 20, 30]

# Float scalar pixel scale factor
spec.meta["_pixel_scale_factor"] = 0.1
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.Jy, spec), [1, 2, 3])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / u.sr, spec), [100, 200, 300])

# conversions with eq pixels
assert_allclose(flux_conversion(values, u.Jy / PIX2, u.Jy, spec), [10, 20, 30])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / PIX2, spec), [10, 20, 30])

# complex translation Jy / sr -> erg / (Angstrom s cm2 sr)
targ = [2.99792458e-12, 1.49896229e-12, 9.99308193e-13] * (u.erg / (u.Angstrom * u.s * u.cm**2 * u.sr)) # noqa: E501
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.erg / (u.Angstrom * u.s * u.cm**2 * u.sr), spec), targ.value) # noqa: E501
# swap sr for pix2 and check conversion
assert_allclose(flux_conversion(values, u.Jy / PIX2, u.erg / (u.Angstrom * u.s * u.cm**2 * PIX2), spec), targ.value) # noqa: E501

# complex translation erg / (Angstrom s cm2 sr) -> Jy / sr
targ = [3.33564095e+13, 2.66851276e+14, 9.00623057e+14] * (u.Jy / u.sr)
assert_allclose(flux_conversion(values, u.erg / (u.Angstrom * u.s * u.cm**2 * u.sr), u.Jy / u.sr, spec), targ.value) # noqa: E501
# swap sr for pix2 and check conversion
assert_allclose(flux_conversion(values, u.erg / (u.Angstrom * u.s * u.cm**2 * PIX2), u.Jy / PIX2, spec), targ.value) # noqa: E501

spectral_values = spec.spectral_axis
spec_unit = u.MJy
eqv = u.spectral_density(spectral_values) + _eqv_pixar_sr(spec.meta["_pixel_scale_factor"])

# test spectrum when target unit in untranslatable unit list
target_values = [5.03411657e-05, 2.01364663e-04, 4.53070491e-04]
expected_units = (u.ph / (u.Hz * u.s * u.cm**2))
for solid_angle in [u.sr, PIX2]:
returned_values, return_units, unit_flag = _indirect_conversion(
values=values, orig_units=(u.MJy),
targ_units=(u.ph / (u.s * u.cm**2 * u.Hz * solid_angle)), # noqa
eqv=eqv,
spec_unit=spec_unit,
indirect_needs_spec_axis=None
)
assert_allclose(returned_values, target_values)
assert (return_units == expected_units)
assert (unit_flag == 'orig')

# test spectrum when original unit in untranslatable unit list
target_values = [1., 2., 3.]
expected_units = (u.ph / (u.Angstrom * u.s * u.cm**2))
returned_values, return_units, unit_flag = _indirect_conversion(
values=values,
orig_units=(u.ph / (u.Angstrom * u.s * u.cm**2 * u.sr)), # noqa
targ_units=(u.MJy), eqv=eqv,
spec_unit=spec_unit, indirect_needs_spec_axis=None # noqa
)
assert_allclose(returned_values, target_values)
assert (return_units == expected_units)
assert (unit_flag == 'targ')

# test the default case where units are translatable
target_values = [10, 20, 30]
expected_units = (u.MJy)
returned_values, return_units, unit_flag = _indirect_conversion(
values=values, orig_units=(u.Jy/u.sr),
targ_units=(u.MJy), eqv=eqv,
spec_unit=spec_unit, indirect_needs_spec_axis=None # noqa
)
assert_allclose(returned_values, target_values)
assert (return_units == expected_units)
assert (unit_flag == 'targ')

# test image viewer data units are untranslatable
target_value = 1.e-18
expected_units = (u.erg / (u.s * u.cm**2 * u.Hz))
returned_values, return_units = _indirect_conversion(
values=1, orig_units=(u.MJy/u.sr),
targ_units=(u.erg / (u.s * u.cm**2 * u.Hz * u.sr)),
eqv=eqv, spec_unit=None, indirect_needs_spec_axis=True
)
assert_allclose(returned_values, target_value)
assert return_units == expected_units

# test image viewer data units are translatable
target_value = 10
expected_units = (u.MJy / u.sr)
returned_values, return_units = _indirect_conversion(
values=10, orig_units=(u.MJy/u.sr), targ_units=(u.Jy/u.sr),
eqv=eqv, spec_unit=None, indirect_needs_spec_axis=True
)
assert_allclose(returned_values, target_value)
assert return_units == expected_units

# Quantity scalar pixel scale factor
spec.meta["_pixel_scale_factor"] = 0.1 * (u.sr / u.pix)
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.Jy, spec), [1, 2, 3])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / u.sr, spec), [100, 200, 300])

# values == 2
values = [10, 20]
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.Jy, spec), [1, 2])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / u.sr, spec), [100, 200])

# float array pixel scale factor
spec.meta["_pixel_scale_factor"] = [0.1, 0.2, 0.3] # min_max = [0.1, 0.3]
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.Jy, spec), [1, 6])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / u.sr, spec), [100, 66.66666666666667])

# Quantity array pixel scale factor
spec.meta["_pixel_scale_factor"] = [0.1, 0.2, 0.3] * (u.sr / u.pix) # min_max = [0.1, 0.3]
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.Jy, spec), [1, 6])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / u.sr, spec), [100, 66.66666666666667])

# values != 2
values = [10, 20, 30]
spec.meta["_pixel_scale_factor"] = [0.1, 0.2, 0.3]
assert_allclose(flux_conversion(values, u.Jy / u.sr, u.Jy, spec=spec), [1, 4, 9])
assert_allclose(flux_conversion(values, u.Jy, u.Jy / u.sr, spec=spec), 100)

# values != 2 but _pixel_scale_factor size mismatch
with pytest.raises(ValueError, match="operands could not be broadcast together"):
spec.meta["_pixel_scale_factor"] = [0.1, 0.2, 0.3, 0.4]
flux_conversion(values, u.Jy / u.sr, u.Jy, spec=spec)

# Other kind of flux conversion unrelated to _pixel_scale_factor.
# The answer was obtained from synphot unit conversion.
spec.meta["_pixel_scale_factor"] = 0.1
targ = [2.99792458e-12, 1.49896229e-12, 9.99308193e-13] * (u.erg / (u.AA * u.cm * u.cm * u.s)) # FLAM # noqa: E501
assert_allclose(flux_conversion(values, u.Jy, targ.unit, spec=spec), targ.value)

# values == 2 (only used spec.spectral_axis[0] for some reason)
values = [10, 20]
targ = [2.99792458e-12, 5.99584916e-12] * (u.erg / (u.AA * u.cm * u.cm * u.s)) # FLAM
assert_allclose(flux_conversion(values, u.Jy, targ.unit, spec=spec), targ.value)
PHOTUTILS_LT_1_12_1 = not minversion(photutils, "1.12.1.dev")


@pytest.mark.parametrize("test_input,expected", [(0, 'a'), (1, 'b'), (25, 'z'), (26, 'aa'),
Expand Down
Loading