Skip to content

Commit

Permalink
fix cosmo unit discovery
Browse files Browse the repository at this point in the history
  • Loading branch information
cbyrohl committed Jan 6, 2025
1 parent 1f722d4 commit 438782a
Show file tree
Hide file tree
Showing 4 changed files with 29 additions and 3 deletions.
3 changes: 3 additions & 0 deletions src/scida/configfiles/simulations.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -121,15 +121,18 @@ data:
- Illustris100
identifiers:
name_contains: "Illustris-1" # can only rely on path name in lack of metadata...
unitfile: [ units/illustris.yaml ] # cosmology explicitly needed for catalogs

Illustris-2:
<<: *illustris_params
identifiers:
name_contains: "Illustris-2"
unitfile: [ units/illustris.yaml ]
Illustris-3:
<<: *illustris_params
identifiers:
name_contains: "Illustris-3"
unitfile: [ units/illustris.yaml ]

EAGLEtype:
identifiers:
Expand Down
2 changes: 2 additions & 0 deletions src/scida/customs/arepo/dataset.py
Original file line number Diff line number Diff line change
Expand Up @@ -175,6 +175,8 @@ def load_catalog(
fileprefix=prfx,
units=self.withunits,
ureg=ureg,
hints=self.hints,
metadata_raw_parent=self._metadata_raw,
)
if "Redshift" in self.catalog.header and "Redshift" in self.header:
z_catalog = self.catalog.header["Redshift"]
Expand Down
17 changes: 16 additions & 1 deletion src/scida/interfaces/mixins/cosmology.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,19 +28,30 @@ def __init__(self, *args, **kwargs):
self._mixins.append(self._mixin_name)
else:
self._mixins = [self._mixin_name]
self.hints["cosmological"] = True

# now call the parent class constructor as we need the metadata_raw
super().__init__(*args, **kwargs)

metadata_raw = self._metadata_raw
c = get_cosmology_from_rawmetadata(metadata_raw)
# sometimes, we want to inherit the cosmology from the parent dataset (e.g. catalogs for snapshots)
if c is None and "metadata_raw_parent" in kwargs:
c = get_cosmology_from_rawmetadata(kwargs["metadata_raw_parent"])
self.cosmology = c
if c is not None:
self.hints["cosmology"] = c
z = get_redshift_from_rawmetadata(metadata_raw)
self.hints["cosmological"] = True
self.redshift = z
self.metadata["redshift"] = self.redshift

if hasattr(self, "ureg"):
ureg = self.ureg
with ignore_warn_redef(ureg):
if c is not None:
ureg.define("h = %s" % str(c.h))
elif "cosmology" in self.hints:
ureg.define("h = %s" % str(self.hints["cosmology"].h))
if z is not None:
a = 1.0 / (1.0 + z)
ureg.define("a = %s" % str(float(a)))
Expand Down Expand Up @@ -97,6 +108,10 @@ def validate(cls, metadata: dict, *args, **kwargs):
if np.isclose(time, 1.0 / (1.0 + redshift)):
# print("Legacy metadata detected for Cosmology Mixin detection.")
valid = True
# sometimes, we have a cosmological runs if we only have the redshift key but no time key (e.g. LGalaxies)
if "/Header" in metadata and "Redshift" in metadata["/Header"]:
if "Time" not in metadata["/Header"]:
valid = True

# in case of SWIFT, we have the "Cosmological run" key in the "Cosmology" group set to 1
if "/Cosmology" in metadata:
Expand Down
10 changes: 8 additions & 2 deletions src/scida/interfaces/mixins/units.py
Original file line number Diff line number Diff line change
Expand Up @@ -262,7 +262,6 @@ def _get_default_units(mode: str, ureg: pint.UnitRegistry) -> dict:
"velocity": ureg.cm / ureg.s,
"time": ureg.s,
}

# common factors in cosmological simulations
if "h" in ureg:
udict["h"] = str_to_unit("h", ureg)
Expand Down Expand Up @@ -294,6 +293,10 @@ def new_unitregistry() -> UnitRegistry:
"""
ureg = UnitRegistry(autoconvert_offset_to_baseunit=True)
ureg.define("unknown = 1.0")
# we remove the hours "h" unit, so we cannot confuse it with the Hubble factor
del ureg._units["h"]
# we remove the year "a" unit, so we cannot confuse it with the scale factor
del ureg._units["a"]
return ureg


Expand Down Expand Up @@ -667,7 +670,10 @@ def check_unit_mismatch(unit, unit_metadata, override=False, path="", logger=log
bool:
Whether the units agree.
"""
without_units = unit == "none" or unit_metadata == "none" # explicitly no units
without_units = isinstance(unit, str) and unit == "none"
without_units |= (
isinstance(unit_metadata, str) and unit_metadata == "none"
) # explicitly no units
if unit is not None and unit_metadata is not None:
msg = None
if without_units and unit == unit_metadata:
Expand Down

0 comments on commit 438782a

Please sign in to comment.