diff --git a/opengate/engines.py b/opengate/engines.py index 41b2970c4..ef6e02ff2 100644 --- a/opengate/engines.py +++ b/opengate/engines.py @@ -5,6 +5,7 @@ import weakref from box import Box from anytree import PreOrderIter +import numpy as np import opengate_core as g4 from .exception import fatal, warning, GateImplementationError @@ -415,6 +416,35 @@ def initialize_optical_material_properties(self): for ( vol ) in self.simulation_engine.simulation.volume_manager.volumes.values(): + if hasattr(vol, "voxel_materials"): + # True if voxelized volume (ImageVolume) + materials = [m[2] for m in vol.voxel_materials] + materials.append(vol.material) + print(materials) + for i in materials: + # Read all materials contained in voxel_materials + for mat in g4.G4Material.GetMaterialTable: + # Loop on all materials created in the simulation + if str(mat.GetName()) == i: + # If material is in voxel_materials + mat_prop = load_optical_properties_from_xml( + self.physics_manager.optical_properties_file, + mat.GetName(), + ) + if mat_prop is not None: + self.g4_optical_material_tables[ + str(mat.GetName()) + ] = create_g4_optical_properties_table(mat_prop) + mat.SetMaterialPropertiesTable( + self.g4_optical_material_tables[ + str(mat.GetName()) + ] + ) + else: + self.simulation_engine.simulation.warn_user( + f"Could not load the optical material properties for material {mat.GetName()} " + f"found in volume {vol.name} from file {self.physics_manager.optical_properties_file}." + ) material_name = vol.g4_material.GetName() material_properties = load_optical_properties_from_xml( self.physics_manager.optical_properties_file, material_name diff --git a/opengate/tests/data b/opengate/tests/data index 0200198b2..be81811c9 160000 --- a/opengate/tests/data +++ b/opengate/tests/data @@ -1 +1 @@ -Subproject commit 0200198b2022f89cff342493a09dff73e5669d6b +Subproject commit be81811c9514bfcbff5014987b587f96ca6774e1 diff --git a/opengate/tests/src/test092_optic_material.py b/opengate/tests/src/test092_optic_material.py new file mode 100755 index 000000000..c4679ee0c --- /dev/null +++ b/opengate/tests/src/test092_optic_material.py @@ -0,0 +1,99 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- +""" +Created on Wed Apr 2 09:48:30 2025 + +@author: fava +""" +import opengate as gate +from opengate.utility import g4_units +import opengate.tests.utility as utility +import opengate_core as g4 +import uproot +import numpy as np +from scipy.spatial.transform import Rotation +import itk + + +if __name__ == "__main__": + + paths = utility.get_default_test_paths(__file__, output_folder="test092") + + # Define the Simulation object + sim = gate.Simulation() + + # Add a material database + sim.volume_manager.add_material_database(paths.data / "GateMaterials.db") + sim.physics_manager.special_physics_constructors.G4OpticalPhysics = True + paths.test_data = paths.data / "test092" + sim.physics_manager.optical_properties_file = paths.test_data / "Materials.xml" + + # Define the units used in the simulation set-up + m = gate.g4_units.m + cm = gate.g4_units.cm + keV = gate.g4_units.keV + eV = gate.g4_units.eV + Bq = gate.g4_units.Bq + mm = gate.g4_units.mm + sec = gate.g4_units.second + gcm3 = gate.g4_units.g_cm3 + deg = g4_units.deg + + sim.visu = False + sim.visu_type = "qt" + sim.number_of_threads = 1 + sim.random_seed = 123456789 + + world = sim.world + world.size = [3 * m, 3 * m, 3 * m] + world.material = "Air" + # world.material = "G4_Galactic" + # world.material = "G4_AIR" + + img = itk.imread(paths.test_data / "vox_volume.mhd") + spacing = img.GetSpacing() + size = np.array(img.GetLargestPossibleRegion().GetSize()) + center = size * spacing / 2 + + patient = sim.add_volume("Image", name="patient") + patient.image = paths.test_data / "vox_volume.mhd" # CT image + # patient.material = "Air" + patient.material = "G4_AIR" + patient.dump_label_image = paths.test_data / "labels_vox_volume.mhd" + patient.translation = [0, 0, 0] + patient.voxel_materials = [[-2, 2, "Muscle"]] + # patient.voxel_materials = [[-500, -49, "Fat"],[-49, 150, "Muscle"]] + + ##################################################################### + + # dose = sim.add_actor("DoseActor", "dose") + # dose.output_filename = paths.output / "imageVolume_photon_unique.mhd" + # dose.attached_to = patient.name + # dose.size = size + # dose.spacing = spacing + # dose.translation = [0, 0, 0] + # dose.edep_uncertainty.active = False + # dose.hit_type = "random" + + ################################################ Optical photon source + source = sim.add_source("GenericSource", "mysource") + source.particle = "opticalphoton" + source.energy.type = "mono" + source.energy.mono = 1.913 * eV + # source.position.type = "box" # plane source + source.position.type = "point" + # source.position.size = [5.5 * cm, 2.6 * cm, 0] + source.direction.type = "momentum" + source.direction.momentum = [0, 0, 1] + source.position.translation = [ + 0, + 0, + -center[2] * mm, + ] # The source is at the border between the image and the exterior + source.n = 100000 / sim.number_of_threads + + stats = sim.add_actor("SimulationStatisticsActor", "Stats") + stats.track_types_flag = True + + sim.run() + print(stats)