How to do Ray Tracing in 3D images, such as CT images?

Support
1

Hi Dear 3DSlicer developers

I want to calculate the length of any ray that propagates from a point source and passes through the CT voxels. In other words, I want to calculate the length of each ray along its path through the voxels. So for each ray I want to calculate two parameters:

  1. The index number of the ray that the ray passes through.

  2. The distance traveled within the voxel.

For a better explanation and understanding, please look at the following figure.

Picture3
Picture31207×1041 76.1 KB

In the above figure, it is assumed that the rays emanated from of a point source are equal to the number of pixels in the imaging system. These rays are lines drawn from the point source to the center of pixel. So if the image matrix is 1024x1024, then there will be 1048576 rays. The figure below is a demonstration of this concept in the sagital slice (2D).

Picture2
Picture21349×1003 60.9 KB

In this particular ray, how can I determine it passes which voxels? What distance has it traveled from each of these voxels? The figure below is the same concept in three dimensions.

Picture1
Picture1818×993 56.5 KB

What is your proposed solution to perform these calculations?
Does VTK have a library to do it?
Please guide me.
Shahrokh

2

Excuse me, I must correct the above message:

… for each ray I want to calculate two parameters:

  1. The index number of the voxel that the ray passes through.
  2. The distance traveled within the voxel.
    …
3

Honestly you should just do this in plastimatch. Use plastimatch drr with the --raytrace-details option. It will create a large file with exactly what you want.

4

I think this problem has been solved.

Take a look at this journal paper: Siddon, Robert L. “Fast calculation of the exact radiological path for a three‐dimensional CT array.” Medical physics 12.2 (1985): 252-255.
and its implementation at : https://code.google.com/archive/p/drrsuite/downloads

There are also two VTK functions that may help:
vtkFixedPointVolumeRayCastMapper
vtkGPUVolumeRayCastMapper

If they do not work for you, you can check some journal papers on Digitally Reconstructed Radiograph (DRR). I saw that they were trying to solve the same problem as yours.

Hope this help.
Manh

5

Hi Dear 3DSlicer and Plastimatch developers,
I apologize for bringing up my question again after 5 years.
Let’s assume that the file stack.mha has the following specifications:

$ plastimatch header stack.mha
Type = short
Planes = 1
Origin = -249.5117 -459.5117 -921.4000
Size = 512 512 72
Spacing = 0.9766 0.9766 3.0000
Direction = 1.0000 0.0000 0.0000 0.0000 1.0000 0.0000 0.0000 0.0000 1.0000
$

Now, I have calculated DRR image and the related ray tracing information using the following command.

$ plastimatch drr --nrm "0 -1 0" --vup "0 0 1" --sid 940 --sad 535 -r "512 512" -z "899.588785047 899.588785047" -c "256 256" -o "-0.4787000 -210.4787 -816.4000" --raytrace-details -I stack.mha -i uniform -O DRRSample -t pfm
I/O time: 0.000524 sec
Total time: 59.2702 secs
$ ll -h
...
-rw-rw-r-- 1 sn sn 1.1M Aug  4 13:33  DRRSample0000.pfm
-rw-rw-r-- 1 sn sn  913 Aug  4 13:33  DRRSample0000.txt
-rw-rw-r-- 1 sn sn 2.6G Aug  4 13:33  -I0000.txt
-rw-rw-r-- 1 sn sn  17M Jul  6 11:48  stack.mha
$ vi -I0000.txt
...

The ray tracing information has been saved in the file -I0000.txt. The contents of this file are extensive, as follows:

51,0,0,0,0,0.732422,0,0
51,1,0,0,0,0.732422,0,0
51,2,0,0,0,0.732422,0,0
51,3,0,0,0,0.732422,0,0
51,4,0,0,0,0.732422,0,0
51,5,0,0,0,0.732422,0,0
51,6,0,0,0,0.732422,0,0
51,7,0,0,0,0.732422,0,0
51,8,0,0,0,0.732422,0,0
51,9,0,0,0,0.732422,0,0
...
455,503,1,0,0,0.732422,0,0
455,504,0,0,0,0.732422,0,0
455,504,1,0,0,0.732422,0,0
455,505,0,0,0,0.732422,0,0
455,505,1,0,0,0.732422,0,0
455,506,0,0,0,0.732422,0,0
455,506,1,0,0,0.732422,0,0
455,507,0,0,0,0.732422,0,0
455,507,1,0,0,0.732422,0,0
455,508,0,0,0,0.732422,0,0
455,508,1,0,0,0.732422,0,0
455,509,0,0,0,0.732422,0,0
455,509,1,0,0,0.732422,0,0
455,510,0,0,0,0.732422,0,0
455,510,1,0,0,0.732422,0,0
455,511,0,0,0,0.732422,0,0
455,511,1,0,0,0.732422,0,0
$

Please explain the eight parameters of each row for me.

I imagine the first parameter corresponds to the Ray ID. If my understanding is correct, why does it start from number 51 and end at 455? Similarly, what are the other parameters?

My ultimate goal is to use this ray tracing information to track each ray passing through a 3D CT matrix. I intend to determine the Hounsfield number of each voxel (from the ray tracing data), then calculate the angle of incidence of the ray with that voxel. Subsequently, based on the Hounsfield number, I will adjust the scattering kernel computed from Monte Carlo simulation according to the calculated angle, and apply it to the voxel position.

Another objective is to utilize this information to identify edge voxels of the RTPlan radiation field.

Please guide me on these matters.
Your sincerely.
Shahrokh