T162 heterogenous sos

simpa/core/simulation_modules/acoustic_forward_module/acoustic_forward_module_k_wave_adapter.py

@@ -17,6 +17,7 @@
 from scipy.spatial.transform import Rotation
 from simpa.core.simulation_modules.acoustic_forward_module import AcousticForwardModelBaseAdapter
 import gc
+from typing import Union
 
 
 class KWaveAdapter(AcousticForwardModelBaseAdapter):
@@ -34,16 +35,16 @@ class KWaveAdapter(AcousticForwardModelBaseAdapter):
         The initial pressure distribution:
             Tags.OPTICAL_MODEL_INITIAL_PRESSURE
         Acoustic tissue properties:
-            Tags.PROPERTY_SPEED_OF_SOUND
-            Tags.PROPERTY_DENSITY
-            Tags.PROPERTY_ALPHA_COEFF
+            Tags.DATA_FIELD_SPEED_OF_SOUND
+            Tags.DATA_FIELD_DENSITY
+            Tags.DATA_FIELD_ALPHA_COEFF
         The digital twin of the imaging device:
             Tags.DIGITAL_DEVICE
         Other parameters:
             Tags.PERFORM_UPSAMPLING
             Tags.SPACING_MM
             Tags.UPSCALE_FACTOR
-            Tags.PROPERTY_ALPHA_POWER
+            Tags.DATA_FIELD_ALPHA_POWER
             Tags.GPU
             Tags.PMLInside
             Tags.PMLAlpha
@@ -127,11 +128,11 @@ def forward_model(self, detection_geometry: DetectionGeometryBase) -> np.ndarray
         return time_series_data
 
     def k_wave_acoustic_forward_model(self, detection_geometry: DetectionGeometryBase,
-                                      speed_of_sound: float, density: float,
+                                      speed_of_sound: Union[float, np.ndarray], density: float,
                                       alpha_coeff: float, initial_pressure: np.ndarray,
                                       optical_path: str = "temporary") -> tuple:
         """
-        Runs the acoustic forward model with the given parameters speed_of_sound (float), density (float),
+        Runs the acoustic forward model with the given parameters speed_of_sound (float or np.ndarray), density (float),
         alpha_coeff (float) for the initial_pressure distribution (numpy array) and a given detection geometry.
         Uses the given optical_path (str) or if none is given a temporary one for saving temporary files.
         Note, that in order to work properly, this function assumes that several settings mentioned above are set.
@@ -165,7 +166,15 @@ def k_wave_acoustic_forward_model(self, detection_geometry: DetectionGeometryBas
         else:
             simulate_2d = False
 
-        data_dict[Tags.DATA_FIELD_SPEED_OF_SOUND] = np.ones_like(initial_pressure) * speed_of_sound
+        self.logger.debug(f"Input SoS has shape {speed_of_sound.shape} and std {np.std(speed_of_sound)}")
+
+        if isinstance(speed_of_sound, float):
+            data_dict[Tags.DATA_FIELD_SPEED_OF_SOUND] = np.ones_like(initial_pressure) * speed_of_sound
+        elif isinstance(speed_of_sound, np.ndarray):
+            data_dict[Tags.DATA_FIELD_SPEED_OF_SOUND] = speed_of_sound
+        else:
+            self.logger.warning("Speed-of-sound has unsuported type (neither float nor np.ndarray)")
+
         data_dict[Tags.DATA_FIELD_DENSITY] = np.ones_like(initial_pressure) * density
         data_dict[Tags.DATA_FIELD_ALPHA_COEFF] = np.ones_like(initial_pressure) * alpha_coeff
         data_dict[Tags.DATA_FIELD_INITIAL_PRESSURE] = initial_pressure
@@ -272,7 +281,7 @@ def k_wave_acoustic_forward_model(self, detection_geometry: DetectionGeometryBas
 
 def perform_k_wave_acoustic_forward_simulation(initial_pressure: np.array,
                                                detection_geometry: DetectionGeometryBase,
-                                               speed_of_sound: float = 1540.0,
+                                               speed_of_sound: Union[float, np.ndarray] = 1540.0,
                                                density: float = 1000.0,
                                                alpha_coeff: float = 0.02,
                                                acoustic_settings: Settings = None,
@@ -288,9 +297,9 @@ def perform_k_wave_acoustic_forward_simulation(initial_pressure: np.array,
                                                spacing_mm: float = 0.5) -> np.array:
     """
     Convenience function for performing a k-Wave acoustic forward simulation using a given detection geometry and
-    initial pressure distribution (numpy array) with the following parameters speed_of_sound (float), density (float),
-    alpha_coeff (float) as well as acoustic_settings (Settings). The acoustic settings may be parsed individually,
-    however, they will be overwritten if they are also set in the acoustic_settings.
+    initial pressure distribution (numpy array) with the following parameters speed_of_sound (float or np.ndarray),
+    density (float), alpha_coeff (float) as well as acoustic_settings (Settings). The acoustic settings may be parsed
+    individually. However, they will be overwritten if they are also set in the acoustic_settings.
 
     :param initial_pressure:
     :param detection_geometry:

simpa/core/simulation_modules/reconstruction_module/__init__.py

@@ -14,6 +14,7 @@
 from simpa.utils import Settings
 from simpa.core.simulation_modules.reconstruction_module.reconstruction_utils import bandpass_filter_with_settings, apply_b_mode
 from simpa.utils.quality_assurance.data_sanity_testing import assert_array_well_defined
+from typing import Union
 
 
 class ReconstructionAdapterBase(SimulationModule):
@@ -89,14 +90,14 @@ def run(self, device):
         self.logger.info("Performing reconstruction...[Done]")
 
 
-def create_reconstruction_settings(speed_of_sound_in_m_per_s: int = 1540, time_spacing_in_s: float = 2.5e-8,
+def create_reconstruction_settings(speed_of_sound_in_m_per_s: Union[int, np.ndarray] = 1540, time_spacing_in_s: float = 2.5e-8,
                                    sensor_spacing_in_mm: float = 0.1,
                                    recon_mode: str = Tags.RECONSTRUCTION_MODE_PRESSURE,
                                    apodization: str = Tags.RECONSTRUCTION_APODIZATION_BOX) -> Settings:
     """
     Function that creates SIMPA settings for reconstruction convenience function.
 
-    :param speed_of_sound_in_m_per_s: (int) speed of sound in medium in meters per second (default: 1540 m/s)
+    :param speed_of_sound_in_m_per_s: (int or np.ndarry) speed of sound in medium in meters per second (default: 1540 m/s)
     :param time_spacing_in_s: (float) time between sampling points in seconds (default: 2.5e-8 s which is equal to 40 MHz)
     :param sensor_spacing_in_mm: (float) space between sensor elements in millimeters (default: 0.1 mm)
     :param recon_mode: SIMPA Tag defining the reconstruction mode - pressure default OR differential

simpa/core/simulation_modules/reconstruction_module/reconstruction_module_delay_and_sum_adapter.py

@@ -10,6 +10,7 @@
     compute_image_dimensions, preparing_reconstruction_and_obtaining_reconstruction_settings
 from simpa.core.device_digital_twins import DetectionGeometryBase
 from simpa.core.simulation_modules.reconstruction_module import create_reconstruction_settings
+from typing import Union
 
 
 class DelayAndSumAdapter(ReconstructionAdapterBase):
@@ -42,8 +43,8 @@ def reconstruction_algorithm(self, time_series_sensor_data, detection_geometry:
 
         values, _ = compute_delay_and_sum_values(time_series_sensor_data, sensor_positions, xdim,
                                                  ydim, zdim, xdim_start, xdim_end, ydim_start, ydim_end, zdim_start, zdim_end, spacing_in_mm, speed_of_sound_in_m_per_s,
-                                                 time_spacing_in_ms, self.logger, torch_device,
-                                                 self.component_settings)
+                                                 time_spacing_in_ms, self.logger, torch_device, self.component_settings,
+                                                 self.global_settings, detection_geometry.device_position_mm)
 
         _sum = torch.sum(values, dim=3)
         counter = torch.count_nonzero(values, dim=3)
@@ -56,7 +57,7 @@ def reconstruction_algorithm(self, time_series_sensor_data, detection_geometry:
 
 def reconstruct_delay_and_sum_pytorch(time_series_sensor_data: np.ndarray,
                                       detection_geometry: DetectionGeometryBase,
-                                      speed_of_sound_in_m_per_s: int = 1540,
+                                      speed_of_sound_in_m_per_s: Union[int, np.ndarray] = 1540,
                                       time_spacing_in_s: float = 2.5e-8,
                                       sensor_spacing_in_mm: float = 0.1,
                                       recon_mode: str = Tags.RECONSTRUCTION_MODE_PRESSURE,
@@ -66,7 +67,7 @@ def reconstruct_delay_and_sum_pytorch(time_series_sensor_data: np.ndarray,
 
     :param time_series_sensor_data: (2D numpy array) sensor data of shape (sensor elements, time steps)
     :param detection_geometry: The DetectionGeometryBase that should be used to reconstruct the given time series data
-    :param speed_of_sound_in_m_per_s: (int) speed of sound in medium in meters per second (default: 1540 m/s)
+    :param speed_of_sound_in_m_per_s: (int or np.ndarray) speed of sound in medium in meters per second (default: 1540 m/s)
     :param time_spacing_in_s: (float) time between sampling points in seconds (default: 2.5e-8 s which is equal to 40 MHz)
     :param sensor_spacing_in_mm: (float) space between sensor elements in millimeters (default: 0.1 mm)
     :param recon_mode: SIMPA Tag defining the reconstruction mode - pressure default OR differential