astropy / astropy

Astronomy and astrophysics core library
https://www.astropy.org
BSD 3-Clause "New" or "Revised" License
4.37k stars 1.74k forks source link

Arithmetic operation on dimensionless mag units #12472

Open olebole opened 2 years ago

olebole commented 2 years ago

Description

Not sure whether this is a bug or just a misunderstanding: I am wondering how dimensionless units work; especially the conversion to linear.

Expected behavior

I would expect that u.mag works equivalently to f.e. u.ABmag when it comes to conversion:

>>> import astropy.units as u
>>> m = 12 * u.ABmag
>>> print(m)
12.0 mag(AB)
>>> print(m.to(u.erg/u.s/u.m**2/u.Hz))
5.754399373371543e-21 erg / (Hz m2 s)
>>> flux = 1e-20 * u.erg/u.s/u.m**2/u.Hz
>>> print(flux.to(u.ABmag))
11.399999999999997 mag(AB)

i.e. u.mag should be convertible to and from a dimensionless unit.

Actual behavior

>>> import astropy.units as u
>>> m = 1 * u.mag
>>> print(m)
1.0 mag
>>> print(m.to(''))
Traceback (most recent call last):
[…]
  File "/usr/lib/python3/dist-packages/astropy/units/core.py", line 1027, in _apply_equivalencies
    raise UnitConversionError(
astropy.units.core.UnitConversionError: 'mag' and '' (dimensionless) are not convertible
>>> f = u.Quantity(10.0)
>>> print(f.to(u.mag))
Traceback (most recent call last):
[…]
  File "/usr/lib/python3/dist-packages/astropy/units/core.py", line 1027, in _apply_equivalencies
    raise UnitConversionError(
astropy.units.core.UnitConversionError: '' (dimensionless) and 'mag' are not convertible

Specific question: what is the "preferred" way to convert to copy to/from dimensionless magnitudes (except explicitly using -2.5*log10(f))?

Steps to Reproduce

See above

System Details

Astropy 5.0rc1

pllim commented 2 years ago

I don't understand your examples. Currently you also cannot do (12 * u.ABmag).to(''), so what is so inconsistent between u.mag and u.ABmag? Also, it can be misleading if they can both blindly convert to unitless, because ABmag has an implied zeropoint and who knows what assumption went into mag because that info is already lost.

What is the use case that would justify the usage of .to('')?

pllim commented 2 years ago

"preferred" way

You can try using https://synphot.readthedocs.io/en/latest/api/synphot.units.convert_flux.html if you want.

olebole commented 2 years ago

"mag" may be the difference between two magnitudes, which would be equivalent to the ratio of fluxes:

>>> m1 = 5*u.ABmag
>>> m2 = 3*u.ABmag
>>> m1-m2
<Magnitude 2. mag>
>>> (m1 - m2).physical  # This works
<Quantity 0.15848932>
>>> (m1 - m2).to('') # This doesn't

Use case for the reverse case: I have a star with 15mag(AB). Another one is twice as bright. How do I calculate its magnitude? I would do it this way:

m1 = 15*u.ABmag
fac = 2.0
m2 = m1 + u.Quantity(fac).to(u.mag)

This does not work, and I didn't find a workaround (except using the plain formula, which contradicts a bit the idea of units).

pllim commented 2 years ago

This seems to work?

>>> fac = 2 * u.dex
>>> m1 + fac
<Magnitude 10. mag(AB)>
olebole commented 2 years ago

Ah - I never heared of dex. What is the difference to magnitudes?

pllim commented 2 years ago

Maybe @mhvk can explain, but when I define fac = 2 * u.ABmag, the unit of m2 is wrong.

pllim commented 2 years ago

Or perhaps https://docs.astropy.org/en/latest/units/logarithmic_units.html can be improved?

olebole commented 2 years ago

I think also with dex the calc is wrong. Normally, it is mag = -2.5*log10(f2/f1), which would for f2/f1=2 give -0.7526. Or, more explicit:

>>> m1 = 15*u.ABmag
>>> fac = 2.0
>>> m2 = (fac * m1.to('erg / (Hz m2 s)')).to(u.ABmag)
>>> print(m2)
14.247425010840047 mag(AB)
pllim commented 2 years ago

Maybe we should just stop using the magnitude system...

mhvk commented 2 years ago

@olebole - to your initial issue, this is a bit of a leftover from the fact that mag existed as a unit before the MagUnit class. Hence, one can have a Quantity with units of mag but also a Magnitude and those have different properties. It is particularly difficult because mag can also be a magnitude difference, which unlike a normal magnitude, can have more operations (e.g., mag/day makes sense for a magnitude difference, not for, say, STmag).

mhvk commented 2 years ago

@olebole - I don't understand why you think the calculation in https://github.com/astropy/astropy/issues/12472#issuecomment-967166632 is wrong: you give a two times higher flux and get a magnitude which is 0.75 brighter.

If you want to calculate a star that is two times brighter, you'd do

m1 + u.Magnitude(2.*u.one)
mhvk commented 2 years ago

I'm confused about one example; I get:

In [10]: m1 = 5*u.ABmag

In [11]: m2 = 3*u.ABmag

In [12]: (m1 - m2).to('')
Out[12]: <Quantity 0.15848932>
mhvk commented 2 years ago

Sorry for the chain of comments, but the following does work:

In [13]: (2*u.one).to(u.mag())
Out[13]: <Magnitude -0.75257499 mag>

It is indeed annoying that u.mag and u.mag()` are different

olebole commented 2 years ago

I see there is a difference between u.mag and u.mag():

>>> 1*u.mag()
<Magnitude 1. mag>
>>> (1*u.mag()).to('')
<Quantity 0.39810717>
>>> 1*u.mag
<Quantity 1. mag>
>>> (1*u.mag).to('')
Traceback (most recent call last):
[…]

This is a bit confusing... Can this be harmonized?