Scattering Module

Direction

Bases: Enum

Source code in pytissueoptics/rayscattering/display/utils/direction.py

class Direction(Enum):
    X_POS = 0
    Y_POS = 1
    Z_POS = 2
    X_NEG = 3
    Y_NEG = 4
    Z_NEG = 5

    def isSameAxisAs(self, other) -> bool:
        return self.value % 3 == other.value % 3

    @property
    def axis(self) -> int:
        """Returns an integer between 0 and 2 representing the x, y, or z axis, ignoring direction sign."""
        return self.value % 3

    @property
    def isNegative(self) -> bool:
        return self.value >= 3

    @property
    def isPositive(self) -> bool:
        return not self.isNegative

    @property
    def sign(self) -> int:
        return 1 if self.isPositive else -1

axis property

Returns an integer between 0 and 2 representing the x, y, or z axis, ignoring direction sign.

EnergyLogger

Bases: Logger

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

class EnergyLogger(Logger):
    _data: dict[InteractionKey, InteractionData]

    def __init__(
        self,
        scene: ScatteringScene,
        filepath: str = None,
        keep3D: bool = True,
        views: Union[ViewGroup, List[View2D]] = ViewGroup.ALL,
        defaultViewEnergyType: EnergyType = EnergyType.DEPOSITION,
        defaultBinSize: Union[float, tuple] = 0.01,
        infiniteLimits=((-5, 5), (-5, 5), (-5, 5)),
    ):
        """
        Log the energy deposited by scattering photons as well as the energy that crossed surfaces. Every interaction
        is linked to a specific solid and surface of the scene when applicable. This `EnergyLogger` has to be given to
        a `Source` instance when propagating.

        :param scene: The scene is required to auto-initialize the 2D views correctly. Even if `keep3D` is True, the
                scene may be used later to correctly create and display the 2D views.
        :param filepath: (Optional) Load the logger from the given file. When saving, the file will be overwritten.
        :param keep3D: (Default to True) If False, logged datapoints are automatically binned to predefined 2D `views`
                and the 3D data is discarded. This allows for a lightweight logger alternative. If True, the 3D data is
                kept and the 2D views are only computed later if displayed.
        :param views: The 2D views to track, particularly used when `keep3D` is set to False. Can be a ViewGroup flag
                or a list of View2D. Default to ViewGroup.ALL which includes the 3 default XYZ projections for the
                whole scene as well as for each solid, and a projection of each surface in the direction of the surface
                normal for both the energy entering and leaving the surface. When `keep3D` is False, ViewGroup.All is
                required to allow calculation of the detailed report for absorbance and transmittance. When `keep3D` is
                True, the 2D views data is only computed when displayed, so it is equivalent to setting `views` to None.
        :param defaultBinSize: The default bin size to use when binning the 3D data to 2D views. In the same physical
                units as the scene. Custom bin sizes can be specified in each View2D.
        :param infiniteLimits: The default limits to use for the 2D views when the scene is infinite (has no solids).
        """
        self._scene = scene
        self._keep3D = keep3D
        self._defaultBinSize = defaultBinSize
        self._infiniteLimits = infiniteLimits
        self._viewFactory = ViewFactory(scene, defaultBinSize, infiniteLimits, energyType=defaultViewEnergyType)

        self._sceneHash = hash(scene)
        self._defaultViews = views
        self._views = self._viewFactory.build(views)
        self._outdatedViews = set()
        self._nDataPointsRemoved = 0

        super().__init__(fromFilepath=filepath)

    def addView(self, view: View2D) -> bool:
        self._viewFactory.build([view])

        if self._viewExists(view):
            return True

        if self.isEmpty:
            self._views.append(view)
            return True

        if self.has3D:
            self._compileViews([view], detectedBy=view.detectedBy)
            self._views.append(view)
            return True

        for i, existingView in enumerate(self._views):
            if view.isContainedBy(existingView):
                view.initDataFrom(existingView)
                if existingView in self._outdatedViews:
                    self._outdatedViews.add(view)
                self._views.append(view)
                return True

        utils.warn(
            f"ERROR: Cannot create view {view.name}. The 3D data was discarded and the required data was not "
            f"found in existing views."
        )
        return False

    def updateView(self, view: View2D):
        if view in self._outdatedViews:
            self._compileViews([view])

    def showView(self, view: View2D = None, viewIndex: int = None, logScale: bool = True, colormap: str = "viridis"):
        assert viewIndex is not None or view is not None, "Either `viewIndex` or `view` must be specified."

        if viewIndex is None:
            created = self.addView(view)
            if not created:
                utils.warn(f"ERROR: Cannot display view {view.name}. Failed to create the view.")
                return
            viewIndex = self._getViewIndex(view)

        view = self.getView(viewIndex)
        self.updateView(view)

        view.show(logScale=logScale, colormap=colormap)

    def listViews(self):
        print("Available views:")
        for i, view in enumerate(self._views):
            print(f"\t{i}: {view.description}")

    def save(self, filepath: str = None):
        if filepath is None and self._filepath is None:
            filepath = self.DEFAULT_LOGGER_PATH
            utils.warn(f"No filepath specified. Saving to {filepath}.")
        elif filepath is None:
            filepath = self._filepath

        with open(filepath, "wb") as file:
            pickle.dump(
                (
                    self._data,
                    self.info,
                    self._labels,
                    self._views,
                    self._defaultViews,
                    self._outdatedViews,
                    self._nDataPointsRemoved,
                    self._sceneHash,
                    self.has3D,
                ),
                file,
            )

    def load(self, filepath: str):
        self._filepath = filepath

        if not os.path.exists(filepath):
            utils.warn(
                "No logger file found at '{}'. No data loaded, but it will create a new file "
                "at this location if the logger is saved later on.".format(filepath)
            )
            return

        with open(filepath, "rb") as file:
            (
                self._data,
                self.info,
                self._labels,
                self._views,
                oldDefaultViews,
                self._outdatedViews,
                self._nDataPointsRemoved,
                oldSceneHash,
                oldHas3D,
            ) = pickle.load(file)

        if oldSceneHash != self._sceneHash:
            utils.warn(
                "WARNING: The scene used to create the logger at '{}' is different from the current "
                "scene. This may corrupt statistics and visualization. Proceed at your own risk.".format(filepath)
            )
        if oldHas3D and not self._keep3D:
            utils.warn(
                "WARNING: The logger at '{}' use to store 3D data, but it was reloaded with keep3D=False. "
                "The 3D data will be compiled to 2D views and discarded.".format(filepath)
            )
        if not oldHas3D and self._keep3D:
            utils.warn(
                "WARNING: The logger at '{}' use to discard 3D data, but it was reloaded with keep3D=True. "
                "This may corrupt the statistics and the 3D visualization. Proceed at your own risk.".format(filepath)
            )
        if self._defaultViews != oldDefaultViews:
            utils.warn(
                "WARNING: Cannot provide new default views to a loaded logger from '{}'."
                "Using only the views from the file.".format(filepath)
            )

    @property
    def views(self) -> List[View2D]:
        return self._views

    def getView(self, index: int) -> View2D:
        if index < 0 or index >= len(self._views):
            raise IndexError(
                f"View index {index} is out of range [0, {len(self._views)}]. Use `.listViews()` to see "
                f"available views."
            )
        return self._views[index]

    def _getViewIndex(self, view: View2D) -> int:
        for i, v in enumerate(self._views):
            if v.isEqualTo(view):
                return i
        raise ValueError(f"View {view.name} not found in the list of views.")

    def _viewExists(self, view: View2D) -> bool:
        return any([view.isEqualTo(v) for v in self._views])

    @property
    def has3D(self) -> bool:
        return self._keep3D

    @property
    def defaultBinSize(self) -> float:
        return self._defaultBinSize

    @property
    def infiniteLimits(self) -> tuple:
        return self._infiniteLimits

    def getSolidLimits(self, solidLabel: str) -> List[List[float]]:
        return self._scene.getSolid(solidLabel).getBoundingBox().xyzLimits

    def logDataPointArray(self, array: np.ndarray, key: InteractionKey):
        """
        Used internally by `Source` when propagating photons. Overwrites the `Logger` method to automatically bin the
        data to 2D views if 3D data is being discarded.
        """
        super().logDataPointArray(array, key)
        self._outdatedViews = set(self._views)

        if not self._keep3D:
            self._compileViews(self._views)
            self._delete3DData()

    def logDataPoint(self, value: float, position: Vector, key: InteractionKey, ID: Optional[int] = None):
        dataPoint = [value, *position.array]
        if ID is not None:
            dataPoint.append(ID)
        self.logDataPointArray(np.array([dataPoint]), key)

    def _compileViews(self, views: List[View2D], detectedBy: Union[str, List[str]] = None):
        if detectedBy is None:
            dataPerInteraction = self._data
            if any(view.detectedBy for view in views):
                utils.warn(
                    "Ignoring the detectedBy property of a view. Can only use detectedBy when adding a view after the "
                    "simulation with keep3D=True."
                )
        else:
            dataPerInteraction = self.getFiltered(detectedBy)._data

        for key, data in dataPerInteraction.items():
            datapointsContainer: Optional[ListArrayContainer] = data.dataPoints
            if datapointsContainer is None or len(datapointsContainer) == 0:
                continue
            for view in views:
                if view.solidLabel and not utils.labelsEqual(view.solidLabel, key.solidLabel):
                    continue
                if view.surfaceLabel and not utils.labelsEqual(view.surfaceLabel, key.surfaceLabel):
                    continue
                if view.surfaceLabel is None and key.surfaceLabel is not None:
                    continue

                data = datapointsContainer.getData()
                if view.energyType == EnergyType.FLUENCE_RATE:
                    data = self._fluenceTransform(key, data)

                view.extractData(data)
        for view in views:
            self._outdatedViews.discard(view)

    def _delete3DData(self):
        self._nDataPointsRemoved += super().nDataPoints
        self._data.clear()

    @property
    def nDataPoints(self) -> int:
        """
        Overwrites the `Logger` method to return the total number of data points logged, including the ones discarded.
        """
        if self._keep3D:
            return super().nDataPoints
        else:
            return self._nDataPointsRemoved

    @property
    def isEmpty(self) -> bool:
        return self.nDataPoints == 0

    def logPoint(self, point: Vector, key: InteractionKey = None):
        raise NotImplementedError("Can only log data points to an EnergyLogger.")

    def logPointArray(self, array: np.ndarray, key: InteractionKey = None):
        raise NotImplementedError("Can only log data points to an EnergyLogger.")

    def logSegment(self, start: Vector, end: Vector, key: InteractionKey = None):
        raise NotImplementedError("Can only log data points to an EnergyLogger.")

    def logSegmentArray(self, array: np.ndarray, key: InteractionKey = None):
        raise NotImplementedError("Can only log data points to an EnergyLogger.")

    def getDataPoints(self, key: InteractionKey, energyType=EnergyType.DEPOSITION) -> np.ndarray:
        """All 3D data points recorded for this InteractionKey (not binned). Array of shape (n, 4) where
        the second axis is (value, x, y, z). The value can be the energy deposited, the fluence rate, or the
        energy that crossed the surface.

        :param key: Filtering the data by solidLabel and surfaceLabel.
        :param energyType: The type of volumetric energy to return when no surfaceLabel is given.

        :return: The data points (value, x, y, z) for the given solidLabel and surfaceLabel. If a surfaceLabel is given,
        the value corresponds to the energy that crossed the surface (positive when in the direction of the normal). If
        only a solidLabel is given, the value corresponds to the volumetric EnergyType at that point.
        """
        if energyType == EnergyType.FLUENCE_RATE:
            return self._getData(DataType.DATA_POINT, key, transform=self._fluenceTransform)

        return self._getData(DataType.DATA_POINT, key)

    def filter(self, detectedBy: Union[str, List[str]]) -> None:
        """Keeps only the data points from photons detected by one of the specified detector(s)."""
        if not self._keep3D:
            utils.warn("Cannot filter a logger that has discarded the 3D data.")
            return

        filteredPhotonIDs = self._getDetectedPhotonIDs(detectedBy)
        self._data = self._getDataForPhotons(filteredPhotonIDs)
        self._outdatedViews = set(self._views)

    def getFiltered(self, detectedBy: Union[str, List[str]]) -> "EnergyLogger":
        """
        Returns a new logger with only data from photons detected by one of the specified detector(s).
        """
        filteredLogger = EnergyLogger(self._scene, views=[])
        filteredPhotonIDs = self._getDetectedPhotonIDs(detectedBy)
        filteredLogger._data = self._getDataForPhotons(filteredPhotonIDs)
        return filteredLogger

    def _getDetectedPhotonIDs(self, detectedBy: Union[str, List[str]]) -> np.ndarray:
        """Helper to get photon IDs detected by one of the specified detector(s)."""
        detector_labels = [detectedBy] if isinstance(detectedBy, str) else detectedBy
        detector_keys = [InteractionKey(label) for label in detector_labels]
        photonIDs = self._getPhotonIDs(detector_keys)
        if len(photonIDs) == 0:
            utils.warn(f"No photons detected by: {detectedBy}")
        return photonIDs

    def _getPhotonIDs(self, key: Union[InteractionKey, List[InteractionKey], None] = None) -> np.ndarray:
        """Get all unique photon IDs that interacted with one of the given interaction key(s)."""
        if isinstance(key, list):
            all_photon_ids = []
            for k in key:
                data = self.getRawDataPoints(k)
                if data is not None and data.shape[1] >= 5:
                    all_photon_ids.append(data[:, 4].astype(np.uint32))

            if len(all_photon_ids) == 0:
                return np.array([], dtype=np.uint32)
            combined = np.concatenate(all_photon_ids)
            return np.unique(combined)

        data = self.getRawDataPoints(key)
        if data is None or data.shape[1] < 5:
            return np.array([], dtype=np.uint32)
        return np.unique(data[:, 4].astype(np.uint32))

    def _getDataForPhotons(self, photonIDs: np.ndarray) -> Dict[InteractionKey, InteractionData]:
        keyToData: Dict[InteractionKey, InteractionData] = {}
        photonIDs = np.asarray(photonIDs, dtype=np.uint32)
        for key, interactionData in self._data.items():
            points: Optional[ListArrayContainer] = interactionData.dataPoints
            if points is None:
                continue
            data = points.getData()
            if data.shape[1] < 5:
                continue
            mask = np.isin(data[:, 4].astype(np.uint32), photonIDs)
            filteredData = data[mask]
            if filteredData.size > 0:
                container = ListArrayContainer()
                container.append(filteredData)
                keyToData[key] = InteractionData(dataPoints=container)
        return keyToData

    def _fluenceTransform(self, key: InteractionKey, data: Optional[np.ndarray]) -> Optional[np.ndarray]:
        # Converts volumetric data to fluence rate when needed.
        if not key.volumetric or data is None:
            return data

        data[:, 0] = data[:, 0] / self._scene.getMaterial(key.solidLabel).mu_a
        return data

    def export(self, exportName: str):
        """
        Export the raw 3D data points to a CSV file, along with the scene information to a JSON file.

        The data file <exportName>.csv will be comma-delimited and will contain the following columns:
        - energy, x, y, z, photon_index, solid_index, surface_index

        Two types of interactions are logged: scattering and surface crossings. In the first case, the energy will be
        the delta energy deposited at the point and the surface index will be -1. In the second case, the energy
        will be the total photon energy when crossing the surface, either as positive if leaving the surface
        (along the normal) or as negative if entering the surface.

        The scene information will be saved in a JSON file named <exportName>.json, which includes details for each solid
        index and surface index, such as their labels, materials, and geometry. The world information is also exported
        as solid index -1.
        """
        if not self.has3D:
            utils.warn("Cannot export data when keep3D is False. No 3D data available.")
            return

        solidLabels = []
        for solid in self._scene.solids:
            if solid.isStack():
                solidLabels.extend(solid.getLayerLabels())
            else:
                solidLabels.append(solid.getLabel())
        solidLabels.sort()

        print("Exporting raw data to file...")
        filepath = f"{exportName}.csv"
        with open(filepath, "w") as file:
            file.write("energy,x,y,z,photon_index,solid_index,surface_index\n")
            self._writeKeyData(file, InteractionKey(WORLD_SOLID_LABEL), -1, -1)
            for i, solidLabel in enumerate(solidLabels):
                self._writeKeyData(file, InteractionKey(solidLabel), i, -1)
                for j, surfaceLabel in enumerate(self._scene.getSurfaceLabels(solidLabel)):
                    self._writeKeyData(file, InteractionKey(solidLabel, surfaceLabel), i, j)
        print(f"Exported data points to {filepath}")

        self._exportSceneInfo(f"{exportName}.json", solidLabels)

    def _writeKeyData(self, file: TextIO, key: InteractionKey, solidIndex: int, surfaceIndex: int):
        if key not in self._data or self._data[key].dataPoints is None:
            return

        dataArray = self._data[key].dataPoints.getData()
        n_rows = dataArray.shape[0]

        output = np.empty((n_rows, 7), dtype=np.float64)
        output[:, :4] = dataArray[:, :4]
        output[:, 4] = dataArray[:, 4].astype(np.uint32)
        output[:, 5] = solidIndex
        output[:, 6] = surfaceIndex

        np.savetxt(
            file,
            output,
            delimiter=",",
            fmt=["%.8e", "%.8e", "%.8e", "%.8e", "%d", "%d", "%d"],
        )

    def _exportSceneInfo(self, filepath: str, solidLabels: List[str]):
        sceneInfo = {}
        material = self._scene.getWorldEnvironment().material
        sceneInfo["-1"] = {"label": "world", "material": material.__dict__ if material else None}
        for i, solidLabel in enumerate(solidLabels):
            material = self._scene.getMaterial(solidLabel)
            solid = self._scene.getSolid(solidLabel)
            surfaces = {}
            for j, surfaceLabel in enumerate(solid.surfaceLabels):
                normals = [s.normal for s in solid.getPolygons(surfaceLabel)[:2]]
                if len(normals) == 1 or normals[0] == normals[1]:
                    normal = normals[0].array
                else:
                    normal = None
                surfaces[j] = {"label": surfaceLabel, "normal": normal}

            sceneInfo[str(i)] = {
                "label": solidLabel,
                "type": solid.__class__.__name__,
                "material": material.__dict__ if material else None,
                "geometry": solid.geometryExport(),
                "surfaces": surfaces,
            }

        with open(filepath, "w") as file:
            json.dump(sceneInfo, file, indent=4)
        print(f"Exported scene information to {filepath}")

nDataPoints property

Overwrites the Logger method to return the total number of data points logged, including the ones discarded.

__init__(scene, filepath=None, keep3D=True, views=ViewGroup.ALL, defaultViewEnergyType=EnergyType.DEPOSITION, defaultBinSize=0.01, infiniteLimits=((-5, 5), (-5, 5), (-5, 5)))

Log the energy deposited by scattering photons as well as the energy that crossed surfaces. Every interaction is linked to a specific solid and surface of the scene when applicable. This EnergyLogger has to be given to a Source instance when propagating.

:param scene: The scene is required to auto-initialize the 2D views correctly. Even if keep3D is True, the scene may be used later to correctly create and display the 2D views. :param filepath: (Optional) Load the logger from the given file. When saving, the file will be overwritten. :param keep3D: (Default to True) If False, logged datapoints are automatically binned to predefined 2D views and the 3D data is discarded. This allows for a lightweight logger alternative. If True, the 3D data is kept and the 2D views are only computed later if displayed. :param views: The 2D views to track, particularly used when keep3D is set to False. Can be a ViewGroup flag or a list of View2D. Default to ViewGroup.ALL which includes the 3 default XYZ projections for the whole scene as well as for each solid, and a projection of each surface in the direction of the surface normal for both the energy entering and leaving the surface. When keep3D is False, ViewGroup.All is required to allow calculation of the detailed report for absorbance and transmittance. When keep3D is True, the 2D views data is only computed when displayed, so it is equivalent to setting views to None. :param defaultBinSize: The default bin size to use when binning the 3D data to 2D views. In the same physical units as the scene. Custom bin sizes can be specified in each View2D. :param infiniteLimits: The default limits to use for the 2D views when the scene is infinite (has no solids).

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

def __init__(
    self,
    scene: ScatteringScene,
    filepath: str = None,
    keep3D: bool = True,
    views: Union[ViewGroup, List[View2D]] = ViewGroup.ALL,
    defaultViewEnergyType: EnergyType = EnergyType.DEPOSITION,
    defaultBinSize: Union[float, tuple] = 0.01,
    infiniteLimits=((-5, 5), (-5, 5), (-5, 5)),
):
    """
    Log the energy deposited by scattering photons as well as the energy that crossed surfaces. Every interaction
    is linked to a specific solid and surface of the scene when applicable. This `EnergyLogger` has to be given to
    a `Source` instance when propagating.

    :param scene: The scene is required to auto-initialize the 2D views correctly. Even if `keep3D` is True, the
            scene may be used later to correctly create and display the 2D views.
    :param filepath: (Optional) Load the logger from the given file. When saving, the file will be overwritten.
    :param keep3D: (Default to True) If False, logged datapoints are automatically binned to predefined 2D `views`
            and the 3D data is discarded. This allows for a lightweight logger alternative. If True, the 3D data is
            kept and the 2D views are only computed later if displayed.
    :param views: The 2D views to track, particularly used when `keep3D` is set to False. Can be a ViewGroup flag
            or a list of View2D. Default to ViewGroup.ALL which includes the 3 default XYZ projections for the
            whole scene as well as for each solid, and a projection of each surface in the direction of the surface
            normal for both the energy entering and leaving the surface. When `keep3D` is False, ViewGroup.All is
            required to allow calculation of the detailed report for absorbance and transmittance. When `keep3D` is
            True, the 2D views data is only computed when displayed, so it is equivalent to setting `views` to None.
    :param defaultBinSize: The default bin size to use when binning the 3D data to 2D views. In the same physical
            units as the scene. Custom bin sizes can be specified in each View2D.
    :param infiniteLimits: The default limits to use for the 2D views when the scene is infinite (has no solids).
    """
    self._scene = scene
    self._keep3D = keep3D
    self._defaultBinSize = defaultBinSize
    self._infiniteLimits = infiniteLimits
    self._viewFactory = ViewFactory(scene, defaultBinSize, infiniteLimits, energyType=defaultViewEnergyType)

    self._sceneHash = hash(scene)
    self._defaultViews = views
    self._views = self._viewFactory.build(views)
    self._outdatedViews = set()
    self._nDataPointsRemoved = 0

    super().__init__(fromFilepath=filepath)

export(exportName)

Export the raw 3D data points to a CSV file, along with the scene information to a JSON file.

The data file .csv will be comma-delimited and will contain the following columns: - energy, x, y, z, photon_index, solid_index, surface_index

Two types of interactions are logged: scattering and surface crossings. In the first case, the energy will be the delta energy deposited at the point and the surface index will be -1. In the second case, the energy will be the total photon energy when crossing the surface, either as positive if leaving the surface (along the normal) or as negative if entering the surface.

The scene information will be saved in a JSON file named .json, which includes details for each solid index and surface index, such as their labels, materials, and geometry. The world information is also exported as solid index -1.

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

def export(self, exportName: str):
    """
    Export the raw 3D data points to a CSV file, along with the scene information to a JSON file.

    The data file <exportName>.csv will be comma-delimited and will contain the following columns:
    - energy, x, y, z, photon_index, solid_index, surface_index

    Two types of interactions are logged: scattering and surface crossings. In the first case, the energy will be
    the delta energy deposited at the point and the surface index will be -1. In the second case, the energy
    will be the total photon energy when crossing the surface, either as positive if leaving the surface
    (along the normal) or as negative if entering the surface.

    The scene information will be saved in a JSON file named <exportName>.json, which includes details for each solid
    index and surface index, such as their labels, materials, and geometry. The world information is also exported
    as solid index -1.
    """
    if not self.has3D:
        utils.warn("Cannot export data when keep3D is False. No 3D data available.")
        return

    solidLabels = []
    for solid in self._scene.solids:
        if solid.isStack():
            solidLabels.extend(solid.getLayerLabels())
        else:
            solidLabels.append(solid.getLabel())
    solidLabels.sort()

    print("Exporting raw data to file...")
    filepath = f"{exportName}.csv"
    with open(filepath, "w") as file:
        file.write("energy,x,y,z,photon_index,solid_index,surface_index\n")
        self._writeKeyData(file, InteractionKey(WORLD_SOLID_LABEL), -1, -1)
        for i, solidLabel in enumerate(solidLabels):
            self._writeKeyData(file, InteractionKey(solidLabel), i, -1)
            for j, surfaceLabel in enumerate(self._scene.getSurfaceLabels(solidLabel)):
                self._writeKeyData(file, InteractionKey(solidLabel, surfaceLabel), i, j)
    print(f"Exported data points to {filepath}")

    self._exportSceneInfo(f"{exportName}.json", solidLabels)

filter(detectedBy)

Keeps only the data points from photons detected by one of the specified detector(s).

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

def filter(self, detectedBy: Union[str, List[str]]) -> None:
    """Keeps only the data points from photons detected by one of the specified detector(s)."""
    if not self._keep3D:
        utils.warn("Cannot filter a logger that has discarded the 3D data.")
        return

    filteredPhotonIDs = self._getDetectedPhotonIDs(detectedBy)
    self._data = self._getDataForPhotons(filteredPhotonIDs)
    self._outdatedViews = set(self._views)

getDataPoints(key, energyType=EnergyType.DEPOSITION)

All 3D data points recorded for this InteractionKey (not binned). Array of shape (n, 4) where the second axis is (value, x, y, z). The value can be the energy deposited, the fluence rate, or the energy that crossed the surface.

:param key: Filtering the data by solidLabel and surfaceLabel. :param energyType: The type of volumetric energy to return when no surfaceLabel is given.

:return: The data points (value, x, y, z) for the given solidLabel and surfaceLabel. If a surfaceLabel is given, the value corresponds to the energy that crossed the surface (positive when in the direction of the normal). If only a solidLabel is given, the value corresponds to the volumetric EnergyType at that point.

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

def getDataPoints(self, key: InteractionKey, energyType=EnergyType.DEPOSITION) -> np.ndarray:
    """All 3D data points recorded for this InteractionKey (not binned). Array of shape (n, 4) where
    the second axis is (value, x, y, z). The value can be the energy deposited, the fluence rate, or the
    energy that crossed the surface.

    :param key: Filtering the data by solidLabel and surfaceLabel.
    :param energyType: The type of volumetric energy to return when no surfaceLabel is given.

    :return: The data points (value, x, y, z) for the given solidLabel and surfaceLabel. If a surfaceLabel is given,
    the value corresponds to the energy that crossed the surface (positive when in the direction of the normal). If
    only a solidLabel is given, the value corresponds to the volumetric EnergyType at that point.
    """
    if energyType == EnergyType.FLUENCE_RATE:
        return self._getData(DataType.DATA_POINT, key, transform=self._fluenceTransform)

    return self._getData(DataType.DATA_POINT, key)

getFiltered(detectedBy)

Returns a new logger with only data from photons detected by one of the specified detector(s).

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

def getFiltered(self, detectedBy: Union[str, List[str]]) -> "EnergyLogger":
    """
    Returns a new logger with only data from photons detected by one of the specified detector(s).
    """
    filteredLogger = EnergyLogger(self._scene, views=[])
    filteredPhotonIDs = self._getDetectedPhotonIDs(detectedBy)
    filteredLogger._data = self._getDataForPhotons(filteredPhotonIDs)
    return filteredLogger

logDataPointArray(array, key)

Used internally by Source when propagating photons. Overwrites the Logger method to automatically bin the data to 2D views if 3D data is being discarded.

Source code in pytissueoptics/rayscattering/energyLogging/energyLogger.py

def logDataPointArray(self, array: np.ndarray, key: InteractionKey):
    """
    Used internally by `Source` when propagating photons. Overwrites the `Logger` method to automatically bin the
    data to 2D views if 3D data is being discarded.
    """
    super().logDataPointArray(array, key)
    self._outdatedViews = set(self._views)

    if not self._keep3D:
        self._compileViews(self._views)
        self._delete3DData()

EnergyType

Bases: Enum

Type of volumetric energy: either as the deposited energy in the solid (absorption) or as the fluence rate.

Source code in pytissueoptics/rayscattering/energyLogging/energyType.py

class EnergyType(Enum):
    """
    Type of volumetric energy: either as the deposited energy in the solid (absorption) or as the fluence rate.
    """

    DEPOSITION = auto()
    FLUENCE_RATE = auto()

PointCloudStyle

3D display configuration for solid and surface point clouds.

Visibility attributes

solidLabel (Optional[str]): Only show the point cloud specific to a single solid. surfaceLabel (Optional[str]): Only show the point cloud specific to a single surface of the solid. showSolidPoints (bool): Show the point clouds of the solids. showSurfacePointsLeaving (bool): Show energy that left the surface (direction with surface normal). showSurfacePointsEntering (bool): Show energy that entered the surface (direction opposite to surface normal). energyType (EnergyType): Type of energy to show for volumetric datapoints (deposition or fluence).

Other attributes

showPointsAsSpheres (bool): Show the points as spheres or as dots. Dots require less memory. pointSize (float): Reference diameter of the points in the point cloud when drawn as spheres. scaleWithValue (bool): Scale the points with their value. A value of 1 yields the pointSize. colormap (str): The name of the colormap to use for the point cloud. reverseColormap (bool): Reverse the colormap. surfacePointSize (float): Same as pointSize, but for the surface points. surfaceScaleWithValue (bool): Same as scaleWithValue but for the surface points. surfaceColormap (str): Same as colormap but for the surface points. surfaceReverseColormap (bool): Same as reverseColormap but for the surface points.

Source code in pytissueoptics/rayscattering/display/viewer.py

class PointCloudStyle:
    """
    3D display configuration for solid and surface point clouds.

    Visibility attributes:
        solidLabel (Optional[str]): Only show the point cloud specific to a single solid.
        surfaceLabel (Optional[str]): Only show the point cloud specific to a single surface of the solid.
        showSolidPoints (bool): Show the point clouds of the solids.
        showSurfacePointsLeaving (bool): Show energy that left the surface (direction with surface normal).
        showSurfacePointsEntering (bool): Show energy that entered the surface (direction opposite to surface normal).
        energyType (EnergyType): Type of energy to show for volumetric datapoints (deposition or fluence).

    Other attributes:
        showPointsAsSpheres (bool): Show the points as spheres or as dots. Dots require less memory.
        pointSize (float): Reference diameter of the points in the point cloud when drawn as spheres.
        scaleWithValue (bool): Scale the points with their value. A value of 1 yields the `pointSize`.
        colormap (str): The name of the colormap to use for the point cloud.
        reverseColormap (bool): Reverse the colormap.
        surfacePointSize (float): Same as `pointSize`, but for the surface points.
        surfaceScaleWithValue (bool): Same as `scaleWithValue` but for the surface points.
        surfaceColormap (str): Same as `colormap` but for the surface points.
        surfaceReverseColormap (bool): Same as `reverseColormap` but for the surface points.
    """

    def __init__(
        self,
        solidLabel: str = None,
        surfaceLabel: str = None,
        showSolidPoints: bool = True,
        showSurfacePointsLeaving: bool = True,
        showSurfacePointsEntering: bool = False,
        energyType=EnergyType.DEPOSITION,
        detectedBy: Union[str, List[str]] = None,
        showPointsAsSpheres: bool = False,
        pointSize: float = 0.15,
        scaleWithValue: bool = True,
        colormap: str = "rainbow",
        reverseColormap: bool = False,
        surfacePointSize: float = 0.01,
        surfaceScaleWithValue: bool = False,
        surfaceColormap: str = None,
        surfaceReverseColormap: bool = None,
    ):
        self.solidLabel = solidLabel
        self.surfaceLabel = surfaceLabel
        self.showSolidPoints = showSolidPoints
        self.showSurfacePointsLeaving = showSurfacePointsLeaving
        self.showSurfacePointsEntering = showSurfacePointsEntering
        self.energyType = energyType
        self.detectedBy = detectedBy
        self.showPointsAsSpheres = showPointsAsSpheres

        self.pointSize = pointSize
        self.scaleWithValue = scaleWithValue
        self.colormap = colormap
        self.reverseColormap = reverseColormap

        self.surfacePointSize = surfacePointSize
        self.surfaceScaleWithValue = surfaceScaleWithValue
        self.surfaceColormap = colormap if surfaceColormap is None else surfaceColormap
        self.surfaceReverseColormap = reverseColormap if surfaceReverseColormap is None else surfaceReverseColormap

ScatteringScene

Bases: Scene

Source code in pytissueoptics/rayscattering/scatteringScene.py

class ScatteringScene(Scene):
    def __init__(self, solids: List[Solid], worldMaterial=ScatteringMaterial(), ignoreIntersections: bool = False):
        super().__init__(solids, worldMaterial=worldMaterial, ignoreIntersections=ignoreIntersections)

    def add(self, solid: Solid, position: Vector = None):
        if solid.isFlat and not solid.isDetector:
            raise Exception(f"Solid '{solid.getLabel()}' is flat. Flat solids must be used with asDetector().")

        polygonSample = solid.getPolygons()[0]
        if not isinstance(polygonSample.insideEnvironment.material, ScatteringMaterial) and not solid.isDetector:
            raise Exception(
                f"Solid '{solid.getLabel()}' has no ScatteringMaterial defined. "
                f"This is required for any RayScatteringScene. "
            )
        super().add(solid, position)

    def show(self, source: Displayable = None, opacity=0.8, colormap="cool", **kwargs):
        viewer = get3DViewer()
        self.addToViewer(viewer, opacity=opacity, colormap=colormap, **kwargs)
        if source:
            source.addToViewer(viewer)
        viewer.show()

    def getEstimatedIPP(self, weightThreshold: float) -> float:
        """
        Get the estimated number of interactions per photon. This gross estimation is done by assuming an infinite
        medium of mean scene albedo. Used as a starting point for the OpenCL kernel optimization.
        """
        materials = self.getMaterials()
        averageAlbedo = sum([mat.getAlbedo() for mat in materials]) / len(materials)
        estimatedIPP = -np.log(weightThreshold) / averageAlbedo
        return estimatedIPP

getEstimatedIPP(weightThreshold)

Get the estimated number of interactions per photon. This gross estimation is done by assuming an infinite medium of mean scene albedo. Used as a starting point for the OpenCL kernel optimization.

Source code in pytissueoptics/rayscattering/scatteringScene.py

def getEstimatedIPP(self, weightThreshold: float) -> float:
    """
    Get the estimated number of interactions per photon. This gross estimation is done by assuming an infinite
    medium of mean scene albedo. Used as a starting point for the OpenCL kernel optimization.
    """
    materials = self.getMaterials()
    averageAlbedo = sum([mat.getAlbedo() for mat in materials]) / len(materials)
    estimatedIPP = -np.log(weightThreshold) / averageAlbedo
    return estimatedIPP

Stats

Source code in pytissueoptics/rayscattering/statistics/statistics.py

class Stats:
    def __init__(self, logger: EnergyLogger):
        self._logger = logger
        self._pointCloudFactory = PointCloudFactory(logger)
        self._extractFromViews = not logger.has3D

        self._photonCount = logger.info["photonCount"]
        self._sourceSolidLabel = logger.info["sourceSolidLabel"]

        self._solidStatsMap = {}

    def report(self, solidLabel: str = None, saveToFile: str = None, verbose=True):
        if solidLabel and solidLabel not in self._logger.getSeenSolidLabels():
            utils.warn(f"WARNING: Cannot compute stats for solid '{solidLabel}' because it was not logged.")
            return

        self._computeStats(solidLabel)

        reportString = self._makeReport(solidLabel=solidLabel)
        if saveToFile:
            self._saveReport(reportString, saveToFile)
        if verbose:
            print(reportString)

    def _computeStats(self, solidLabel: str = None):
        solidLabels = [solidLabel]
        if solidLabel is None or utils.labelsEqual(solidLabel, WORLD_SOLID_LABEL):
            solidLabels = self._logger.getSeenSolidLabels()

        for solidLabel in solidLabels:
            if solidLabel == WORLD_SOLID_LABEL:
                continue
            try:
                absorbance = self.getAbsorbance(solidLabel)
            except ZeroDivisionError:
                utils.warn("WARNING: No energy input for solid '{}'".format(solidLabel))
                absorbance = None
            self._solidStatsMap[solidLabel] = SolidStats(
                absorbance,
                self.getAbsorbance(solidLabel, useTotalEnergy=True),
                self.getTransmittance(solidLabel),
                self._getSurfaceStats(solidLabel),
            )

    def _makeReport(self, solidLabel: str = None, reportString: str = ""):
        if solidLabel:
            if solidLabel == WORLD_SOLID_LABEL:
                reportString += self._reportWorld(solidLabel)
            else:
                reportString += self._reportSolid(solidLabel)
        else:
            for solidLabel in self._logger.getSeenSolidLabels():
                reportString = self._makeReport(solidLabel, reportString)
        return reportString

    def _reportWorld(self, worldLabel: str):
        totalSolidEnergy = sum([solidStats.totalAbsorbance for solidStats in self._solidStatsMap.values()])
        reportString = "Report of '{}'\n".format(worldLabel)
        reportString += "  Absorbed {:.2f}% of total power\n".format(100 - totalSolidEnergy)
        return reportString

    def _reportSolid(self, solidLabel: str):
        solidStats = self._solidStatsMap[solidLabel]
        reportString = "Report of solid '{}'\n".format(solidLabel)

        if solidStats.absorbance == math.inf:
            # Detectors don't log energy input.
            reportString += "  Detected {:.2f}% of total power\n".format(solidStats.totalAbsorbance)
            return reportString

        if solidStats.absorbance is None:
            reportString += "  Absorbance: N/A ({:.2f}% of total power)\n".format(solidStats.totalAbsorbance)
            reportString += "  Absorbance + Transmittance: N/A\n"
            return reportString

        reportString += "  Absorbance: {:.2f}% ({:.2f}% of total power)\n".format(
            solidStats.absorbance, solidStats.totalAbsorbance
        )
        reportString += "  Absorbance + Transmittance: {:.1f}%\n".format(
            solidStats.absorbance + solidStats.transmittance
        )

        for surfaceLabel, surfaceStats in solidStats.surfaces.items():
            reportString += "    Transmittance at '{}': {:.1f}%\n".format(surfaceLabel, surfaceStats.transmittance)

        return reportString

    def getAbsorbance(self, solidLabel: str, useTotalEnergy=False) -> float:
        if self._extractFromViews:
            return self._getAbsorbanceFromViews(solidLabel, useTotalEnergy)
        points = self._getPointCloud(solidLabel).solidPoints
        energyInput = self.getEnergyInput(solidLabel) if not useTotalEnergy else self.getPhotonCount()
        return 100 * self._sumEnergy(points) / energyInput if energyInput else math.inf

    def _getAbsorbanceFromViews(self, solidLabel: str, useTotalEnergy=False) -> float:
        energyInput = self.getEnergyInput(solidLabel) if not useTotalEnergy else self.getPhotonCount()
        absorbedEnergy = self._getAbsorbedEnergyFromViews(solidLabel)
        return 100 * absorbedEnergy / energyInput if energyInput else math.inf

    def _getAbsorbedEnergyFromViews(self, solidLabel: str) -> float:
        for view in self._logger.views:
            if not isinstance(view, View2DProjection):
                continue
            if not utils.labelsEqual(solidLabel, view.solidLabel):
                continue
            if not self._viewContainsSolid(view, solidLabel):
                continue
            return view.getSum()

        raise Exception(
            f"Could not extract absorbance for solid '{solidLabel}'. The 3D data was discarded and "
            f"no stored 2D view corresponds to this solid."
        )

    def _viewContainsSolid(self, view, solidLabel: str) -> bool:
        solidLimits = self._logger.getSolidLimits(solidLabel)
        requiredLimitsU = solidLimits[view.axisU]
        requiredLimitsV = solidLimits[view.axisV]
        if min(view.limitsU) > min(requiredLimitsU) or max(view.limitsU) < max(requiredLimitsU):
            return False
        if min(view.limitsV) > min(requiredLimitsV) or max(view.limitsV) < max(requiredLimitsV):
            return False
        return True

    def getPhotonCount(self) -> int:
        return self._photonCount

    def getEnergyInput(self, solidLabel: str = None) -> float:
        if solidLabel is None:
            return self.getPhotonCount()
        if self._extractFromViews:
            return self._getEnergyInputFromViews(solidLabel)
        points = self._getPointCloudOfSurfaces(solidLabel).enteringSurfacePoints
        energy = self._sumEnergy(points)

        if utils.labelsEqual(self._sourceSolidLabel, solidLabel):
            energy += self.getPhotonCount()
        return energy

    def _getEnergyInputFromViews(self, solidLabel: str) -> float:
        return self._getEnergyCrossingSolidFromViews(solidLabel, leaving=False)

    def _getEnergyLeavingFromViews(self, solidLabel: str):
        return self._getEnergyCrossingSolidFromViews(solidLabel, leaving=True)

    def _getEnergyCrossingSolidFromViews(self, solidLabel: str, leaving: bool) -> float:
        energy = 0
        for surfaceLabel in self._logger.getSeenSurfaceLabels(solidLabel):
            energy += self._getSurfaceEnergyFromViews(solidLabel, surfaceLabel, leaving=leaving)

        if utils.labelsEqual(self._sourceSolidLabel, solidLabel) and not leaving:
            energy += self.getPhotonCount()
        return energy

    def _getSurfaceEnergyFromViews(self, solidLabel: str, surfaceLabel: str, leaving: bool) -> float:
        for view in self._logger.views:
            if not utils.labelsEqual(solidLabel, view.solidLabel):
                continue
            if not utils.labelsEqual(surfaceLabel, view.surfaceLabel):
                continue
            if view.surfaceEnergyLeaving != leaving:
                continue
            if not self._viewContainsSolid(view, solidLabel):
                continue
            return view.getSum()
        raise Exception(
            f"Could not extract energy {['entering', 'leaving'][leaving]} surface '{surfaceLabel}' "
            f"of solid '{solidLabel}'. The 3D data was discarded and no stored 2D view corresponds "
            f"to this surface."
        )

    def _getSurfaceStats(self, solidLabel: str) -> Dict[str, SurfaceStats]:
        stats = {}
        for surfaceLabel in self._logger.getSeenSurfaceLabels(solidLabel):
            stats[surfaceLabel] = SurfaceStats(self.getTransmittance(solidLabel, surfaceLabel))
        return stats

    def getTransmittance(self, solidLabel: str, surfaceLabel: str = None, useTotalEnergy=False):
        """Uses local energy input for the desired solid by default. Specify 'useTotalEnergy' = True
        to compare instead with total input energy of the scene."""
        if self._extractFromViews:
            return self._getTransmittanceFromViews(solidLabel, surfaceLabel, useTotalEnergy)

        if surfaceLabel is None:
            points = self._getPointCloudOfSurfaces(solidLabel).leavingSurfacePoints
        else:
            points = self._getPointCloud(solidLabel, surfaceLabel).leavingSurfacePoints

        energyInput = self.getEnergyInput(solidLabel) if not useTotalEnergy else self.getPhotonCount()
        return 100 * self._sumEnergy(points) / energyInput if energyInput else math.inf

    def _getTransmittanceFromViews(self, solidLabel: str, surfaceLabel: str = None, useTotalEnergy=False):
        if surfaceLabel is None:
            energyLeaving = self._getEnergyLeavingFromViews(solidLabel)
        else:
            energyLeaving = self._getSurfaceEnergyFromViews(solidLabel, surfaceLabel, leaving=True)

        energyInput = self.getEnergyInput(solidLabel) if not useTotalEnergy else self.getPhotonCount()
        return 100 * energyLeaving / energyInput if energyInput else math.inf

    @staticmethod
    def _sumEnergy(points: np.ndarray):
        return np.abs(np.sum(points[:, 0])) if points is not None else 0

    def _getPointCloud(self, solidLabel: str = None, surfaceLabel: str = None) -> PointCloud:
        return self._pointCloudFactory.getPointCloud(solidLabel, surfaceLabel)

    def _getPointCloudOfSurfaces(self, solidLabel: str = None) -> PointCloud:
        return self._pointCloudFactory.getPointCloudOfSurfaces(solidLabel)

    @staticmethod
    def _saveReport(report: str, filepath: str = None):
        if filepath is None:
            filepath = "simulation_report"
            utils.warn(f"WARNING: No filepath specified. Saving to {filepath}.")
        i = 0
        filename, extension = filepath.split(".")
        if extension == "":
            extension = "txt"
        if os.path.exists(filepath):
            while os.path.exists("{}_{}.{}".format(filepath, i, extension)):
                i += 1
            filename = "{}_{}".format(filepath, i)
        filepath = "{}.{}".format(filename, extension)
        with open(filepath, "wb") as file:
            file.write(report.encode("utf-8"))
            file.close()

getTransmittance(solidLabel, surfaceLabel=None, useTotalEnergy=False)

Uses local energy input for the desired solid by default. Specify 'useTotalEnergy' = True to compare instead with total input energy of the scene.

Source code in pytissueoptics/rayscattering/statistics/statistics.py

def getTransmittance(self, solidLabel: str, surfaceLabel: str = None, useTotalEnergy=False):
    """Uses local energy input for the desired solid by default. Specify 'useTotalEnergy' = True
    to compare instead with total input energy of the scene."""
    if self._extractFromViews:
        return self._getTransmittanceFromViews(solidLabel, surfaceLabel, useTotalEnergy)

    if surfaceLabel is None:
        points = self._getPointCloudOfSurfaces(solidLabel).leavingSurfacePoints
    else:
        points = self._getPointCloud(solidLabel, surfaceLabel).leavingSurfacePoints

    energyInput = self.getEnergyInput(solidLabel) if not useTotalEnergy else self.getPhotonCount()
    return 100 * self._sumEnergy(points) / energyInput if energyInput else math.inf

ViewGroup

Bases: Flag

A ViewGroup is a bit Flag representing a group of 2D views. They can be combined with the | operator (bitwise OR).

It has two different usages

1. Specify which default views to create when initializing an EnergyLogger, particularly when keep3D is False.
2. Specify which views to show during visualization.

In the first use case, the default views created are three 2D projections of the whole scene (SCENE) or of each solid (SOLIDS) in every XYZ direction. Except for surface groups, where the default views created include a single 2D projection in the direction of the surface normal. SURFACES_ENTERING specifies the energy that entered the surface (energy direction opposite to the surface normal).

Source code in pytissueoptics/rayscattering/display/views/view2D.py

class ViewGroup(Flag):
    """
    A ViewGroup is a bit Flag representing a group of 2D views. They can be combined with the `|` operator (bitwise OR).

    It has two different usages:
        1. Specify which default views to create when initializing an EnergyLogger, particularly when `keep3D` is False.
        2. Specify which views to show during visualization.

    In the first use case, the default views created are three 2D projections of the whole scene (`SCENE`) or of each solid (`SOLIDS`) in every XYZ direction.
     Except for surface groups, where the default views created include a single 2D projection in the direction of the surface normal.
     `SURFACES_ENTERING` specifies the energy that entered the surface (energy direction opposite to the surface normal).
    """

    SCENE = 1
    SOLIDS = 2
    SURFACES_ENTERING = 4
    SURFACES_LEAVING = 8
    SURFACES = SURFACES_ENTERING | SURFACES_LEAVING
    ALL = SCENE | SOLIDS | SURFACES

Visibility

Bases: Flag

A Visibility is a bit Flag representing what to show inside a 3D visualization. They can be combined with the | operator (bitwise OR). AUTO will automatically switch to DEFAULT_3D if 3D data is present, else DEFAULT_2D.

Source code in pytissueoptics/rayscattering/display/viewer.py

class Visibility(Flag):
    """
    A Visibility is a bit Flag representing what to show inside a 3D visualization. They can be combined with the `|`
    operator (bitwise OR). `AUTO` will automatically switch to DEFAULT_3D if 3D data is present, else DEFAULT_2D.
    """

    SCENE = 1
    SOURCE = 2
    POINT_CLOUD = 4
    VIEWS = 8
    DEFAULT_3D = SCENE | SOURCE | POINT_CLOUD
    DEFAULT_2D = SCENE | SOURCE | VIEWS
    AUTO = 0