| Home > Examples |
Loading...
Searching...
No Matches
The DICOM application has been originally developed by the Geant4 users:
Louis Archambault,(1)Luc Beaulieu, (2)Vincent Hubert-Tremblay.
And it has been deeply reviewed by Pedro Arce in December 2007.
Very small changes by Stephane Chauvie in January 2008.
Stephane Chauvie, Oct 2009: changed Physics list; changes in DICOM read.
Stephane Chauvie and Andrea Armando; June 2010 adapted for reading whatever DICOM file
Jonathan Madsen, Nov 2013: updated DICOM to utilize multithreading now available in Geant4.10
A new way to read DICOM files has been implemented since release 10.3, to avoid the often problems found by users when reading DICOM files. It can also read RT structures in DICOM format as well as RT plans. This utility uses the DCMTK (http://dicom.offis.de/dcmtk.php.en).
A new DICOM Digital Head included by S. Guatelli (susan.nosp@m.na@u.nosp@m.ow.ed.nosp@m.u.au) and V. Giacometti. Available since Geant4 10.4. The Digital model is documented in: Giacometti, V., Guatelli, S., Bazalova-Carter, M., Rosenfeld, A.B., Schulte, R.W., "Development of a high resolution voxelised head phantom for medical physics applications", (2017) Physica Medica, 33, pp. 182-188.
This example serves first to convert a DICOM file to a simple ASCII file, where the Hounsfield numbers are converted to materials and densities so that it can be used by GEANT4. It serves also to create a GEANT4 geometry based on the DICOM file information using the G4PhantomParameterisation.
You can find the phantom reproduced in the image PhantomCT.jpg. In the application the phantom is placed on a table.
You have to download and install DCMTK, see http://dicom.offis.de/dcmtk.php.en . At least you need version 3.6.1 or newer; the important thing is that you make sure it contains the dcmrt package. Then define the enviromental variable DCMTK_BASE_DIR to point to the directory where you have installed it. And to run it, you have to add to the enviromental variable LD_LIBRARY_PATH the path ${DCMTK_BASE_DIR}/lib (i.e. export LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${DCMTK_BASE_DIR}/lib )
Then you have to set the enviromental variable DICOM_USE_DCMTK to 1.
#DICOM_USE_DCMTK := true
'cmake -DGeant4_DIR=/path_to_geant4_install/lib/Geant4-x.x.x \
-DDCMTK_DIR=/path_to_dcmtk_install /path/to/DICOM/source'
then makeDICOM run.mac
DICOMthe file vis.mac is read in order to visualise the phantom with OpenGL, DAWN or VRML
The old version of "Data.dat" is found in "Data.dat.old", when the project is configured with DICOM_USE_DCMTK=OFF, "Data.dat.old" is copied into the binary directory at "Data.dat".
The new version of "Data.dat" is found in "Data.dat.new", when the project is configured with DICOM_USE_DCMTK=ON, "Data.dat.new" is copied into the binary directory as "Data.dat".
The file Data.dat has the following information
In case you want to convert DICOM files to text files, it must have the following lines:
As for the previous version, a Data.dat file has to be defined to manage the conversion options. The format of this file is though quite different from the previous version. The format of this file is based on tags (similary to the ASCII geometry files). The following tags should be used:
Where "level" is the number of voxels that will be merged into one in the X and Y dimen- sions. The Hounsfield numbers of the voxels merged are averaged to give the resulting value for the new voxel.
Example:
:COMPRESSION 4
4 X 4 voxels will be merged, so that the number of voxels in X and Y dimensions will be reduced by a factor 4
These are the list of files (one line per file) in DICOM format that will be treated. They can be of modality CT, RTSTRUCT or RTPLAN (the code will automatically detect its modality and treat it correspondingly).
Example:
:FILE 1.dcm :FILE 2.dcm :FILE 3.dcm
These sets of value pairs build the calibration curve (linearly interpolating between them). In other words, each Hounsfield number is given a material density using a function that is built interpolating between this list of value pairs.
Example:
:CT2D -5000 0. :CT2D -1000 0.01 :CT2D -400 0.602 :CT2D 300 1.145 :CT2D 2000 1.856
This serves for the Hounsfield number to material name conversion. The voxels with a Hounsfield number between 0. and the first upper bound will be assigned to the first material, those with a Hounsfiled number between the first upper bound and the second upper bound will be assigned to the second material, etc.
Example:
:MATE G4_AIR -800 :MATE G4_LUNG_ICRP -145 :MATE G4_ADIPOSE_TISSUE_ICRP -60 :MATE G4_WATER 0
Alternatively to the use of :MATE, you can use the :MATE_DENS
This serves for the material density to material name conversion. The voxels with a density between 0. and the first upper bound will be assigned to the first material, those with a density between the first upper bound and the second upper bound will be assigned to the second material, etc.
Example:
:MATE_DENS G4_AIR 0.207 :MATE_DENS G4_LUNG_ICRP 0.919 :MATE_DENS G4_ADIPOSE_TISSUE_ICRP 0.979 :MATE_DENS G4_WATER 1.01
We recommend the use of :MATE instead of :MATE_DENS as this is the way is used more often in the literature.
Name of output file containing the DICOM information in ASCII format
After reading the name of files from Data.dat, if a file .dcm is found, then it looks for the corresponding .g4dcm file and if not found creates it. Each file corresponds to a Z slice. The Z slices will be merged at runtime to form a unique patient volume; therefore the different slices have to be contiguous in Z.
The DICOM images pixel values represent CT (Hounsfield) numbers and they should be converted, first, to a given density and then to a material type. The relation between CT number and density is more or less linear. The file CT2Density.dat contains the calibration curve to convert CT (Hounsfield) number to physical density The assignment of material densities to materials is done following the information from the file Data.dat (see below). In this case we have used:
#####################################################
# Density Range Material #
#---------------------------------------------------#
# mg/cm3 - #
#---------------------------------------------------#
# [ 0. , 0.207 ) Air #
# [ 0.207 , 0.481 ) Lungs (inhale) #
# [ 0.481 , 0.919 ) Lungs (exhale) #
# [ 0.919 , 0.979 ) Adipose #
# [ 0.979 , 1.004 ) Breast #
# [ 1.004 , 1.043 ) Phantom #
# [ 1.043 , 1.109 ) Liver #
# [ 1.109 , 1.113 ) Muscle #
# [ 1.113 , 1.496 ) Trabecular Bone#
# [ 1.496 , 1.654 ] Dense Bone #
#####################################################
Data taken from the International Commission on Radiation Units and measurements (ICRU) report 46 was used to build the materials (lung, liver, breast, bones, ...)
In the class DicomDetectorConstruction, it is defined a density interval
G4double densityDiff = 0.1;
This means that the voxels of each material will be grouped in density intervals of 0.1 g/cm3 and a new material will be created for each group of voxels.
The file Colormap.dat defines the colour that will be assigned to the voxels of each material.
The DICOM files are converted to a simple text format. You may create your own file with the following format (see e.g. 14196616.g4dcm):
As commented before the DICOM files (.dcm) are assumed to describe one Z slice per file, and therefore the GEANT4 text files (.g4dcm) created from them have also one unique Z slice per file. Nevertheless if you create your own .g4dcm file you may include as many Z slices as desired. In any case you have to respect the rule that the Z slices must be contiguous.
The same information is also used to fill a file in binary format, that contains the same information as the text format. Its name ends in .g4dcmb, instead of .g4dcm .
There are four possible ways in GEANT4 to treat the navigation in regular voxelised volumes:
Obsolete option:
You can select among the four options in the following way:
patient_phys->SetRegularStructureId(0);
G4PVParameterised* patient_phys
= new G4PVParameterised("Patient",voxel_logic,container_logic,
kXAxis, nVoxelX*nVoxelY*nVoxelZ, param);
by G4PVParameterised * patient_phys
= new G4PVParameterised("Patient",voxel_logic,container_logic,
kUndefined, nVoxelX*nVoxelY*nVoxelZ, param);
Note also you must not set the enviromental variable DICOM_NESTED_PARAM.
As mentioned above the regular navigation has the option to keip voxel frontiers when two voxels share the same material, what can make the CPU time several times smaller. But this option makes that all energy deposited is computed in the last voxel, instead of distributing it along the voxels traversed. To properly calculate the dose in each voxel the G4PSDoseDeposit_RegNav scorer can be used.
It takes into account the fact that, when the particle travels through the voxels it looses energy and therefore the energy lost per length (dEdx) is bigger and also the effect of the multiple scattering is bigger. The algorithm to make this correction is an iterative one, as the step length increase due multiple scattering (that converts the geometrical step length in what we will call the true step length) and the energy loss are correlated. It works in the folloing way: first the total true step length is distributed among the voxels proportionally to their geometrical step length; with these values it is calculated one voxel after another the value of dEdx and then the value of the kinetic energy at the entrance of each voxel; with these values it is calculated the geometrical to true step corrections due to multiple scattering for each voxel; finally these new values are used to recalculate the energy lost in each voxel. It has been demonstrated for dose in a water phantom and in a real phantom that the two-step iteration described is enough to reproduce the dose calcualted when no skipping of voxel frontiers is done.
This scorer is implemented in this examples if the regular navigation option is chosen. It is triggered at the method RegularDicomDetectorConstruction::ConstructPhantom() by the call
SetScorer(voxel_logic);
For the nested parameterisation the geometry comprises replicas in X and Y which are then parameterised in Z. This means that to get the correct voxel idendification the replica depth has to be taken into account. The G4PSDoseDeposit3D scorers uses a fixed algorithm to calculate the voxel ID, according to the number of voxels in each axes and the associated replica depth.
G4PSDoseDeposit3D("DoseDeposit", fNoVoxelsZ, fNoVoxelsY, fNoVoxelsX, 0, 2, 1)
contains the number of voxels at the top (Z) level (0) and then two daughter levels down for the Y-voxels and one depth down for X.
dicom.out is produced running the macro file run.mac. It has 2 columns: the first is the number of voxel (ordered in x,y,z) and the second the dose there deposited (in Gy) It is produced, as an example, with a compression value of 32
It is possible to create a partial phantom, that is the intersection of a phantom with a volume. You may define the volume with the command
/dicom/intersectWithUserVolume 0. 0. 0. 45.*deg 0. 0. TUBE 0. 150. 100.
where the first three arguments are its position, its second three arguments are the rotation around the global X, Y and Z axis and the rest of the parameters are the same that you use to build a solid using the ASCII geometry format
Alternatively you can intersect the phantom with an existing Geant4 volume with the command
/dicom/intersectWithG4Volume VOLUME_NAME
The job will create an ASCII file names "phantom.g4pdcm" containing the partial phantom. To read this file all what is needed is to set the enviromental variable DICOM_PARTIAL_PARAM to 1
The Geant4 drivers are not meant for visualizing millions of voxel and visualising the DICOM geometries can be very computationally demanding. The users may want to visualise each DICOM slice separately or use higher compression values when visualising a part of DICOM project.
1.9.8