@pllim - thanks! Should this go in a astropy.photometry package (into which we will later add aperture and PSF photometry)? If we add in things like filters etc it seems all that should belong together?

Very exciting! One general comment: my guess would be that soon enough we will want to be able to convert between different magnitude systems, i.e., something like

m = Magnitude(xxx, system='ST')
m.AB
m.SDSS(transformation='fuku,...')

This will likely require a magnitude superclass, which contains a magnitude in a given system (just like Time contains a time in a given scale, and SkyCoordinates will soon contain an object holding the coordinates in a given reference frame). As far as I can see from a quick glance, nothing in what you've done would stop this, except that perhaps the Magnitude name itself should be reserved for the superclass rather than the base class.

p.s. The units.py file has on l29:

# TODO: Support functional units, e.g. log(x), ln(x)

Maybe @mdboom can comment whether he had specific plans, as obviously that would be useful here.

I cannot find astropy.photometry in the current master. Is it another active PR?

@pllim - I meant I wonder whether this PR should create astropy.photometry. Maybe @eteq has some ideas?

Or maybe astropy.units?

I'd wondered about units, as in practice people will mostly use quantities (which is strangely hidden in the documentation). I'd certainly only not move magnitudes to that group unless we rename it. But perhaps think first about the content rather than the location!

@astrofrog - Will naming it astropy.photometry confuse users with the affiliated package photutils? If not, I am okay with renaming this package if no one else speaks against it. I would rather not put this into astropy.units, which I think should only contain building blocks for more complex unit/quantity. Furthermore, units is deemed stable but not this package.

@mhvk - I could rename Magnitude to GenericMagnitude or BaseMagnitude, whichever is preferred. Thanks for the comment on ABMAG, I'll update the doc in the next round.

Am I right that this only applies to spectral flux values? Do you have a plan for magnitudes in passbands?

@kbarbary - Not in this PR. But the class is written to be expandable to include other magnitude systems.

@mhvk I like your idea of a system argument, but I think it should be tied to a registry of MagSystem classes underneath. MagSystem and derived classes would be responsible for converting between flux (both spectral and integrated) and magnitude. Sample inheritance diagram:

MagSystem --> ABMagSystem
          --> STMagSystem
          --> SpectralMagSystem (e.g. Vega, BD+17)

The user would see the same thing:

m = Magnitude(xxx, system='ab')
m = Magnitude(xxx, system='st')
m = Magnitude(xxx, system='vega')

But underneath there would be a lookup registry of systems, like this:

_magsystems = {'ab': ABMagSystem(),
               'st': STMagSystem(),
               'vega': SpectralMagSystem(<spectrum of vega>),
               'bd17': SpectralMagSystem(<spectrum of bd+17>)}

Here is a sample API for MagSystem:

>>> ab = ABMagSystem()

# Spectral flux of a magnitude 0 object at specific wavelength
>>> ab.zpflux(4000.*u.AA)
<Quantity xxx erg / (cm2 s A)>

# Integrated flux of a magnitude 0 object in some_bandpass:
# Once calculated for a given bandpass, this value is cached within the MagSystem instance, so it is very fast.
>>> ab.zpbandflux(<some_bandpass>)
<Quantity xxx ph / (cm2 s)>

The zpflux and zpbandflux methods above are the important ones. The rest is just log and division stuff that will be identical for every magnitude system. You could also have convenience methods that do the math for you, or this functionality could be in the Magnitude class. If it were here, it might look like:

# Convert flux to magnitude. These would return either a Magnitude object or a float.
>>> ab.flux_to_mag(1. * (u.erg / u.s / u.cm**2 / u.AA), 4000.*u.AA)
>>> ab.band_flux_to_mag(1. * (u.ph * u.s / u.cm**2), <some_bandpass>)

# Convert magnitude to flux. These would return Quantities.
>>> ab.mag_to_flux(..., 4000.*u.AA)
>>> ab.band_mag_to_flux(..., <some_bandpass>)

It is very important that this "magnitude system" functionality be available outside of the Magnitude class, for cases where one is measuring a negative flux. Very often, I measure a negative flux on a CCD that is given a physical meaning by (1) a zeropoint (float) and (2) associated magnitude system for that zeropoint. Converting the negative flux value to a Magnitude is either not possible or is conceptually undesirable. For this purpose, I've been using similar MagSystem classes and a registry of built-in magnitude systems in sncosmo for some time now and it has worked well for me.

@kbarbary - When spectrum or passband is involved, you need classes for them, which is out of the scope of this sub-package but will hopefully be handled by http://github.com/spacetelescope/pysynphot (in active development and not ready for review yet).

@pllim Absolutely this depends on having Spectrum and Bandpass classes, which I sort of left as "blanks" in my examples. My intent was to show how they would be integrated in the future, not necessarily to suggest that they be part of this PR. For the current PR, one could simply define ABMagSystem and STMagSystem and only define the method zpflux. This could be done without reference to any spectrum or bandpass.

I think it is important to think now about how the design would work with both spectra and bandpasses.

@pllim BTW, what is the relation between pysynphot development and specutils? They seem to have overlapping planned functionality.

@kbarbary - As regard to pysynphot development, my priority is to have it working for STScI first. Then I'll coordinate with specutils developers. Let's take this discussion off this PR though -- feel free to contact me privately.

We are looking at that now. I think we will try to have it use the model framework much more directly in the astropy version we are developing, so there ought to be a very strong overlap of spectral models with what is in pysynphot. One could almost think of it now as extending models to handle spectral issues.

I'm also going to try to get some discussion going on what spectral analysis priorities there are to help guide our other work in this area.

@pllim BTW, what is the relation between pysynphot development and specutils? They seem to have overlapping planned functionality.

Real annoying point here: in other packages, we opt for verbosity. Can you change Mag to Magnitude wherever relevant? Also why all caps? ABMAG -> ABMagnitude

@kbarbary Also also, I'm against string systems -- why not pass in the class itself?

import astropy.whatever as yy
m = yy.Magnitude(xxx, system=yy.ABMagnitude)

@adrn - don't agree with the argument against string systems -- I like how those can match attributes, as in
ABmag(...).stmag

@adrn I'm imagining that the system keyword could accept either a MagSystem instance or a string corresponding to a specific MagSystem instance in the registry:

m = Magnitude(xxx, system='vega')
m = Magnitude(xxx, system=SpectralMagSystem(<spectrum of Vega>))
m = Magnitude(xxx, system='ab')
m = Magnitude(xxx, system=ABMagSystem())

I guess I'm wondering why you need both? Shouldn't there be "one -- and preferably one -- way to do it"? It's a difference of one line:

vega = SpectralMagnitudeSystem(<spectrum of Vega>)
m = Magnitude(xxx, system=vega)

@adrn In your example, you'd have to load the spectrum of vega first, so it is actually more than one line. Using the string argument/registry, this would be taken care of behind the scenes. I admit that it is simpler for the AB system, where it would just be

ab = ABMagSystem()
m = Magnitude(xxx, system=ab)

The other thing is that there need only be a single instance of any given magnitude system. (Well there could be multiple instances, but it will only lead to inefficient code (they wouldn't share cached zeropoint vales for example). The registry would take care of this, so that multiple calls of Magnitude(xxx, system='ab') use the same instance of ABMagSystem. If you didn't use strings, you'd have to keep track of your one ab = ABMagSystem() instance throughout your code. Not hard, but an added headache.

In short, the preferable way to do it would be using strings, but for flexibility it would also accept MagSystem instances, for custom magnitude systems.

A possible problem with Magnitudes as we've been writing them: in order to convert between arbitrary magnitude systems, you'd need to know the wavelength at which the spectral flux density is observed. E.g.,

m = Magnitude(1.*u.erg / u.s/u.cm**2/u.AA, system='ab')
m.vega # ??? can't convert from AB mags to Vega mags without knowing wavelength

@kbarbary - I think your one-but-last point about system (as well as my earlier one) brings up a more general worry about having different magnitude subclasses. Indeed, this is becoming eerily like Coordinates, where the design started with different FK4, ICRS, etc, subclasses, but it has now been realised that this doesn't work, as, e.g., ICRS needs to carry information about equinox and epoch for it to be able to convert to FK4.

I think that like for coordinates, it will be good to keep two worries separate: how simply to transform between flux and magnitude (this PR, as @pllim stressed), and how to deal with magnitude systems (a future PR; @pllim probably has ideas given her work on pysynphot). I guess in other words, there should be a basic magnitude which allows

mag(flambda)
mag(fnu)
mag.flambda
mag.fnu

(and photon rates as well...)

You are right, though, that this means the magnitude does need information about the wavelength, so it can convert one type of flux to the other (using u.spectral_density(wavelength)). The problem, of course, is that which wavelength one uses depends on the actual spectral shape, detailed band-pass, etc. So perhaps, at this point, this is the best we can do...

@pllim: it would be very helpful to have an idea about what "system" magnitude API you are envisioning for synphot; do you have a link? Also, separately, to what extent is the current PR useful for pysynphot? Is it not too simplistic?

@mhvk

Currently in synphot (https://github.com/spacetelescope/pysynphot/), there is no distinct magnitude system API. They are defined as different variants of u.mag in synphot/units.py and the actual conversion takes place in synphot/spectrum.py under convert_fluxes() function, which actually does a two-pass conversion that centers around PHOTLAM (otherwise I would have to define a plethora of equivalencies for all possible combinations).

As for my vision, I am hoping to move the simpler flux conversions to spectral_density() equivalencies (#1560). And have a simple mag-flux conversion via this PR. I purposefully excluded OBMAG and VEGAMAG from this PR for now (although they are supported by synphot) because they are complicated; one involves primary mirror area and pixel widths, the other depends on which Vega spectrum one prefers. And then I am hoping to use these pieces in synphot to handle the actual conversion involving wavelengths, passband, spectrum, etc.

I am certainly open to suggestions. There seem to be a lot of good ideas so far in the comments already.

I'll have some more specific comments here later, but just to address the location:

I agree with @astrofrog that astropy.magnitudes is not good, but I pause a bit about astropy.photometry, because that makes me think it actually does photometry. But in my comments (that I will make later), I'm going to say this should be combined with a Flux class to make it useful. And "photometry" is technically the measurement of fluxes... So is there some synonym to "photometry" that won't make people think it actually does photometry?

@aconley - my definitive source on everything tells me irradiance equals radiative flux [1]. So, still W/m2?
[1] http://en.wikipedia.org/wiki/Irradiance

@mhvk - Yes, yet another synonym, again in W m^-2. There's a welter of different systems here, and astronomers and radiometers often have different names for the same thing. That's why saying 'flux' is so confusing -- it is used to mean different things by different people in different communities. To a radiometer it should mean 'radiant flux', which is Luminosity to an astronomer. To a photometer, it's luminous flux (in lumens). To an astronomer, flux seems to mean 'whatever I'm talking about right now that has anything whatsoever to do with radiometry or photometry.'

Flux density, as far as I know, is unique (although it has a synonym, spectral irradiance, for radiometers, and
sometimes people put a spectral in front (spectral flux density)).

But, anyways, in practice, telescope pipelines return either 1) magnitudes, 2) counts, 3) photons, or 4) (spectral) flux densities. So those are the cases that are important to handle. 2 and 3 are special because they aren't universal -- you need more information about your telescope/detector/atmosphere, encapsulated as a zeropoint, to go anywhere with them. But you always know exactly how to go between, say, AB mags and a flux density.

Irradiance/radiative flux is proportional to the thing you need for a broadband magnitude, but magnitude systems are cleverly arranged so that you don't need to know the constant of proportionality -- meaning you don't need to know the true system throughput, just it's wavelength dependence. The conventions around broadband flux density, as returned by telescope pipelines, are also set up so that you don't need to know the true throughput -- which is good, because you never know it!

@mhvk - I agree with @mhvk that there's a lot of valuable discussion here. But I think it's valuable discussion mostly about a different, future PR (possibly in an affiliated package). There's a ton of valuable stuff to do here, but all I think should be in this actual PR is classes to represent the magnitude-to-physical quantity ratio conversion. This gives users something specific to use when they want to represent magnitudes or anything that's flux-like (+).

Now to address @mhvk's

In particular, while Magnitude is subclassed from Quantity, no real use is made of it, while I think it could be extremely useful. In particular, actually using the unit could have real benefits.

I can see the use here, and nothing of what I propose is incompatible with adding that later. My point here (and I think @astrofrog was saying the same?) is that sometimes you don't know the system. But you still want a way to represent the mags and fluxes, and switch between them. I promise I would use that literally today if it existed, and I know of at least one person in an office down the hall who asked me if such a thing exists in astropy (without any prompting).

So I would say @pllim should remove the AB and STMAG part and any specific reference to systems, and leave the Magnitude and add the Flux (+) class in this PR. This is setting the groundwork for all the more in-depth ideas discussed here. Then, @pllim, you can just move the AB and STMAG classes to pysynphot for now. @astrofrog @aconley @kbarbary and @mhvk, does that sound like a way to move forward to you?

If so, then the remaining question is the zeropoint attribute. I think it really ought to be in Magnitude, because it would certainly make my life easier, and if so, it also makes sense to carry it in Flux (+). Remember that it also can just be None, which would mean there's not a valid way to convert between the two. Also remember that in any more advanced system, i.e., VegaMag, you would definitely want to override to_flux, and that could be much more complicated, but that option is lost if we mandate zeropoint as an argument to to_flux.

(+) While in the original conversation we were talking about flux density (energy / time /area / spectral_unit), in fact what I'm saying above is that this generalizes to "fluxish" meaning any of the things astronomers sometimes use magnitudes to represent the ratio of. That is sometimes [energy or counts], [energy or counts] / time, [energy or counts] / time / area, [energy or counts] / area, [energy or counts] / time / area / spectral_unit (I agree with @aconley that this is the one that has the unique "flux density" name). I don't care too much about the actual name, Flux was just meant as a colloquialism used in the field. We could also call it Fluxish if that makes everyone feel better. Or if someone has a shorter synonym for ThingThatMagsAreARatioOf, I'm all ears :)

@eteq - Indeed, one sometimes does not know the units and one still would want to transfer. So, in my proposed API, one would use dimensionless as the underlying unit or even define a new onet. My larger point, however, remains, that the API could be much better, and as such I'm not in favour of merging this.

I'll comb through this next week and re-implement Magnitude and add Flux (Fluxish? FluxLike? FluxYou? 😃 ). It does not have to go into 0.3. If the 6-month release time starts soon, pysynphot won't be released until the next astropy release cycle anyway.

I also want to point out that pysynphot already has its own units and conversion system. It does not need this PR to work. The idea is to move some of the more general units/conversion codes over to astropy so others can benefit from them. So if you think this is not useful to your work, I do not have to do this.

Just so we all know, this is definitely not going to make 0.3, as we decided on a feature freeze today in the interest of getting it out soon.

@mhvk

So, in my proposed API, one would use dimensionless as the underlying unit or even define a new one.

I don't understand how what @pllim has (especially with the modifications I suggest) is in conflict with your suggestion here. I'm saying we start with the simplest thing, and later on we can add support for additional systems. As far as I can tell, this could be done with units just as you suggest after this PR is merged. (I should add, I don't think we should use units to represent different magnitude systems for a few different reasons... But I think that belongs in a separate discussion, not this PR.)

My larger point, however, remains, that the API could be much better, and as such I'm not in favour of merging this.

I am advocating we use this as a starting point. I think it's terribly unwise to say we're going to wait until we've figured out the whole photometry system scheme. Long PRs with a lot of functionality in one have historically been much more work for everyone (and more prone to bugs), so I'm saying we should get the groundwork in place now, because even just that would be valuable right now for doing science in python.

(and 👍 to FluxYou ... ok, not really, but still funny, @pllim :)

EDIT: @mhvk - one thing that might clear things up: are you suggesting having only Magnitude, and not Flux? I see how that would be incompatible, but I'm not sure if that's actually what you're saying.

@eteq - you and @pllim certainly convinced me that we should keep photometric systems out of this PR, so let's indeed focus on a general magnitude.

For those, too, I feel we need to think carefully about API (am I turning into @demitri?). My particular problems with the current API, in order of increasing generality:

1. It introduces only two methods (from_flux and to_flux), both of which I find problematic:
   • Aesthetically, I find these very non-(astro)pythonesque (we were finally getting rid of the to_string stuff!)
   • from_flux is unnecessary (a quantity with the right units should just be converted; see my API sketch).
   • to_flux may be good to have (though hopefully mostly unnecessary with a good equivalency), but at least it should be flux -- it is just a different representation, not a conversion.
2. unit should be used; e.g., reddening and extinction are truly dimensionless, but measurements are not. And, hey, we might as well have AB, ST, instrumental, and "some unknown flux" magnitudes (there is an UnrecognizedUnit class that would be suitable just for that; we may also need a unit that has known physical type with fixed but unknown scale -- this corresponds to the undefined zeropoint case).
3. Basic operations should be defined: what happens for m1+m2, m1*m2??, m1.mean(), etc. This requires units, in that the outcome of, e.g., m1-m2 depends on the units of m1 and m2 (think extinction vs a comparison of two sources in the same image). Here, it may be better to contain a Quantity rather than subclass from it.
4. Most of the above problems would be true for other "functional units" such as Lupton's magnitude, dB, dex, etc., and thus to, e.g., logg, gains/losses in radio astronomy instrumentation, [Fe/H]. (and here also reality includes usage like dB/m...) It would be good if the work is general enough that it can be used there.

Of course, I also realise one has to start somewhere, but it seems to me at least 1--3 should be addressed. I also think making use of unit and equivalencies can solve many of these [1], but whether it is easy remains to be seen. If so, it would certainly solve (4) almost implicitly.

To me, it would seem the best way to proceed is first construct a requirement list for magnitudes.

[1] E.g., one can map m1-m2 to Magnitude(m1.flux/m2.flux); with units defined, this would do at least extinction and relative brightnesses right.

p.s. On Flux, strangely, I'm not sure the two are that coupled; I think it is totally fine to construct a Flux class first without worrying about magnitudes at all.

One further thing that probably is obvious to most, but struck me again only now: one often thinks of magnitudes as dimensionless, since they are -2.5*log(f1/f_ref), but that is not right: really they are fluxes in units of f_ref, of which one just happens to have taken -2.5*log(). This is no more unit-less thanlog(g.cgs.value) is unitless. One happy consequence is that the units of the fluxes described by AB magnitudes are just Unit('3630 Jy'). Of course, this is much more obvious for instrumental magnitudes, where the units are u.count/u.s (and more bizarre for vega, where they are Unit('vega')). Hence, the magnitude class does in fact not need a zeropoint; the unit suffices.

This would also suggest the Flux class does not need a zero point either, but rather the abiltity to have rather odd units, such as Unit('Vega').

@pllim, @eteq - I owe you both an apology -- I clearly had not looked carefully enough at the PR, it is much closer to what I envisioned than I thought. Clearly, I made too hurried a judgment. But at least this also meant it was relatively easy to implement all I thought would be nice to have, and hence I did a (very long) PR to your branch (pllim#1, or https://github.com/mhvk/astropy/tree/mag-math). It includes updated documentation and adapted tests (but not the many tests that will need to be added). With it, one can now do:

EDIT -- with bright daylight came a better way to represent magnitudes; now changed below

In [1]: import astropy.units as u; import astropy.magnitudes as m
In [2]: mag = -5. * m.mag

In [3]: mag, mag.flux
Out[3]: (<Magnitude -5.0 mag>, <Quantity 100.0 >)

In [4]: mag = -5. * m.ST

In [5]: mag, mag.flux, mag.to(u.erg/u.cm**2/u.s/u.AA)
Out[5]: 
(<Magnitude -5.0 mag(ST)>,
 <Quantity 100.0 ST>,
 <Quantity 3.6307805477010025e-07 erg / (Angstrom cm2 s)>)

In [6]: mag = -5. * m.AB

In [7]: mag, mag.flux, mag.to(u.kJy)
Out[7]: 
(<Magnitude -5.0 mag(AB)>,
 <Quantity 100.0 AB>,
 <Quantity 363.07805477010027 kJy>)

In [8]: mag + 5.*m.mag
Out[8]: <Magnitude 0.0 mag(AB)>

In [9]: mag - 1.*m.AB    # change to dimensionless automatically!
Out[9]: <Magnitude -6.0 mag>

In [10]: mag = m.Magnitude(1000.*u.count)

In [11]: mag, mag.flux
Out[11]: (<Magnitude -7.5 mag(ct)>, <Quantity 1000.0 ct>)

In [12]: zp = 0.*m.ST - m.Magnitude(1e10*u.count)

In [13]: mag += zp

In [14]: mag, mag.flux, mag.to(u.erg/u.cm**2/u.s/u.AA)
Out[14]: 
(<Magnitude 17.5 mag(ST)>,
 <Quantity 1e-07 ST>,
 <Quantity 3.6307805477010027e-16 erg / (Angstrom cm2 s)>)

By introducing u.AB and u.ST flux density units, it became much more obvious how to display magnitudes: just do <Magnitude <value> mag(<unit)> (see edited entry above). Even better, this made the PR cleaner too!

I thought it might be handy to provide a direct link to the documentation, so one can see how things should work:
http://www.astro.utoronto.ca/~mhvk/astropy/magnitudes/index.html
(Note: not all in-text links are right -- I'm still rather inexperienced with sphinx/doctext)

I (mostly) like what you've done here @mhvk! However, how will things work when we want to introduce magnitude systems like 'Vega'? I don't think we will be able to define a u.Vega, will we? (Note that I'm not saying that we should implement spectral magnitude systems now, but the API should be able to handle them when we do.)

My solution is to simply replace your MagUnits with MagSystems. The behavior would look identical:

>>> mag = m.Magnitude(u.Quantity(3.63e-9, u.erg / u.cm ** 2 / u.s / u.AA),
...                   system=m.ST)
>>> mag                                 
<Magnitude 0.0002334... ST>

But m.ST, m.AB would be instances of MagSystem or derived classes. These would be fairly simple classes that can tell you the zeropoint flux of the system (flux corresponding to mag = 0). Making them full fledged classes rather than trying to make them look like Units will give us the flexibility to add things like Vega in the future.

@kbarbary - Thanks! I think a MagSystem may work, but it would definitely be nice to have it at least be sub-classed from Unit. My hope with the current scheme is eventually to do things like,

dm = u.Unit('abs', 4.*np.pi*(10.*u.pc)**2)
DMmag = MagUnit(dm**0.5)
(10.*m.ST - av - 5*DMmag).to(u.erg/u.s/u.AA)

Note that in the unit-based scheme, one could still have Vega-based fluxes, but they might be in 'Vega' units, which are not directly convertible to anything else, except when one applied one's favourite Vega system equivalency. (This all makes use of unit's equivalency system, which in principle allows for arbitrary functions to convert from one to the other -- even a full pysynphot calculation ;-)

p.s. Just as a heads-up: I'm now working on my point (4), making this work for all logarithmic units rather than just magnitudes. The API is still similar, though the implementation somewhat different. The framework would be fairly easily extensible to other functional units.

@mhvk - I need to consider the matter a bit more, I think, but there one thing in your API I find confusing

In [12]: zp = 0._m.ST - m.Magnitude(1e10_u.count)

This seem opaque as to what's going on, and this is the use case I thought the API @pllim and I were talking about does with nicely; you have to actually type zeropoint=1e10 * u.count, rather than having to always carry around the zeropoint and remember to subtract it out. Or am I mis-interpreting how you're saying this should be used?

In [12]: zp = 0._m.ST - m.Magnitude(1e10_u.count)

This seem opaque as to what's going on, and this is the use case I thought the API @pllim and I were talking about does with nicely; you have to actually type zeropoint=1e10 * u.count, rather than having to always carry around the zeropoint and remember to subtract it out. Or am I mis-interpreting how you're saying this should be used?

I think my example may have been somewhat misleading. In it, I tried to show how my proposal could be used as part of photometric calibration, where at some point one has an observed count rate of some reference star with a known magnitude in the ST system. Then, the zero point of the calibration is:

zp = <refmag> * m.ST - m.Magnitude(<observed-ref-counts>)

Here, since the subtraction implies that units are divided, zp has the perhaps somewhat odd units of ST/ct. But if you add this to another instrumental magnitude,

starmag = m.Magnitude(<observed-star-counts>) + zp

then the multiplication of the physical units means that starmag has units of ST (or whatever other flux unit used for the reference system).

Returning now to the question of how to store system zeropoints (rather than do calibration), the difference from what @pllim and you suggested is that in my proposal the zero point is simply defined by the physical unit, rather than being a log value corresponding to it. So, for the m.ST example, rather thatn store -2.5*np.log10(3.63e-9)=21.1 and then still having to check that values given are in erg/cm^2/s/A, my logic has one just use a new unit, Unit(3.63e-9*erg/cm**2/s/AA), to which input quantities are transformed before taking the log. In other words, I store the flux corresponding to mag=0 instead of the magnitude corresponding to 1 * physical_unit.

@pllim - do you think this is superceded by #1894 (assuming it gets merged)? Or does this include some features that are not in #1894 that you want?

@eteq, you can close this. I can always build on the other PR as needed. Thanks.

Ok, thanks for the work on this, anyway, @pllim !

@pllim - thanks also from me!

You're welcome :)