I cannot adjust the Hounsfield value of water in basic simulations

[juanenofbit]

I am performing a basic simulation and reconstruction experiment to check if the HU of water and air are set to 0 and -1000.

Using the most of the default configuration parameters, simulating a water cylinder 30 cm in diameter. To make it as simple as possible, no noise, monochromatic (so no beam-hardening correction is needed)

But I am getting a mean value in water is HU=4.4, and in air -992, why is this bias? Is there a simulation example where you can check the Hounsfield units in the reconstructions?

These are the some of the parameters changed to get a simplified simulation:

ct.physics.energyCount=1  #monochromatic
ct.protocol.bowtie = ""
ct.physics.enableQuantumNoise = 0
ct.physics.enableElectronicNoise = 0
ct.physics.monochromatic = 50 #  for monoenergetic specify the energy in keV
ct.protocol.spectrumFilename = []

Reconstruction attenuation from the retrieved value at 50keV:

from catsim.pyfiles.GetMu import GetMu
ct.recon.mu = GetMu('water', 50)
ct.recon.mu = ct.recon.mu[0]/10  # in mm^-1

The value is 0.022698049247264863 Then:

import reconstruction.pyfiles.recon as recon
ct.recon.kernelType='standard'
ct.recon.fov = 400
ct.recon.imageSize = 800

ct.do_Recon=True
ct.waitForKeypress=False
ct.recon.unit = 'HU'
ct.recon.huOffset = -1000
ct.recon.reconType = 'fdk_equiAngle'

recon.recon(ct)

I have noticed this commented line in the recon.py file, could this be the problem?

# A hack until the previous line is fixed.
#imageVolume3D = fdk_equiAngle(cfg, prep)
imageVolume3D = scaleReconData(cfg, imageVolume3D)

Thank you in advance!

[zhangjy-ge]

Hi, I was not able to reproduce your errors using a water phantom and monochromatic spectrum. Could you try to post your complete configurations (including phantom) so that I can run on my local machine? that will be easier for me to debug.

[juanenofbit]

Hi,

This is my complete configuration:

import catsim as xc ###------------ import catsim-xcist  version 0.1.8 
ct = xc.CatSim() # default parameters

# Modified parameters
ct.phantom.callback = 'Phantom_Analytic'              # name of function that reads and models phantom
ct.phantom.projectorCallback = 'C_Projector_Analytic' # name of function that performs projection through phantom
ct.phantom.filename = 'Cylinder_30cm.ppm'  

ct.physics.monochromatic = 50 #  for monoenergetic specify the energy in keV
ct.physics.energyCount = 1  #monochromatic
ct.physics.enableQuantumNoise    = 0
ct.physics.enableElectronicNoise = 0

ct.protocol.viewsPerRotation = 1440  # total numbers of view per rotation
ct.protocol.viewCount = 1440      
ct.protocol.stopViewId = ct.protocol.startViewId + ct.protocol.viewCount - 1 # index of the last view in the scan

ct.protocol.bowtie = []
ct.protocol.spectrumFilename = []
ct.protocol.mA = 500 
ct.protocol.flatFilter = []

ct.scanner.sid = 540.0   # source-to-iso distance (in mm)
ct.scanner.sdd = 950.0   # source-to-detector distance (in mm)
ct.scanner.detectionGain = 15.0  

ct.scanner.detectorColsPerMod = 1   # number of detector columns per module (4)
ct.scanner.detectorRowsPerMod = 1   # number of detector rows per module
ct.scanner.detectorColOffset = 0.0  # detector column offset relative to centered position (in detector columns)
ct.scanner.detectorColCount  = 900  # total number of detector columns
ct.scanner.detectorRowCount  = 1    # total number of detector rows
ct.resultsName ='Cylinder_30cm_monochromatic_github'

Running the simulation and reconstruct the sinogram

ct.run_all()

import reconstruction.pyfiles.recon as recon
from catsim.pyfiles.GetMu import GetMu
ct.recon.mu = GetMu('water', 50) # ct.physics.monochromatic = 50
ct.recon.mu = ct.recon.mu[0]/10

ct.recon.kernelType='standard'
ct.recon.fov        = 400
ct.recon.imageSize  = 800

ct.do_Recon=True
ct.waitForKeypress=False
ct.recon.unit = 'HU'
ct.recon.huOffset = -1000
ct.recon.reconType = 'fdk_equiAngle'

recon.recon(ct)

THen I load the reconstructed image: im = recon.loadImageVolume(ct)

The mean and std values of the regions inside the cylinder (avoiding border effects) and around the cylinder (air), as in the attached image are

Mean and std of water (red ROI): 
[4.5076714, 0.19668283]
Mean and std of air (yellow ROI): 
[-992.41547, 0.8617506]

The phantom file (Cylinder_30cm.ppm) has the following parameters:

materialList = {'water'};
object.center(1,:) = [0.000000 0.000000 0.000000];
object.half_axes(1,:) = [150.000000 150.000000 50.000000];
object.euler_angs(1,:) = [0.000000 0.000000 0.000000];
object.density(1) = 1.000000;
object.type(1) = 2;
object.material(1) = 1;
object.axial_lims(1,:) = [0 0];
object.shape(1) = 0;
object.clip{1} = [];

Cylinder_30cm.zip

[zhangjy-ge]

Just wanted to give you an update - we suspect this is a bug in our reconstruction code and we are still working on it. Thanks for bringing this up to us