add mcst unit support (#184)

* add mcst test case

.github/actions/get-testdata-all/action.yml

@@ -103,3 +103,8 @@ runs:
         outputfolder: ${{ inputs.outputfolder }}
         filename: minimal_flamingo_snapshot.hdf5
         url: https://heibox.uni-heidelberg.de/f/8749a46fdfc74577b29a/?dl=1
+    - uses: ./.github/actions/get-testdata
+      with:
+        outputfolder: ${{ inputs.outputfolder }}
+        filename: minimal_mcst_arepo_mcstmodel.hdf5
+        url: https://heibox.uni-heidelberg.de/f/802d5c8cbee1483a983f/?dl=1

src/scida/configfiles/units/gadget_cosmological.yaml

@@ -43,17 +43,28 @@ fields:
     GFM_MetalsTagged:
     GFM_WindDMVelDisp: km / s
     GFM_WindHostHaloMass: code_mass
+    GrackleTemperature: K
+    GrackleCoolTime: (kpc / h) / (km / s)  # double checked that physical kpc
     HighResGasMass: code_mass
+    HeIIIMassFraction:  # no units
+    HeIIMassFraction:  # no units
+    HeIMassFraction:  # no units
+    HIIMassFraction:  # no units
+    HIMassFraction:  # no units
     InternalEnergy: (km / s) ^ 2
     InternalEnergyOld: (km / s)^2
     Machnumber:
     MagneticField: (h / a^2) * (code_pressure)^(1/2)
     MagneticFieldDivergence: (h^3 / a^2) * code_mass^(1/2) * (km/s) * (ckpc)^(-5/2)
     # we guess the alternative has the same units as the original field...
     MagneticFieldDivergenceAlternative: (h^3 / a^2) * code_mass^(1/2) * (km/s) * (ckpc)^(-5/2)
-    NeutralHydrogenAbundance:
+    Metallicity:  # no units
+    NeutralHydrogenAbundance:  # no units
     # TODO: Pressure?
+    RadiationEnergyDensity: code_mass * (km / s)^2 / (kpc / h)^3  # double checked that physical kpc
     SmoothingLength: code_length
+    StromgrenRadius: code_length
+    StromgrenSourceID: none
     SubfindDMDensity: code_mass / code_length^3
     Temperature: K
     TimebinHydro:
@@ -62,6 +73,8 @@ fields:
     Acceleration: unknown  # TODO: From MTNG
   PartType1:
     SubfindDMDensity: code_mass / code_length^3
+  PartType2:
+    Metallicity: # field existence is a bug in MCST model, needs to be removed.
   PartType3:
     ParentID: none
     TracerID: none
@@ -75,7 +88,13 @@ fields:
     GFM_MetalsTagged:
     GFM_StellarFormationTime: # expressed as signed scale factor
     GFM_StellarPhotometrics:  # not supported mag for now
-    StellarHsml:
+    IMFMass: code_mass
+    IonisingPhotonRate1e49: 1e49 / s
+    LowMass: code_mass
+    Metallicity: # no units
+    StellarArray: none
+    StellarHsml: code_length
+    StellarLuminosity: 1e10 * Msun * (km^3 / s^3) / kpc
     SubfindDMDensity: code_mass / code_length^3
   PartType5:
     BH_BPressure: (h^4 / a^4) * (1e10 Msun) * (km / s)^2 / ckpc^3

tests/test_discovery.py

@@ -15,7 +15,7 @@
 from scida.io import load_metadata
 from scida.misc import check_config_for_dataset
 from scida.series import DatasetSeries
-from tests.testdata_properties import require_testdata_path
+from tests.testdata_properties import get_testdata_yaml_config, require_testdata_path
 
 
 def return_intended_stype(name) -> Type[DatasetSeries]:
@@ -42,7 +42,7 @@ def return_intended_dstype(name, simconf=False) -> Type[Dataset]:
             return MTNGArepoSnapshot
     elif any([k in name for k in ["tng-cluster", "tngcluster"]]):
         return TNGClusterSnapshot
-    elif any(k in name for k in ["tng", "illustris", "auriga"]):
+    elif any(k in name for k in ["tng", "illustris", "auriga", "mcst"]):
         if "group" in name:
             return ArepoCatalog
         else:
@@ -65,6 +65,19 @@ def return_intended_dstype(name, simconf=False) -> Type[Dataset]:
     elif any(k in name for k in ["rockstar"]):
         return RockstarCatalog
 
+    # alternative check testdata.yaml metadata
+    tdprops = get_testdata_yaml_config()
+    dname = name.split("/")[-1]
+    for k in tdprops:
+        # remove extension if present
+        dname_without_ext = dname
+        if dname.endswith(".hdf5"):
+            dname_without_ext = dname[:-5]
+        fname = tdprops[k].get("fn", "")
+        if str(k).lower() == dname_without_ext or fname == dname:
+            if "areposnapshot" in tdprops[k]["types"]:
+                return ArepoSnapshot
+
     print("No intended type found for %s" % name)
     return Dataset
     # raise ValueError("Have not specified intended type for %s" % name)

tests/testdata.yaml

@@ -209,11 +209,17 @@ testdata:
     fn: minimal_TNG50-1_group_z3.hdf5
 
   TNG50-1_group_z0_minimal:
-      types:
-          - interface
-          - areposnapshot
-          - areposnapshot_withcatalog|G|1|2
-      fn: minimal_TNG50-1_group_z0.hdf5
+    types:
+      - interface
+      - areposnapshot
+      - areposnapshot_withcatalog|G|1|2
+    fn: minimal_TNG50-1_group_z0.hdf5
+
+  MCST_snapshot_z3_minimal:
+    types:
+      - interface
+      - areposnapshot
+    fn: minimal_mcst_arepo_mcstmodel.hdf5
 
   FIRE2_snapshot_z1_minimal:
     types:

tests/testdata_properties.py

@@ -57,6 +57,10 @@ def add_testdata_entry(name, types=None, marks=None, fn=None):
         print(name, properties)
 
 
+def get_testdata_yaml_config():
+    return testdata_properties
+
+
 testdataskip = testdataskip.split()
 testdata_local = []
 for k, v in testdatadict.items():