pytomography.io.SPECT.helpers
#
This module contains helper function for inputting/outputting/interpretting SPECT/CT data in the SIMIND and DICOM file formats. A considerable amount of these functions have to do with creating attenuation maps from CT data files for attenuation correction in SPECT imaging.
Module Contents#
Functions#
|
|
|
Function used for curve fitting of linear attenuation coefficient vs. photon energy curves from NIST. It's given by the functional form \(f(x) = c_1e^{-d_1\sqrt{x}} + c_2e^{-d_2\sqrt{x}}\). It was chosen purely because it gave good fit results. |
|
Function used to convert between Hounsfield Units at an effective CT energy and linear attenuation coefficient at a given SPECT radionuclide energy. It consists of two distinct linear curves in regions \(HU<0\) and \(HU \geq 0\). |
|
Given a list of seperate DICOM files, opens them up and stacks them together into a single CT image. |
Return energy and linear attenuation data from NIST datafiles of mass attenuation coefficients between 50keV and 511keV. |
|
|
Gets linear attenuation corresponding to a given energy in in the data from |
|
Obtains the Hounsfield Units of some material at a given energy |
|
Gets Housnfield Unit vs. linear attenuation coefficient for air, water, and cortical bone data points |
|
Converts hounsfield units to linear attenuation coefficient |
|
Obtains the Hounsfield Unit corresponding to cortical bone from a CT scan. |
Finds the effective CT energy that gives a cortical bone Hounsfield Unit value corresponding to |
|
|
Obtains the HU to mu conversion function that converts CT data to the required linear attenuation value in units of 1/cm required for attenuation correction in SPECT/PET imaging. |
Attributes#
- pytomography.io.SPECT.helpers.compute_TEW(projection_lower, projection_upper, width_lower, width_upper, width_peak)[source]#
- pytomography.io.SPECT.helpers.dual_sqrt_exponential(energy, c1, c2, d1, d2)[source]#
Function used for curve fitting of linear attenuation coefficient vs. photon energy curves from NIST. It’s given by the functional form \(f(x) = c_1e^{-d_1\sqrt{x}} + c_2e^{-d_2\sqrt{x}}\). It was chosen purely because it gave good fit results.
- Parameters:
energy (float) – Energy of photon
c1 (float) – Fit parameter 1
c2 (float) – Fit parameter 2
d1 (float) – Fit parameter 3
d2 (float) – Fit parameter 4
- Returns:
_description_
- Return type:
float
- pytomography.io.SPECT.helpers.bilinear_transform(HU, a1, a2, b1, b2)[source]#
Function used to convert between Hounsfield Units at an effective CT energy and linear attenuation coefficient at a given SPECT radionuclide energy. It consists of two distinct linear curves in regions \(HU<0\) and \(HU \geq 0\).
- Parameters:
HU (float) – Hounsfield units at CT energy
a1 (float) – Fit parameter 1
a2 (float) – Fit parameter 2
b1 (float) – Fit parameter 3
b2 (float) – Fit parameter 4
- Returns:
Linear attenuation coefficient at SPECT energy
- Return type:
float
- pytomography.io.SPECT.helpers.open_CT_file(files_CT, return_max_slice_loc=False)[source]#
Given a list of seperate DICOM files, opens them up and stacks them together into a single CT image.
- Parameters:
files_CT (Sequence[str]) – List of DICOM files corresponding to a particular scan
return_max_slice_loc (bool, optional) – Whether or not to return the max slice location of all the CT files. Needed when applying affine transformation to align CT and SPECT images. Defaults to False.
- Returns:
CT scan in units of Hounsfield Units at the effective CT energy.
- Return type:
np.array
- pytomography.io.SPECT.helpers.get_E_mu_data_from_datasheet(file)[source]#
Return energy and linear attenuation data from NIST datafiles of mass attenuation coefficients between 50keV and 511keV.
- Parameters:
file (str) – Location of NIST data file. Corresponds to a particular element/material.
- Returns:
Energy and linear attenuation values.
- Return type:
tuple[np.array, np.array]
- pytomography.io.SPECT.helpers.get_mu_from_spectrum_interp(file, energy)[source]#
Gets linear attenuation corresponding to a given energy in in the data from
file
.- Parameters:
file (str) – Filepath of the mu-energy data
energy (float) – Energy at which mu is computed
- Returns:
Linear attenuation coefficient (in 1/cm) at the desired energies.
- Return type:
np.array
- pytomography.io.SPECT.helpers.get_HU_from_spectrum_interp(file, energy)[source]#
Obtains the Hounsfield Units of some material at a given energy
- Parameters:
file (str) – Filepath of material
energy (float) – Energy at which HU is desired
- Returns:
HU at the desired energies.
- Return type:
np.array
- pytomography.io.SPECT.helpers.get_HU_mu_curve(files_CT, CT_kvp, E_SPECT)[source]#
Gets Housnfield Unit vs. linear attenuation coefficient for air, water, and cortical bone data points
- Parameters:
files_CT (Sequence[str]) – Filepaths of all CT slices
CT_kvp (float) – Value of kVp for the CT scan
E_SPECT (float) – Photopeak energy of the SPECT scan
- Return type:
tuple[np.array, np.array]
- pytomography.io.SPECT.helpers.HU_to_mu(HU, E, p_water_opt, p_air_opt)[source]#
Converts hounsfield units to linear attenuation coefficient
- Parameters:
HU (float) – Hounsfield Unit value
E (float) – Effective CT energy
p_water_opt (Sequence[float]) – Optimal fit parameters for mu vs. E data for water
p_air_opt (Sequence[float]) – Optimal fit parameters for mu vs. E data for air
- Returns:
_description_
- Return type:
_type_
- pytomography.io.SPECT.helpers.get_HU_corticalbone(files_CT)[source]#
Obtains the Hounsfield Unit corresponding to cortical bone from a CT scan.
- Parameters:
files_CT (Sequence[str]) – CT data files
- Returns:
Hounsfield unit of bone. If not found, then returns
None
.- Return type:
float | None
- pytomography.io.SPECT.helpers.get_ECT_from_corticalbone_HU(HU)[source]#
Finds the effective CT energy that gives a cortical bone Hounsfield Unit value corresponding to
HU
.- Parameters:
HU (float) – Hounsfield Unit of Cortical bone at effective CT energy
- Returns:
Effective CT energy
- Return type:
float
- pytomography.io.SPECT.helpers.get_HU2mu_conversion(files_CT, CT_kvp, E_SPECT)[source]#
Obtains the HU to mu conversion function that converts CT data to the required linear attenuation value in units of 1/cm required for attenuation correction in SPECT/PET imaging.
- Parameters:
files_CT (Sequence[str]) – CT data files
CT_kvp (float) – kVp value for CT scan
E_SPECT (float) – Energy of photopeak in SPECT scan
- Returns:
Conversion function from HU to mu.
- Return type:
function