Having issues with photon timing and refractive index

[bcleary1120]

I've been using MCMatlab in order to generate data for time-domain fluorescence imaging by placing individual light sources at all points within the desired "tumor" object and summing/convolving each of their individual signals to get a set of impulse signals at each desired point. While doing this, I noticed that the time-of-flight from the sources to the detector was not aligning with the predicted TOF using c (~0.3 mm/ps) and the set refractive index (1.33). Instead, the time difference between light sources placed 1mm apart gives a value for c that is ~1 mm/ps. I am unsure how to resolve this issue since c is defined in the code as 29,979,245,800 cm/s. Additionally, I have noticed that changing the refractive index for the material (mediaProperties.n) does not affect this time difference at all. I've attached a normalized time-domain plot to show what I'm seeing. The first inclusion is a point source at 1mm below the surface and the second is an identical point source 2mm below the surface with their peak intensities marked in ps. I am using the July 2023 release of MCMatlab.

[ankrh]

Hi, thanks for the feedback!

One thing to be aware of is that if you have matchedInterfaces = true, then all refractive indices are set to 1, even for the purposes of time tagging. Some of your discrepancies might stem from that.

I also noticed a bug where the plotting function would not plot results from a light collector that is time-resolved but not spatially resolved, although it's unlikely to have anything to do with your issues.

I have made a new commit that fixes this bug and adds a warning when the user uses time tagging and has matchedInterfaces = true.

Try unpacking and running the attached test.m file on the new commit and see how the expected photon arrival time matches the measured arrival time.

test.zip

If your code still gives wrong results then I can take a look at it if you can share a minimal working example model file.

[bcleary1120]

For some reason, whenever I switch the matchedInterfaces flag to false, the data turns to zero at all points. I am entirely unsure as to why this is happening, but I have attached my code here.

TDFI_v4.zip

[ankrh]

Your problem is that only light that hits a boundary at a location where the refractive index is 1 is allowed to "escape". If the refractive index is not 1 then the photon is "killed" and cannot be counted on the detector. There is a warning for this but it will only be checked for when the user runs the plot() function, which your script doesn't :-). I should probably move these kinds of warnings to the end of the runMonteCarlo() function rather than the plot() function.

The reason that photons hitting a boundary where n != 1 are killed is that in order to know if the photon would hit the detector, we need to know the correct angle it's propagating at, which means the refraction out of the medium must be considered.

You can fix this by adding a thin layer of air on top of your model, just as I did in my test.m file.

[bcleary1120]

This fixed the code running, though I am still unsure as to why my timings seem off. Running this with just a one-way (source to detector) simulation, my difference in arrival times between light sources spaced 1mm apart is only 0.65ps, rather than ~4.4ps. Do you know why this might be? I could have easily made an error in my code, but I'm not sure where. I greatly appreciate all of your help on this.

[ankrh]

I think I found it.

In this line:
[x,y,z] = posFunc(i,j,k,inc,xMin,yMin,zMin);
you're incrementing x, y and z by the same increment even though model.G.dz is not the same as model.G.dx or model.G.dy.

A better way of getting x, y, z values of voxels is this:

for i = 1:nX
x = model.G.x(i);
for j = 1:nY
y = model.G.y(j);
for k = 1:nZ
z = model.G.z(k);
...

[bcleary1120]

Thank you so much! This fixed it completely.