Scattering
To do scattering, no new classes are made, only functions relying on all the previous classes are used.
Convolutional Scattering
This type of scatter uses the same method as the scatter presented in the XCIST paper, and the kernel for the convolution has been retrieved from XCIST GitHub.
- convolutional_scatter(xray_source, photon_count, detector, sfactor)
This method, takes in the photon count across all energy bins and applies a convolution to this array, the result is then added to the original photon count array to create the final count.
- Parameters:
xray_source (
source) – The xray_source used to create the photon count array.photon_count (
double) – An array of dimensions [energy_bins, ny_pix, nz_pix, nrotation], generated by ray tracing the source through the phantom.detector (
detector) – The detector used to create the photon count array.sfactor (
double) – A scatter factor used to customise the amount of scatter in the image.
- Returns:
scatter – The additional scatter to be added to the photon count array.
- Return type:
nbins x ny x nz x nrot double
To retrieve the scatter kernel, the user can use the following function.
- get_scatter_kernel()
Since the data is retrieved from the file taken from the xcist repository, the data needs to be reshaped to column major order, as the data is stored in row major order.
- Returns:
scatter_kernel – The scatter kernel used in the convolutional scatter.
- Return type:
65x49 double
To normalise the photon count, we calculate what the count would be if there was no phantom, an air scan. This function is not limited to the convolutional scatter, it can be used by the user to calculate the air scan, if they so wish.
- air_scan(xray_source, detector)
This function assumes that the ray length provided by the detector array is the exact length between the source and the detector, therefore calculates the total attenuation without performing any ray tracing, using \(\mu_a\,l\), where \(\mu_a\) is the linear attenuation coefficient of air and \(l\) is the length of the ray.
- Parameters:
- Returns:
photon_count – An array of dimensions [energy_bins, ny_pix, nz_pix, nrotation], representing the number of photons that would be detected if there was no phantom.
- Return type:
nbins x ny x nz x nrot double
Deterministic Scatter
This scatter method is still under development and may be subject to change. The method is takes a more deterministic approach to scatter than true monte carlo methods, where we rely on the rays generated for ray tracing, and then sample the scatter from these rays.
We sample every ray at 0.1 mean free paths, and then trace 100*sfactor rays from this point, randomly choosing the direction of the scatter, and the energy of the scattered ray. We then calculate if this ray will hit the detector, and if it does, we add the number of photons that will hit this pixel to the scatter array.
We determine the number of photons that will hit the pixel by using the attenuation of the up to this sample point, the probability of the ray scattering at this point, the number of rays being sampled from this point, and the attenuation of the ray from this point to the detector.
There are several functions involved in this method, but the main function is the following.
- monte_carlo_scatter(xray_source, phantom, detector_obj, sfactor)
This function is very similar to
compute_sinogram(), with the only real difference being that this function is used to calculate the scatter, and not the sinogram. This also means that the energy that the rays are scattered with is not the same as the energy that the rays are generated with, and so we need to keep track of this.- Parameters:
xray_source (
source) – The source of xrays to be used to create the scatter.phantom (
voxel_array) – The phantom that the rays will be scattered through.detector_obj (
detector) – The detector that the rays will be scattered to.sfactor (
double) – A scatter factor used to customise the amount of rays traced and potentially scattered.
- Returns:
scatter – The scatter array to be added to the photon count array.
- Return type:
nbins x ny x nz x nrot double
- compton_dist(nrjs)
- Parameters:
nrjs (
1xN double) – A list of energies to be used to sample the scatter.- Returns:
thetas – A list of thetas produced by sampling the distribution using the Geant4 physics reference manual. Used alongside the
compton_scatter()function. The produced thetas should be randomly sampled to produce the angles of scatter.- Return type:
1xN double
- compton_scatter(direction, inrj, thetas)
This function calculates the new direction and energy of the ray after a Compton scatter. The sampling of the new direction is done using the formulae from the Geant4 physics reference manual and assisted by the Geant4 source code.
This function uses some code directly from CLHEP to calculate the vector transformation in order to retrieve the new direction of the ray.
- Parameters:
direction (
3x1 double) – The initial direction of the ray.inrj (
double) – The initial energy of the ray.thetas (
1xN double) – A list of thetas to be used to sample the scatter.
- Returns:
direction (
3x1 double) - The new direction of the ray after scattering.nrj (
double) - The new energy of the ray after scattering.