Merge pull request #285 from MC-kit/devel

convert gamma rates to photon/s

pyproject.toml

@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "xpypact"
-version = "0.5.2"
+version = "0.5.3"
 description = "\"Python workflow framework for FISPACT.\""
 authors = ["dvp <[email protected]>"]
 license = "MIT"

src/xpypact/dao/duckdb/implementation.py

@@ -6,6 +6,8 @@
 from dataclasses import dataclass
 from pathlib import Path
 
+import numpy as np
+
 import msgspec as ms
 
 from xpypact.dao import DataAccessInterface
@@ -89,7 +91,6 @@ def save(self, inventory: Inventory, material_id=1, case_id=1) -> None:
         _save_nuclides(cursor, inventory)
         _save_time_steps(cursor, inventory, material_id, case_id)
         _save_time_step_nuclides(cursor, inventory, material_id, case_id)
-        _save_gbins(cursor, inventory)
         _save_gamma(cursor, inventory, material_id, case_id)
 
     def load_rundata(self) -> db.DuckDBPyRelation:
@@ -268,27 +269,26 @@ def _save_time_step_nuclides(
 
 
 # noinspection SqlNoDataSourceInspection
-def _save_gbins(cursor: db.DuckDBPyConnection, inventory: Inventory) -> None:
+def _save_gamma(cursor: db.DuckDBPyConnection, inventory: Inventory, material_id=1, case_id=1):
     gs = inventory[0].gamma_spectrum
     if gs is None:
-        return
+        return  # pragma: no coverage
     sql = """
         insert into gbins values(?, ?);
     """
-    cursor.executemany(sql, enumerate(gs.boundaries))
-
-
-# noinspection SqlNoDataSourceInspection
-def _save_gamma(cursor: db.DuckDBPyConnection, inventory: Inventory, material_id=1, case_id=1):
+    boundaries = np.asarray(gs.boundaries, dtype=float)
+    cursor.executemany(sql, enumerate(boundaries))
     sql = """
         insert into timestep_gamma values(?, ?, ?, ?, ?);
     """
+    mids = 0.5 * (boundaries[:-1] + boundaries[1:])
     cursor.executemany(
         sql,
         (
-            (material_id, case_id, t.number, x[0] + 1, x[1][1])  # timestep number  # g  # rate
+            (material_id, case_id, t.number, x[0] + 1, x[1])  # use gbins index for upper boundary
             for t in inventory.inventory_data
             if t.gamma_spectrum
-            for x in enumerate(zip(t.gamma_spectrum.boundaries[1:], t.gamma_spectrum.values))
+            for x in enumerate(np.asarray(t.gamma_spectrum.values, dtype=float) / mids)
+            # convert rate MeV/s -> photon/s
         ),
     )

tests/test_duckdb_dao.py

@@ -59,6 +59,8 @@ def test_save(inventory_with_gamma) -> None:
             ],
         )
         assert not time_step_nuclides.loc[2, 290630].empty
+        gbins = dao.load_gbins().df().set_index("g")
+        assert gbins.loc[0].boundary == pytest.approx(1e-11)
         gamma = dao.load_gamma().df()
         assert not gamma.empty
         gamma = gamma.set_index(["time_step_number", "g"])