astropy / astropy

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

Equivalencies and log-scale units #5339

Open bwinkel opened 8 years ago

bwinkel commented 8 years ago

Equivalencies seem to work only into one direction with log-scale units:

from astropy import units as apu
from math import log10, sqrt

# Note, engineers define the dB(uV/m) to operate on |E|^2 (to make dB-scale equivalent to power)
dB_uV_m = apu.dB(apu.uV ** 2 / apu.m ** 2)

a = 10 * apu.uV / apu.m
b = 20 * dB_uV_m

# Without equivalencies, the following cannot work (of course)
a.to(dB_uV_m)
b.to(apu.uV / apu.m)
# both raise:
# UnitConversionError: 'uV2 / m2' and 'uV / m' (electrical field strength)
# are not convertible

E_field_equivalency = [(
    apu.uV / apu.m,
    dB_uV_m,
    lambda x: 10. * log10(x ** 2),
    lambda x: sqrt(10 ** (x / 10.))
    )]

# Now, the following works:
a.to(dB_uV_m, equivalencies=E_field_equivalency)
# <Decibel 20.0 dB(uV2 / m2)>

# But the other direction is still not ok:
b.to(apu.uV / apu.m, equivalencies=E_field_equivalency)
# UnitConversionError: 'uV2 / m2' and 'uV / m' (electrical field strength)
# are not convertible

# Interestingly, if I use 
apu.add_enabled_equivalencies(E_field_equivalency)
# the Exception that is raised, is a different one:
a.to(dB_uV_m)
# <Decibel 20.0 dB(uV2 / m2)>
b.to(apu.uV / apu.m)
# AttributeError: 'DecibelUnit' object has no attribute '_is_equivalent'
mhvk commented 8 years ago

@bwinkel - thanks for reporting! The underlying trouble is that when you define an equivalency, it should not be directly with the logarithmic unit itself, since that one is "built in" to the logarithmic unit, but rather between the physical units. So, following your script, it works if one does:

from astropy import units as apu
# importing from numpy means you'll be able to handle arrays
from numpy import log10, sqrt
dB_uV_m = apu.dB(apu.uV ** 2 / apu.m ** 2)
a = 10 * apu.uV / apu.m
b = 20 * dB_uV_m
E_field_equivalency = [(
    apu.uV / apu.m,
    apu.uV**2 / apu.m**2,
    lambda x: x ** 2,
    lambda x: sqrt(x)
    )]
b.to(apu.uV / apu.m, equivalencies=E_field_equivalency)
# <Quantity 10.0 uV / m>

I realise the above is somewhat confusing; it is related to the fact that, generally, we cannot "stack" equivalencies -- the only exception is that for logarithmic units we do try to combine the equivalency built-in when you make the unit, and any single other one (#5017). If you can think of a good place to put this in the documentation, do let me know!

mhvk commented 8 years ago

Separately, not being an engineer: I had vaguely read about versions of dB where not a difference of 10 but of 20 dB corresponds to a factor 10, but I couldn't think how to implement this sensibly (and didn't really understand it either...). If this is a very common occurrence, perhaps we should have an additional unit that does this by default (i.e., where you'd have some dB_power(u.V/u.m); maybe there is a standard notation??). Or was your method of defining it as dB(V^2/m^2) a very natural one for you? Even so, it may be possible to have an extra argument that tells the magnitude that the unit is a power, and that thus the default conversion should include the squaring/square-rooting.

bwinkel commented 8 years ago

Thanks a lot for the quick solution, which works nicely. I'm not quite sure, if it could have any strange side effects (the equivalency between E and E^2 may not always be what one naively expects). For my application it will probably be fine, though.

Perhaps, one could indeed add an optional parameter to the "dex" (or dB) unit, which allows additional scaling. After all, the dB is 0.1 * dex, mag is -0.4 * dex, and in my case dB_uV_m (= dB_power) = 0.2 * dex. If I understand the astropy source code correctly, there are quite a few places, which would need to be changed to implement something like the db_power, right?

And to answer the remaining questions: (1) No, the dB(V^2/m^2) feels not natural to me, (2) I've not heard of a version of dB where 20 dB is a factor of ten, but perhaps this was referring to the fact that that you can write E [dB_uV_m] = 10 log_10 E[uV^2 / m^2] or E [dB_uV_m] = 20 log_10 E[uV / m].