amuse icon indicating copy to clipboard operation
amuse copied to clipboard

Astropy units compatibility / conversion?

Open rieder opened this issue 8 years ago • 23 comments

Since Astropy has a unit systems different from the one in AMUSE while it's quite plausible both packages are used in the same project, it seems practical to make sure these unit systems are in some way compatible. The implementations are somewhat different (e.g. Astropy uses the multiplication operator '*' to assign units, where Amuse uses '|', among others), but not completely.

I think it would be useful to check:

  • if we use the same definition of units (we really should)
  • if we can make an easy conversion module ('astropy_to_amuse' and VV) to replace one unit type with the other for better compatibility
  • if it would make any sense at all to switch unit systems (this would probably be a lot of work, so I doubt it's worthwhile unless we do a big overhaul)

A conversion module should be relatively straightforward, using a table of Amuse units and their Astropy equivalent.

rieder avatar Oct 18 '17 07:10 rieder

Simple concept for a conversion script (hardcoded only for meters):

import astropy.units as astropyunits
from amuse.units import units as amuseunits


def amuse_equivalent_unit(
        astropyunit,
        ):
    # Check if the named unit exists directly in Amuse
    # TODO: read from a lookup table
    if astropyunit == astropyunits.m:
        return amuseunits.m
    # TODO: Throw Astropy-style exception if the unit can't be found
    return -1


def astropy_equivalent_unit(
        amuseunit,
        ):
    # Check if the named unit exists directly in Astropy
    # TODO: read from a lookup table
    if amuseunit == amuseunits.m:
        return astropyunits.m
    # TODO: Throw AMUSE-style exception if the unit can't be found
    return -1


def astropy_to_amuse(
        astropyquantity,
        ):
    astropyunit = astropyquantity.unit
    value = astropyquantity.value
    amuseunit = amuse_equivalent_unit(astropyunit)
    amusequantity = value | amuseunit
    return amusequantity


def amuse_to_astropy(
        amusequantity,
        ):
    amuseunit = amusequantity.unit
    value = amusequantity.value_in(amuseunit)
    astropyunit = astropy_equivalent_unit(amuseunit)
    astropyquantity = value * astropyunit
    return astropyquantity

rieder avatar Oct 18 '17 07:10 rieder

I think switching unit systems would not work (mainly for the nbody units etc)..in any case it would be another heavy prerequisite (remember that we also want the framework to be more domain agnostic)...what is doable is to maybe accept astropy units as input, but that might need some cleanup/refactoring of the code for the data stores to minimize the adaptations...step two would be to output the same unit system as was used in the input...

ipelupessy avatar Oct 18 '17 14:10 ipelupessy

I agree switching over probably won't work, at least it's not quite worth the huge effort.

rieder avatar Oct 19 '17 03:10 rieder

👋 I was linked to this thread by Nora L. I'm a deputy maintainer of astropy.units, so this is something I care about! I think there are some things we can do on both of our ends to make the typical user experience a bit more painless (for example, see discussion here about interoperability between astropy.units and unyt). In particular, it would be nice to think about whether we can agree on a few high-level API decisions that would allow us to, e.g., initialize our respective unit objects with one from the (or any) other units package. For example, ideally astropyunits.Unit(amuseunits.some_unit) should work. Is this something you would be open to discussing at some point?

Also cc @mhvk, the lead maintainer of astropy.units.

adrn avatar Jul 31 '19 15:07 adrn

Thanks for responding @adrn! Yes, I would be very happy to discuss this at some point. @ipelupessy should also be included in that discussion as he has more knowledge about the AMUSE unit system than I do.

rieder avatar Jul 31 '19 15:07 rieder

also available...

ipelupessy avatar Aug 01 '19 09:08 ipelupessy

Am on vacation, so reading e-mail randomly, but units at least seems relatively straightforward as it would seem one should be able to go via the string representation, perhaps accessible via some standard dunder method, like __unit__ (for amuse, I'd think this would even work for all the IAU units, but one could envision also having a standard way to get the SI version).

Quantities may be harder. Does amuse support arrays of values with a unit?

mhvk avatar Aug 01 '19 12:08 mhvk

Am on vacation, so reading e-mail randomly, but units at least seems relatively straightforward as it would seem one should be able to go via the string representation, perhaps accessible via some standard dunder method, like unit (for amuse, I'd think this would even work for all the IAU units, but one could envision also having a standard way to get the SI version).

One problem we'd have to fix is that AMUSE in some cases uses slightly different values for IAU units than Astropy does (that's #171). But for the SI route, this should indeed probably be straightforward.

Quantities may be harder. Does amuse support arrays of values with a unit?

Yes, that's supported.

rieder avatar Aug 01 '19 12:08 rieder

Assigned this to the AMUSE 13 milestone, as I hope that by that version we can have some conversion process in place from the AMUSE side, perhaps similar to the method we currently use for representation of quantities (some_quantity.in_(units.parsec), some_quantity.value_in(units.parsec), some_quantity.value_in(units.parsec)). Not sure if we can use these methods directly (that would be ideal), but at least we could have a method that works similarly (some_quantity.in_external_units(astropyunits.parsec)) (this looks very ugly, need to think of something better).

rieder avatar Aug 20 '19 14:08 rieder

Reassigning to AMUSE 14.

rieder avatar Dec 14 '19 15:12 rieder

Here's a working conversion function from AMUSE -> Astropy. Two things to note:

  • it always returns to base (SI) units (not really a problem)
  • the values can/will change slightly, since AMUSE and Astropy use slightly different constants (#171) (this is a problem)
from amuse.units import units, constants
from amuse.units import si
from astropy import units as apunits
from astropy import constants as apconstants

def to_astropy(quantity):
    # Find the SI base of the unit
    unit = quantity.unit
    unit_bases = unit.base

    # Find the quantity's value in base units
    value = quantity.value_in(unit.base_unit())
    
    # Find the corresponding base unit in astropy
    apquantity = value
    for base_unit in unit_bases:
        if base_unit[1] == si.m:
            apquantity = apquantity * apunits.m**base_unit[0]
        elif base_unit[1] == si.kg:
            apquantity = apquantity * apunits.kg**base_unit[0]
        elif base_unit[1] == si.s:
            apquantity = apquantity * apunits.s**base_unit[0]
        elif base_unit[1] == si.A:
            apquantity = apquantity * apunits.A**base_unit[0]
        elif base_unit[1] == si.K:
            apquantity = apquantity * apunits.K**base_unit[0]
        elif base_unit[1] == si.mol:
            apquantity = apquantity * apunits.mol**base_unit[0]
        elif base_unit[1] == si.cd:
            apquantity = apquantity * apunits.cd**base_unit[0]

    return apquantity

rieder avatar Dec 14 '19 15:12 rieder

  • the values can/will change slightly, since AMUSE and Astropy use slightly different constants (#171) (this is a problem)

Or perhaps this is not so much a problem (since in the SI base AMUSE and Astropy agree)? @ipelupessy what do you think?

rieder avatar Dec 14 '19 15:12 rieder

I think this can be added to the units; about the values: do you have an example of the difference?

ipelupessy avatar Dec 15 '19 18:12 ipelupessy

Looks OK, though my sense is that going through the string representation is simpler (and not much slower). Faster would be to use that you implement bases quite similarly from how astropy does it, and feed bases and powers directly to a u.CompositeUnit.

Also, if your units are hashable (ours are), then one could just create a dict mapping one set of SI units to another.

mhvk avatar Dec 17 '19 01:12 mhvk

Ok, here are to/from functions that can be added to AMUSE to convert between AMUSE/Astropy:

def from_astropy(ap_quantity):
    "Convert a unit from Astropy to AMUSE"

    # Find SI bases of the unit
    si_bases = ap_quantity.si.unit.bases
    si_powers = ap_quantity.si.unit.powers
    si_units = list(zip(si_powers, si_bases))

    # Find the quantity's value in base units
    si_value = ap_quantity.si.value

    # Reconstruct the quantity in AMUSE units
    amuse_quantity = si_value
    for base_unit in si_units:
        if base_unit[1] == apu.m:
            amuse_quantity = amuse_quantity * (1 | units.m**base_unit[0])
        elif base_unit[1] == apu.kg:
            amuse_quantity = amuse_quantity * (1 | units.kg**base_unit[0])
        elif base_unit[1] == apu.s:
            amuse_quantity = amuse_quantity * (1 | units.s**base_unit[0])
        elif base_unit[1] == apu.A:
            amuse_quantity = amuse_quantity * (1 | units.A**base_unit[0])
        elif base_unit[1] == apu.K:
            amuse_quantity = amuse_quantity * (1 | units.K**base_unit[0])
        elif base_unit[1] == apu.mol:
            amuse_quantity = amuse_quantity * (1 | units.mol**base_unit[0])
        elif base_unit[1] == apu.cd:
            amuse_quantity = amuse_quantity * (1 | units.cd**base_unit[0])

    return amuse_quantity


def to_astropy(quantity):
    "Convert a unit from AMUSE to Astropy"

    # Find the SI bases of the unit
    unit = quantity.unit
    unit_bases = unit.base

    # Find the quantity's value in base units
    value = quantity.value_in(unit.base_unit())
    
    # Reconstruct the quantity in Astropy units
    ap_quantity = value
    for base_unit in unit_bases:
        if base_unit[1] == m:
            ap_quantity = ap_quantity * apu.m**base_unit[0]
        elif base_unit[1] == kg:
            ap_quantity = ap_quantity * apu.kg**base_unit[0]
        elif base_unit[1] == s:
            ap_quantity = ap_quantity * apu.s**base_unit[0]
        elif base_unit[1] == A:
            ap_quantity = ap_quantity * apu.A**base_unit[0]
        elif base_unit[1] == K:
            ap_quantity = ap_quantity * apu.K**base_unit[0]
        elif base_unit[1] == mol:
            ap_quantity = ap_quantity * apu.mol**base_unit[0]
        elif base_unit[1] == cd:
            ap_quantity = ap_quantity * apu.cd**base_unit[0]

    return ap_quantity

rieder avatar Apr 23 '20 12:04 rieder

I suggest I add them to a new 'amuse.units.astropy' file

rieder avatar Apr 23 '20 12:04 rieder

This is another Astropy attempt at adding compatibility: https://github.com/astropy/astropy/pull/10205

rieder avatar Jun 02 '20 11:06 rieder

Closing this issue in favour of #687 which probably requires a broader solution.

rieder avatar Oct 27 '20 11:10 rieder

When astropy is installed, AMUSE quantities now have an "as_astropy_quantity" method that converts to astropy.

rieder avatar Nov 11 '20 10:11 rieder

What was the conclusion, though, for things like Msun? Just convert to Msun or go through kg and take into account possibly different definitions?

mhvk avatar Nov 11 '20 13:11 mhvk

What was the conclusion, though, for things like Msun? Just convert to Msun or go through kg and take into account possibly different definitions?

I solved this by going through SI units like kg. Those will always be the same so we won't get any conversion errors. Currently, an AMUSE quantity will be converted to an Astropy quantity with an SI unit (it doesn't get converted back to the equivalent unit). We may add this later but it seems of less importance.

rieder avatar Nov 12 '20 12:11 rieder

OK, always going through SI will at least be consistent, though I think it will also lead to confusion. But probably unsolvable... Though with the IAU definitions of those types of numbers as constants, it will hopefully eventually go away.

mhvk avatar Nov 12 '20 12:11 mhvk

This issue has been automatically marked as stale because it has not had recent activity. It will be closed in 14 days if no further activity occurs. Thank you for your contributions.

stale[bot] avatar Mar 04 '22 15:03 stale[bot]