Add CDR to the by component processing

f4e_radwaste/data_formats/data_mass.py

@@ -41,6 +41,20 @@ def get_mass_from_cells(self, cell_ids: List[int]) -> float:
  filtered_dataframe = self.get_filtered_dataframe(cells=cell_ids)
  return filtered_dataframe[KEY_MASS_GRAMS].sum()
 
+ def calculate_material_id_proportions(
+ self, cell_ids: List[List[int]]
+ ) -> List[pd.Series]:
+ mat_id_proportions = []
+
+ for comp_cell_ids in cell_ids:
+ df = self.get_filtered_dataframe(cells=comp_cell_ids)
+ masses_by_material = df.groupby(KEY_MATERIAL).sum()
+ mat_id_proportion = masses_by_material / masses_by_material.sum()
+ mat_id_proportion = mat_id_proportion[KEY_MASS_GRAMS]
+ mat_id_proportions.append(mat_id_proportion)
+
+ return mat_id_proportions
+
  @property
  def materials(self) -> np.ndarray:
  return self._dataframe.index.unique(level=KEY_MATERIAL).values

f4e_radwaste/main.py

@@ -2,6 +2,7 @@
 import shutil
 from pathlib import Path
 
+from f4e_radwaste.post_processing.calculate_dose_rates import DoseCalculator
 from f4e_radwaste.post_processing.folder_paths import FolderPaths
 from f4e_radwaste.post_processing.input_data import InputData
 from f4e_radwaste.post_processing.post_processing import (
@@ -14,6 +15,7 @@
  isotope_criteria_file,
  filter_cells_file,
 )
+from f4e_radwaste.readers.aux_material_file import read_element_mixes_of_materials
 from f4e_radwaste.readers.component_ids_file import (
  get_component_ids_from_folder,
  get_relevant_cells_from_components,
@@ -65,13 +67,22 @@ def by_component_process(input_folder_path: Path):
  folder_paths = get_folder_paths(input_folder_path)
  input_data = load_input_data_from_folder(folder_paths.input_files)
  component_ids = get_component_ids_from_folder(folder_paths.input_files)
+ dose_calculator = DoseCalculator(
+ dose_1_m_factors=read_dose_1_m_factors(PATH_TO_DOSE_FACTORS_FILE),
+ cdr_factors=read_contact_dose_rate_factors(PATH_TO_DOSE_FACTORS_FILE),
+ element_mix_by_material_id=read_element_mixes_of_materials(
+ folder_paths.input_files
+ ),
+ )
 
  # Filter in only the cells that appear in the components for performance reasons
  cells_to_include = get_relevant_cells_from_components(component_ids)
  input_data.apply_filter_include_cells(cells_to_include)
 
  # Process and save the data grouped by component in CSV
- process_input_data_by_components(input_data, folder_paths, component_ids)
+ process_input_data_by_components(
+ input_data, folder_paths, component_ids, dose_calculator
+ )
 
 
 def get_folder_paths(input_folder_path: Path) -> FolderPaths:
@@ -98,17 +109,13 @@ def load_input_data_from_folder(folder_path: Path) -> InputData:
  data_absolute_activity = dgs_file.read_file(folder_path / FILENAME_DGS_DATA)
  data_mesh_info = mesh_info_file.read_file(folder_path / FILENAME_MESHINFO)
  isotope_criteria = isotope_criteria_file.read_file(PATH_TO_CRITERIA_FILE)
- dose_1_m_factors = read_dose_1_m_factors(PATH_TO_DOSE_FACTORS_FILE)
- cdr_factors = read_contact_dose_rate_factors(PATH_TO_DOSE_FACTORS_FILE)
 
  return InputData(
  data_absolute_activity,
  data_mesh_info,
  isotope_criteria,
- dose_1_m_factors,
- cdr_factors,
  )
 
 
 if __name__ == "__main__":
- standard_process(Path(r"D:\WORK\tryingSimple\tests\old_data\ivvs_cart"))
+ by_component_process(Path(r"D:\WORK\tryingSimple\tests\old_data\ivvs_cart"))

f4e_radwaste/post_processing/calculate_dose_rates.py

@@ -1,59 +1,78 @@
 from dataclasses import dataclass
-from typing import Dict
+from typing import Dict, List
 
 import pandas as pd
 
-from f4e_radwaste.constants import KEY_DOSE_1_METER, KEY_CDR
+from f4e_radwaste.constants import (
+ KEY_DOSE_1_METER,
+ KEY_CDR,
+)
 from f4e_radwaste.data_formats.data_mesh_activity import DataMeshActivity
 
 
 @dataclass
 class DoseCalculator:
  dose_1_m_factors: pd.Series
  cdr_factors: pd.DataFrame
- material_mixes_by_id: Dict[int, pd.Series]
+ element_mix_by_material_id: Dict[int, pd.Series]
 
  def calculate_doses(
- self,
- data_mesh_activity: DataMeshActivity,
- material_proportion: Dict[int, float],
+ self, comp_activity: DataMeshActivity, cdr_factor_columns: List[pd.Series]
  ) -> DataMeshActivity:
- activity_df = data_mesh_activity.get_filtered_dataframe()
+ activity_df = comp_activity.get_filtered_dataframe()
 
  dose_1m_column = (activity_df * self.dose_1_m_factors).sum(axis=1)
 
- cdr_column = self._calculate_cdr_column(activity_df, material_proportion)
+ cdr_column = self._calculate_cdr_values(activity_df, cdr_factor_columns)
 
- updated_df = data_mesh_activity.get_dataframe_with_added_columns(
+ updated_df = comp_activity.get_dataframe_with_added_columns(
  {KEY_DOSE_1_METER: dose_1m_column, KEY_CDR: cdr_column}
  )
  return DataMeshActivity(updated_df)
 
- def _calculate_cdr_column(self, activity_df, material_proportion) -> pd.Series:
-  element_mix = self._calculate_element_mix_from_material_id_proportion(
-  material_proportion
-  )
- cdr_factors_mix = self._calculate_cdr_factors_for_element_mix(element_mix)
+ @staticmethod
+ def _calculate_cdr_values(
+ activity_df: pd.DataFrame, cdr_factor_columns: List[pd.Series]
+ ) -> pd.Series:
+ cdr_values = []
 
- return (activity_df * cdr_factors_mix).sum(axis=1)
+ for (_, row), cdr_factors in zip(activity_df.iterrows(), cdr_factor_columns):
+ cdr_values.append((row * cdr_factors).sum())
 
- def _calculate_element_mix_from_material_id_proportion(
- self, materials_proportion: Dict[int, float]
- ) -> pd.Series:
- proportioned_mixes = []
+ return pd.Series(index=activity_df.index, data=cdr_values)
 
-  material_ids = materials_proportion.keys()
- for material_id in material_ids:
-  material_mix = self.material_mixes_by_id[material_id]
-  proportion = materials_proportion[material_id]
+ def calculate_cdr_factors_list(
+ self, material_id_proportions: List[pd.Series]
+ ) -> List[pd.Series]:
+ element_mixes = self._calculate_element_mixes(material_id_proportions)
 
- proportioned_mixes.append(material_mix * proportion)
+ cdr_factors = []
+ for element_mix in element_mixes:
+ cdr_factors.append((self.cdr_factors * element_mix).sum(axis=1))
 
- element_mix = pd.concat(proportioned_mixes, axis=1)
- element_mix = element_mix.sum(axis=1)
- return element_mix
+ return cdr_factors
 
- def _calculate_cdr_factors_for_element_mix(
- self, element_mix: pd.Series
- ) -> pd.Series:
- return (self.cdr_factors * element_mix).sum(axis=1)
+ def _calculate_element_mixes(self, material_id_proportions: List[pd.Series]):
+ element_mixes = []
+
+ for mat_id_proportion in material_id_proportions:
+ proportioned_mixes = []
+
+ for mat_id, proportion in mat_id_proportion.items():
+ if mat_id not in self.element_mix_by_material_id:
+ continue
+
+ # noinspection PyTypeChecker
+ proportioned_mixes.append(
+ self.element_mix_by_material_id[mat_id] * proportion
+ )
+
+ if len(proportioned_mixes) == 0:
+ element_mixes.append(pd.Series())
+ continue
+
+ element_mix = pd.concat(proportioned_mixes, axis=1)
+ element_mix = element_mix.sum(axis=1)
+ element_mixes.append(element_mix)
+
+ return element_mixes

f4e_radwaste/post_processing/components_info.py

@@ -0,0 +1,20 @@
+from typing import List
+
+from f4e_radwaste.data_formats.data_mass import DataMass
+from f4e_radwaste.post_processing.calculate_dose_rates import DoseCalculator
+
+
+class ComponentsInfo:
+ def __init__(
+ self,
+ component_ids: List[List],
+ data_mass: DataMass,
+ dose_calculator: DoseCalculator,
+ ):
+ self.names, self.cell_ids = zip(*component_ids)
+
+ mat_id_proportions = data_mass.calculate_material_id_proportions(self.cell_ids)
+
+ self.cdr_factors = dose_calculator.calculate_cdr_factors_list(
+ material_id_proportions=mat_id_proportions
+ )

f4e_radwaste/post_processing/input_data.py

@@ -15,8 +15,10 @@
 from f4e_radwaste.data_formats.data_mesh_activity import DataMeshActivity
 from f4e_radwaste.data_formats.data_mesh_info import DataMeshInfo
 from f4e_radwaste.helpers import format_time_seconds_to_str
+from f4e_radwaste.post_processing.calculate_dose_rates import DoseCalculator
 from f4e_radwaste.post_processing.classify_waste import classify_waste
 from f4e_radwaste.post_processing.component_output import ComponentOutput
+from f4e_radwaste.post_processing.components_info import ComponentsInfo
 from f4e_radwaste.post_processing.folder_paths import FolderPaths
 from f4e_radwaste.post_processing.mesh_ouput import MeshOutput
 
@@ -89,14 +91,30 @@ def get_mesh_activity_by_time_and_materials(
  return DataMeshActivity(voxel_activity_dataframe)
 
  def get_component_output_by_time_and_ids(
- self, decay_time: float, component_ids: List[List]
+ self,
+ decay_time: float,
+ component_ids: List[List],
+ dose_calculator: DoseCalculator,
  ) -> ComponentOutput:
  comp_mesh_activity = self.get_component_mesh_activity_by_time_and_ids(
  decay_time=decay_time, component_ids=component_ids
  )
 
  comp_mesh_activity = classify_waste(comp_mesh_activity, self.isotope_criteria)
 
+ ComponentsInfo(component_ids, self.data_mesh_info.data_mass, dose_calculator)
+ dose_calculator.calculate_doses(
+ comp_activity=comp_mesh_activity,
+ cdr_factor_columns=)
+
+ # data_mass = self.data_mesh_info.data_mass
+ # material_proportions = data_mass.calculate_material_proportions(
+ # component_ids=component_ids
+ # )
+ # comp_mesh_activity = dose_calculator.calculate_doses(
+ # comp_mesh_activity, material_proportions
+ # )
+
  return ComponentOutput(
  name=f"{format_time_seconds_to_str(decay_time)}_by_component",
  data_mesh_activity=comp_mesh_activity,

f4e_radwaste/post_processing/post_processing.py

@@ -1,5 +1,6 @@
 from typing import List
 
+from f4e_radwaste.post_processing.calculate_dose_rates import DoseCalculator
 from f4e_radwaste.post_processing.folder_paths import FolderPaths
 from f4e_radwaste.post_processing.input_data import InputData
 
@@ -24,13 +25,18 @@ def process_input_data_by_material(input_data: InputData, folder_paths: FolderPa
 
 
 def process_input_data_by_components(
- input_data: InputData, folder_paths: FolderPaths, component_ids: List[List]
+ input_data: InputData,
+ folder_paths: FolderPaths,
+ component_ids: List[List],
+ dose_calculator: DoseCalculator,
 ):
  decay_times = input_data.data_absolute_activity.decay_times
 
  for decay_time in decay_times:
  component_output = input_data.get_component_output_by_time_and_ids(
- decay_time=decay_time, component_ids=component_ids
+ decay_time=decay_time,
+ component_ids=component_ids,
+ dose_calculator=dose_calculator,
  )
 
  component_output.save(folder_paths)

tests/test_data_formats/test_data_mass.py

@@ -78,3 +78,16 @@ def test_get_mass_from_cells(self):
 
  def test_materials(self):
  np.testing.assert_array_equal(self.data_mass.materials, np.array([10, 20, 40]))
+
+ def test_calculate_material_id_proportions(self):
+ result = self.data_mass.calculate_material_id_proportions(
+ cell_ids=[[11, 12], [14]]
+ )
+
+ first_comp_mat_10_proportion = result[0][10]
+ last_comp_mat_40_proportion = result[-1][40]
+
+ self.assertAlmostEqual(
+ (2.34 + 1.09) / (2.34 + 1.09 + 3.13), first_comp_mat_10_proportion
+ )
+ self.assertAlmostEqual(1, last_comp_mat_40_proportion)