corrct.physics package
Submodules
corrct.physics.attenuation module
Incident beam and emidded radiation attenuation support.
@author: Nicola VIGANÒ, Computational Imaging group, CWI, The Netherlands, and ESRF - The European Synchrotron, Grenoble, France
- class corrct.physics.attenuation.AttenuationVolume(incident_local: ~numpy.ndarray[~typing.Any, ~numpy.dtype[~numpy.floating]] | None, emitted_local: ~numpy.ndarray[~typing.Any, ~numpy.dtype[~numpy.floating]] | None, angles_rot_rad: ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes], angles_det_rad: ~numpy.ndarray[~typing.Any, ~numpy.dtype[~numpy.floating]] | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes] = 1.5707963267948966, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float32'>)[source]
Bases:
object
Attenuation volume computation class.
- angles_det_rad: ndarray[Any, dtype[floating]]
- angles_rot_rad: ndarray[Any, dtype[floating]]
- compute_maps(use_multithreading: bool = True, verbose: bool = True) None [source]
Compute the correction maps for each angle.
- Parameters:
use_multithreading (bool, optional) – Use multi-threading for computing the attenuation maps. The default is True.
verbose (bool, optional) – Show verbose output. The default is True.
- dtype: dtype[Any] | None | type[Any] | _SupportsDType[dtype[Any]] | str | tuple[Any, int] | tuple[Any, SupportsIndex | Sequence[SupportsIndex]] | list[Any] | _DTypeDict | tuple[Any, Any]
- emitted_local: ndarray[Any, dtype[floating]] | None
- get_maps(roi: _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes] | None = None, rot_ind: int | slice | Sequence[int] | ndarray[Any, dtype[integer]] | None = None, det_ind: int | slice | Sequence[int] | ndarray[Any, dtype[integer]] | None = None) ndarray[Any, dtype[_ScalarType_co]] [source]
Return the attenuation maps.
- Parameters:
roi (ArrayLike, optional) – The region-of-interest to select. The default is None.
rot_ind (int, optional) – A specific rotation index, if only one is to be desired. The default is None.
det_ind (int, optional) – A specific detector index, if only one is to be desired. The default is None.
- Returns:
The attenuation maps.
- Return type:
NDArray
- get_projector_args(roi: _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes] | None = None, rot_ind: int | slice | Sequence[int] | ndarray[Any, dtype[integer]] | None = None, det_ind: int | slice | Sequence[int] | ndarray[Any, dtype[integer]] | None = None) dict[str, ndarray[Any, dtype[_ScalarType_co]]] [source]
Return the projector arguments.
- Parameters:
roi (ArrayLike, optional) – The region-of-interest to select. The default is None.
rot_ind (int, optional) – A specific rotation index, if only one is to be desired. The default is None.
det_ind (int, optional) – A specific detector index, if only one is to be desired. The default is None.
- Returns:
A dictionary containing the attenuation maps and the detector angle.
- Return type:
dict[str, NDArray]
- incident_local: ndarray[Any, dtype[floating]] | None
- maps: ndarray[Any, dtype[_ScalarType_co]]
- plot_map(ax: Axes, rot_ind: int, det_ind: int = 0, slice_ind: int | None = None, axes: Sequence[int] | ndarray[Any, dtype[integer]] = (-2, -1)) Sequence[float] [source]
Plot the requested attenuation map.
- Parameters:
ax (matplotlib axes) – The axes where to plot.
rot_ind (int) – Rotation angle index.
det_ind (int, optional) – Detector angle index. The default is 0.
slice_ind (Optional[int], optional) – Volume slice index (for 3D volumes). The default is None.
axes (Sequence[int] | NDArray, optional) – Axes of the slice. The default is (-2, -1).
- Returns:
The extent of the axes plot (min-max coords).
- Return type:
Sequence[float]
- Raises:
ValueError – In case a slice index is not passed for a 3D volume.
- vol_shape_zyx: ndarray[Any, dtype[_ScalarType_co]]
- corrct.physics.attenuation.get_linear_attenuation_coefficient(compound: str | dict, energy_keV: float, pixel_size_um: float, density: float | None = None) float [source]
Compute the linear attenuation coefficient for given compound, energy, and pixel size.
- Parameters:
compound (Union[str, dict]) – The compound for which we compute the linear attenuation coefficient
energy_keV (float) – The energy of the photons
pixel_size_um (float) – The pixel size in microns
density (Union[float, None], optional) – The density of the compound (if different from the default value), by default None
- Returns:
The linear attenuation coefficient
- Return type:
float
- corrct.physics.attenuation.plot_emission_line_attenuation(compound: str | dict, thickness_um: float, mean_energy_keV: float, fwhm_keV: float, line_shape: str = 'lorentzian', num_points: int = 201) None [source]
Plot spectral attenuation of a given line.
- Parameters:
compound (Union[str, dict]) – Compound to consider
thickness_um (float) – Thickness of the compound (in microns)
mean_energy_keV (float) – Average energy of the line
fwhm_keV (float) – Full-width half-maximum of the line
line_shape (str, optional) – Shape of the line, by default “lorentzian”. Options are: “gaussian” | “lorentzian” | “sech**2”.
num_points (int, optional) – number of discretization points, by default 201
- Raises:
ValueError – When an unsupported line is chosen.
- corrct.physics.attenuation.plot_transmittance_decay(compounds: str | dict | Sequence[str | dict], mean_energy_keV: float, thickness_range_um: tuple[float, float, int] = (0.0, 10.0, 101)) None [source]
Plot transmittance decay curve(s) for the given compound(s) at a given energy and thickness range.
- Parameters:
compounds (str | dict | Sequence[str | dict]) – The compound(s) description
mean_energy_keV (float) – The mean photon energy
thickness_range_um (tuple[float, float, int], optional) – The thickness range as (start, end, num_points), by default (0.0, 10.0, 101)
corrct.physics.materials module
Materials support functions and classes.
@author: Nicola VIGANÒ, ESRF - The European Synchrotron, Grenoble, France, and CEA-IRIG, Grenoble, France
- class corrct.physics.materials.VolumeMaterial(materials_fractions: Sequence[ndarray[Any, dtype[_ScalarType_co]]], materials_composition: Sequence, voxel_size_cm: float, dtype: dtype[Any] | None | type[Any] | _SupportsDType[dtype[Any]] | str | tuple[Any, int] | tuple[Any, SupportsIndex | Sequence[SupportsIndex]] | list[Any] | _DTypeDict | tuple[Any, Any] = None, verbose: bool = False)[source]
Bases:
object
VolumeMaterial class, that can be used for predicting fluorescence and Compton productions, attenuation, etc.
- Parameters:
material_fractions (Sequence) – Concentration fractions of each material for each voxel.
material_compounds (Sequence) – Compound description of each material.
voxel_size_cm (float) – Voxel size in cm.
dtype (DTypeLike, optional) – Data type of the produced data. The default is None.
- Raises:
ValueError – Raised in case of incorrect parameters.
- get_attenuation(energy_keV: float) ndarray[Any, dtype[_ScalarType_co]] [source]
Compute the local attenuation for each voxel.
- Parameters:
energy_keV (float) – The X-ray energy.
- Returns:
The computed local attenuation volume.
- Return type:
NDArray
- get_compton_scattering(energy_in_keV: float, angle_rad: float | None = None, detector: DetectorXRF | None = None) tuple[float, ndarray[Any, dtype[_ScalarType_co]]] [source]
Compute the local Compton scattering.
- Parameters:
energy_in_keV (float) – Incoming beam energy.
angle_rad (float, optional) – The detector angle, with respect to incoming beam direction. The default is None.
detector (DetectorXRF, optional) – The detector object. The default is None.
- Raises:
ValueError – In case neither of angle_rad nor detector have been passed.
- Returns:
energy_out_keV (float) – The energy of the Compton radiation received by the detector.
cmptn_prod (NDArray) – Local production of Compton radiation.
Either angle_rad or detector need to be supplied.
- get_element_mass_fraction(element: str | int) ndarray[Any, dtype[_ScalarType_co]] [source]
Compute the local element mass fraction through out all the materials.
- Parameters:
element (str | int) – The element to look for.
- Returns:
mass_fraction – The local mass fraction in each voxel.
- Return type:
NDArray
- get_fluo_production(element: str | int, energy_in_keV: float, fluo_lines: str | FluoLine | Sequence[FluoLine], detector: DetectorXRF | None = None) tuple[float, ndarray[Any, dtype[_ScalarType_co]]] [source]
Compute the local fluorescence production, for the given line of the given element.
Using Eq. (1) from: - T. Schoonjans et al., “The xraylib library for X-ray-matter interactions. Recent developments,” Spectrochim. Acta
Part B At. Spectrosc., vol. 66, no. 11-12, pp. 776-784, Nov. 2011, doi: 10.1016/j.sab.2011.09.011.
- Parameters:
element (str | int) – The element to consider.
energy_in_keV (float) – The incoming X-ray beam energy.
fluo_lines (str | FluoLine | Sequence[FluoLine]) – The fluorescence line to consider.
detector (DetectorXRF, optional) – The detector geometry. The default is None.
- Returns:
energy_out_keV (float) – The emitted fluorescence energy.
el_prod (NDArray) – The local fluorescence production in each voxel.
corrct.physics.phase module
Phase contrast support functions.
@author: Nicola VIGANÒ, CEA-IRIG, Grenoble, France
- corrct.physics.phase.apply_propagation_filter(data_wvu: ndarray[Any, dtype[_ScalarType_co]], pix_size_um: float, dist_um: float, wlength_um: float, delta_beta: float, filter_type: str = 'tie') ndarray[Any, dtype[_ScalarType_co]] [source]
Apply a requested propagation filter to an image or a stack of images.
- Parameters:
data_wvu (NDArray) – The Image or stack of images
pix_size_um (float) – Pixel size of the detector (in microns)
dist_um (float) – Propagation distance (in microns)
wlength_um (float) – Wavelength of the radiation (in microns)
delta_beta (float) – Delta-over-beta value for the given material
filter_type (str, optional) – Type of the filter, by default “tie”. Options: “ctf” (Contrast Transfer Function) | “tie” (Transport of Intensity Equation)
- Returns:
The filtered image
- Return type:
NDArray
- corrct.physics.phase.get_delta_beta(cmp_name: str, energy_keV: float, density: float | None = None) float [source]
Compute the delta-over-beta parameter for a specific compound.
- Parameters:
cmp_name (str) – Molar composition of a compound
energy_keV (float) – Energy at which the d/b value should be computed
density (Optional[float], optional) – Density of the compound, by default None
- Returns:
The computed delta-over-beta.
- Return type:
float
- corrct.physics.phase.get_delta_beta_curves(compounds: Sequence[str], energy_keV_range: tuple[float, float, int] = (1.0, 800.0, 500), plot: bool = True) Sequence[ndarray[Any, dtype[_ScalarType_co]]] [source]
Compute the delta-over-beta curves for the listed compounds, in the requested energy range.
- Parameters:
compounds (Sequence[str]) – Sequence of compounds.
energy_keV_range (tuple[float, float, int], optional) – The energy range in keV, by default (1.0, 800.0, 500)
plot (bool, optional) – Whether to plot the result, by default True
- Returns:
The computed delta-over-beta valueas.
- Return type:
Sequence[NDArray]
- corrct.physics.phase.get_propagation_filter(img_shape: ~collections.abc.Sequence[int] | ~numpy.ndarray[~typing.Any, ~numpy.dtype[~numpy._typing._array_like._ScalarType_co]], pix_size_um: float, dist_um: float, wlength_um: float, delta_beta: float, filter_type: str = 'ctf', use_rfft: bool = False, plot_result: bool = False, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float32'>) tuple[ndarray[Any, dtype[_ScalarType_co]], ndarray[Any, dtype[_ScalarType_co]]] [source]
Compute the phase contrast propagation filter for the given parameters.
- Parameters:
img_shape (Union[Sequence[int], NDArray]) – Shape of the target image
pix_size_um (float) – Pixel size of the detector (in microns)
dist_um (float) – Propagation distance (in microns)
wlength_um (float) – Wavelength of the radiation (in microns)
delta_beta (float) – Delta-over-beta value for the given material
filter_type (str, optional) – Type of the filter, by default “ctf”. Options: “ctf” (Contrast Transfer Function) | “tie” (Transport of Intensity Equation)
use_rfft (bool, optional) – Whether to use the rfft, by default False
plot_result (bool, optional) – Whether to plot the result, by default False
dtype (DTypeLike, optional) – Data type to use, by default np.float32
- Returns:
The Fourier-space and real-space filters, respectively
- Return type:
tuple[NDArray, NDArray]
- Raises:
ValueError – When choosing an incorrect filter type
- corrct.physics.phase.plot_filter_responses(filter_length: int, pix_size_um: float, dist_um: float, wlength_um: float, delta_beta: float, domain: str = 'fourier') tuple [source]
Plot frequency response of the wave propagation.
- Parameters:
filter_length (int) – Length of the filter in pixels
pix_size_um (float) – Pixel-wise of the detector in microns
dist_um (float) – Distance of the detector from sample in microns
wlength_um (float) – Wavelength of the wave in microns
delta_beta (float) – Radio between the refraction index decrement and the absorption coefficient.
domain (str) – Whether to plot Fourier or direct-space responses. By default, “Fourier”.
- Returns:
Figure and axes of the plot
- Return type:
tuple
corrct.physics.xraylib_helper module
xraylib handling functions.
@author: Nicola VIGANÒ, CEA-IRIG, Grenoble, France
- corrct.physics.xraylib_helper.get_compound(cmp_name: str, density: float | None = None) dict [source]
Build a compound from the compound composition string.
- Parameters:
cmp_name (str) – Compund name / composition.
density (float, optional) – The density of the compound. If not provided it will be approximated from the composing elements. The default is None.
- Returns:
cmp – The compound structure.
- Return type:
dict
- corrct.physics.xraylib_helper.get_compound_cross_section(compound: dict, mean_energy_keV: float) float [source]
Compute a compound cross section for the given incoming photon energy.
- Parameters:
compound (dict) – The compound structure (as returned by get_compound)
mean_energy_keV (float) – The average photon energy
- Returns:
The computed cross-section
- Return type:
float
- corrct.physics.xraylib_helper.get_element_number(element: str | int) int [source]
Return the element number from the symbol.
- Parameters:
element (str | int) – The element symbol (or number, which won’t be converted).
- Returns:
The corresponding element number.
- Return type:
int
- corrct.physics.xraylib_helper.get_element_number_and_symbol(element: str | int) tuple[str, int] [source]
Return both the element symbol and number from either the symbol or the number.
- Parameters:
element (str | int) – The element symbol (or number, which won’t be converted).
- Returns:
The element symbol and number.
- Return type:
tuple[str, int]
corrct.physics.xrf module
XRF support functions.
@author: Nicola VIGANÒ, ESRF - The European Synchrotron, Grenoble, France, and CEA-IRIG, Grenoble, France
- class corrct.physics.xrf.DetectorXRF(surface_mm2: float, distance_mm: float | ndarray[Any, dtype[_ScalarType_co]], angle_rad: float = 1.5707963267948966)[source]
Bases:
object
Simple XRF detector model.
- angle_rad: float = 1.5707963267948966
- distance_mm: float | ndarray[Any, dtype[_ScalarType_co]]
- property solid_angle_sr: float | ndarray[Any, dtype[_ScalarType_co]]
Compute the solid angle covered by the detector.
- Returns:
The computed solid angle of the detector geometry.
- Return type:
float | NDArray
- surface_mm2: float
- class corrct.physics.xrf.FluoLine(name: str, indx: int)[source]
Bases:
object
Fluorescence line description class.
- indx: int
- name: str
- class corrct.physics.xrf.LinesSiegbahn[source]
Bases:
object
Siegbahn fluorescence lines collection class.
- static get_lines(line: str) Sequence[FluoLine] [source]
Return the list of xraylib line macro definitions for the requested family.
- Parameters:
line (str) – The requested line. It can be a whole shell (transition to that shell), or sub-shells.
- Returns:
List of corresponding lines.
- Return type:
Sequence
- lines = [FluoLine(name='KA1', indx=-3), FluoLine(name='KA2', indx=-2), FluoLine(name='KA3', indx=-1), FluoLine(name='KB1', indx=-6), FluoLine(name='KB2', indx=-11), FluoLine(name='KB3', indx=-5), FluoLine(name='KB4', indx=-13), FluoLine(name='KB5', indx=-8), FluoLine(name='LA1', indx=-90), FluoLine(name='LA2', indx=-89), FluoLine(name='LB1', indx=-63), FluoLine(name='LB2', indx=-95), FluoLine(name='LB3', indx=-34), FluoLine(name='LB4', indx=-33), FluoLine(name='LB5', indx=-102), FluoLine(name='LB6', indx=-91), FluoLine(name='LB7', indx=-98), FluoLine(name='LB9', indx=-36), FluoLine(name='LB10', indx=-35), FluoLine(name='LB15', indx=-94), FluoLine(name='LB17', indx=-62), FluoLine(name='LG1', indx=-68), FluoLine(name='LG2', indx=-38), FluoLine(name='LG3', indx=-39), FluoLine(name='LG4', indx=-47), FluoLine(name='LG5', indx=-65), FluoLine(name='LG6', indx=-75), FluoLine(name='LG8', indx=-72), FluoLine(name='LE', indx=-60), FluoLine(name='LH', indx=-60), FluoLine(name='LL', indx=-86), FluoLine(name='LS', indx=-88), FluoLine(name='LT', indx=-87), FluoLine(name='LU', indx=-96), FluoLine(name='LV', indx=-70), FluoLine(name='MA1', indx=-207), FluoLine(name='MA2', indx=-206), FluoLine(name='MB', indx=-187), FluoLine(name='MG', indx=-165)]
- corrct.physics.xrf.get_energy(element: str | int, lines: str | FluoLine | Sequence[FluoLine], *, compute_average: Literal[False] = False, verbose: bool = False) ndarray[Any, dtype[_ScalarType_co]] [source]
- corrct.physics.xrf.get_energy(element: str | int, lines: str | FluoLine | Sequence[FluoLine], *, compute_average: Literal[True] = True, verbose: bool = False) float
Return the energy(ies) of the requested line for the given element.
- Parameters:
element (Union[str, int]) – The requested element.
line (str) – The requested line. It can be a whole shell (transition to that shell), or sub-shells.
compute_average (bool, optional) – Weighted averaging the lines, using the radiation rate. The default is False.
- Returns:
energy_keV – Either the average energy or the list of different energies.
- Return type:
Union[float, NDArray]
- corrct.physics.xrf.get_radiation_rate(element: str | int, lines: str | FluoLine | Sequence[FluoLine], verbose: bool = False) ndarray[Any, dtype[_ScalarType_co]] [source]
Return the radiation rates of the requested lines for the given element.
- Parameters:
- Returns:
The list of radiation rates
- Return type:
NDArray
Module contents
Physics module.
- corrct.physics.convert_energy_to_wlength(energy_keV: float | ndarray[Any, dtype[_ScalarType_co]]) float | ndarray[Any, dtype[_ScalarType_co]] [source]
Convert from energy in keV to wavelength in meters.
- Parameters:
energy_keV (float | NDArray) – The energy in keV
- Returns:
The wavelength in m
- Return type:
float | NDArray
- corrct.physics.convert_wlength_to_energy(wlength_m: float | ndarray[Any, dtype[_ScalarType_co]]) float | ndarray[Any, dtype[_ScalarType_co]] [source]
Convert wavelength in meters to energy in keV.
- Parameters:
wlength_m (float | NDArray) – The wavelength in meter
- Returns:
The energy in keV
- Return type:
float | NDArray
- corrct.physics.pencil_beam_profile(voxel_size_um: float, beam_fwhm_um: float, profile_size: int = 1, beam_shape: str = 'gaussian', verbose: bool = False) ndarray[Any, dtype[_ScalarType_co]] [source]
Compute the pencil beam integration point spread function.
- Parameters:
voxel_size_um (float) – The integration length.
beam_fwhm_um (float) – The beam FWHM.
profile_size (int, optional) – The number of pixels of the PSF. The default is 1.
verbose (bool, optional) – Whether to print verbose information. The default is False.
- Returns:
The beam PSF.
- Return type:
NDArray