ray tracing, indices and path length registration

Dear all

Thank you for reading my following question and I greatly appreciate any of your help.

Problem: Given a grid (N*N) and a non-straight ray with known discrete coordinates, how to register the indices and path length in the pixels crossed by the ray? This is on a plane, what if it is in 3d space? Are there any existing MATLAB script to carry out this function?

Thank you very much.

Sincerely

Albert

[Jan]

What is a non-straight ray?

[Albert]

Thank you Jan.

What I meant by non-straight ray could be zigzag line or curve. The ray is represented by many dots with known coordinates.

[Corey Silva]

Is this related to your question?

https://www.mathworks.com/matlabcentral/fileexchange/45355-optometrika

[Albert]

Thank you Corey, although it's not what I need. I am doing muon tomography. We have already found way to register muon most likely trajectory represented by discrete points. Now I am focusing on a function which could register the indexes and segment length in these voxels/pixels intersected by the trajectory. Here is a function

https://www.mathworks.com/matlabcentral/fileexchange/56527-fast-raytracing-through-a-3d-grid

It's quite close to what I need, but it doesn't register segment length in these voxels/pixels. I am modifying this function now.

Thank you

[Alexander Jones]

Hello, I am trying to do some microwave tomography, and am in need of a very similar thing - I was wondering if you found a solution to the problem of finding path length, and if so if you could please point me in the right direction.

[Albert]

Hi Jones,

To be clear, are you trying to trace a straight or curve? If it is for straight, please see https://www.mathworks.com/matlabcentral/fileexchange/28496-tomobox?focused=5159748&tab=function

Anyhow you may want to have a look at this document for the theory behind the code, see here (much efficient than calculating the interaction) http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.42.3443&rep=rep1&type=pdf

If you are looking at curve rather than straight ray, as Anton Semechko mentioned, the curve ray may be represented by segments, for each segment you can use the above code to calculate index and path length, but you need truncate the result before you concatenate them together for one ray. Because the above code calculates the index and path length in each voxel for straight line across the whole image reconstruction volume.

The above code can handle 2D and 3D.

Good luck.

Answer 1

[Anton Semechko]

For 2D reconstruction, the problem can be approached follows:

1) Let's say you have a set of K rays {Rk}, where Rk is k-th ray parameterized by the equation Rk(t)=Xk+t*Dk, where Xk is the point of origin and Dk is the direction of travel.

2) Given a M-by-M grid (G) sampled on a plane passing through point P with unit normal N, find points of intersection between G and all of the rays {Rk}.

3) Find transformation that aligns N with the z-axis, and horizonal and vertical grid lines of G with the x- and y-axes, respectively.

4) Apply the transformation from (3) to the points of intersection in (2).

5) Bin data from (4) to get the intensity image corresponding to the number of hits in a given bin (i.e., grid cell/pixel).

This process can be easily generalized to get a volumetric image using a 3D grid. This is because 3D grid is just a stack of 2D grids.

[Albert]

Thank you very much, Anton Semechko. You are right. Actually I did something similar to what you described in your step 2, but it was relatively slow, especially when the number of rays is very large. Hopefully, we are on the same page, but looks like we are not. I was trying to trace a curve through an image reconstruction volume, not a straight line. For a straight ray, it will be much easier to carry out. It is already done, please see here: https://www.mathworks.com/matlabcentral/fileexchange/28496-tomobox?focused=5159748&tab=function

[Anton Semechko]

Yeah, I missed the part where you said that the rays aren't straight. In that case, rays can be (i) approximated as piecewise linear curves composed of several straight line segments, (ii) points intersection can be approximated by finding line segments that intersect with the plane containing the grid, (iii) if the rays have analytical representation, points of intersection can be refined using a root-finding algorithm.