pytomography.io.CT.attenuation_map#

This module is used to create attenuation maps from CT images required for SPECT/PET attenuation correction.

Module Contents#

Functions#

bilinear_transform(HU, a1, a2, b1, b2)

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\).

get_HU_from_spectrum_interp(file, energy)

Obtains the Hounsfield Units of some material at a given energy

get_HU_mu_curve(files_CT, CT_kvp, E_SPECT)

Gets Housnfield Unit vs. linear attenuation coefficient for air, water, and cortical bone data points

HU_to_mu(HU, E, p_water_opt, p_air_opt)

Converts hounsfield units to linear attenuation coefficient

get_HU_corticalbone(files_CT)

Obtains the Hounsfield Unit corresponding to cortical bone from a CT scan.

get_ECT_from_corticalbone_HU(HU)

Finds the effective CT energy that gives a cortical bone Hounsfield Unit value corresponding to HU.

get_HU2mu_conversion(files_CT, CT_kvp, E_SPECT)

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#

module_path

DATADIR

FILE_WATER

FILE_AIR

FILE_CBONE

pytomography.io.CT.attenuation_map.module_path[source]#
pytomography.io.CT.attenuation_map.DATADIR[source]#
pytomography.io.CT.attenuation_map.FILE_WATER[source]#
pytomography.io.CT.attenuation_map.FILE_AIR[source]#
pytomography.io.CT.attenuation_map.FILE_CBONE[source]#
pytomography.io.CT.attenuation_map.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.CT.attenuation_map.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.CT.attenuation_map.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.CT.attenuation_map.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.CT.attenuation_map.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.CT.attenuation_map.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.CT.attenuation_map.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