Widefield src and detector inconsistencies

[alejandro martinez]

Hello again Dr. Qianqian, 

I'm working with widefield source and detector in a reflectance configuration and have encountered several inconsistencies when running the experiments.

Depending on the order of the calls to mccaddsrc and mmcadddet the mesh completely changes, even sometimes removing the -1 elem of the light source. This causes a variety of bugs from random hangs, to incorrect detected photosn.

Below I present four simple experiments on a 2 layer cube using a widefield source and detector to explain the issues I've encountered. 

Common matlab code to generate the cubes:

voxelCube = ones(20, 20, 20);

voxelCube = uint8(voxelCube);

voxelCube(:,:,1:10) = 2;

[node, elem, face] = v2m(voxelCube, [], 1, 1, 'cgalmesh');

• CASE 1: Source further than detector.

1. • CASE 1.1: First call adddet and then addsrc:

srcPos = [5, 5, 35]; % FURTHER

detPos = [10,10,25]; % CLOSER

srcStruct = struct( ...

'srctype', 'planar', ...

'srcpos', srcPos, ...

'srcdir', [ 0,0,-1], ...

'srcparam1', [10 0 0 0], ...

'srcparam2', [0 10 0 0]);

detStruct = struct( ...

'srctype', 'disk', ...

'srcpos', detPos, ...

'srcdir', [ 0,0,-1], ...

'srcparam1', [20 0 0 0], ...

'srcparam2', [0 0 0 0]);

[nodes, elem] = mmcadddet(nodes, elem, detStruct);

[nodes, elem] = mmcaddsrc(nodes, elem, srcStruct);

elem(:, 1:4) = meshreorient(nodes, elem(:, 1:4));

eid = tsearchn(nodes, elem(:,1:4), srcPos);

savemmcmesh(fileNameBase, nodes, elem);

plotmesh(nodes, elem, 'y>8');

un = unique(elem(:,5));

We can correctly see cube (yellows), remesh scaffolding (light blue), source (blue) and detector (dark blue).

1. • CASE 1.2: First call addsrc and then adddet: (Just changed the order in the following lines)

[nodes, elem] = mmcaddsrc(nodes, elem, srcStruct);

[nodes, elem] = mmcadddet(nodes, elem, detStruct);

Detector is moved to the side.

• CASE 2: Detector further than source

1. • CASE 2.1: First call adddet and then addsrc:

srcPos = [5, 5, 25]; % CLOSER

detPos = [10,10,35]; % FURTHER

In this case, the unique elements show [-2;0;1;2], light source (-1) is gone, all the other show [-2;-1;0;1;2],

1. • CASE 2.2: First call addsrc and then adddet:

[nodes, elem] = mmcaddsrc(nodes, elem, srcStruct);

[nodes, elem] = mmcadddet(nodes, elem, detStruct);

This causes a range of inconsistencies also when running the simulations. As I detail below:

• CASE 1.1: No detected photons. Fluence looks ok.

retrieving flux ... saved 4232477 trajectory positions, total: 4232477 detected 0 photons, total: 0 transfer complete: 1909 ms saving data to file ... saving data complete : 1914 ms simulated 100000 photons (100000) with 1 devices (ray-tet 18890768) MCX simulation speed: 187.97 photon/ms total simulated energy: 100000.00 absorbed: 70.26094%

• CASE 1.2: Simulation just hangs forever, when activating -M P and cancelling progress stays at 99%, had to reduce number of photons to 10000. Still hang is not deterministic, sometimes it hangs forever and sometimes it finishes successfully, seems like the number of photons increases probability of hang.(maybe as there are more photons they reach some bad element or something that causes the hang, just speculation). This is the output from one of the successful runs.

retrieving flux ... saved 417725 trajectory positions, total: 417725 detected 522 photons, total: 522 transfer complete: 1729 ms saving data to file ... saving data complete : 1733 ms simulated 10000 photons (10000) with 1 devices (ray-tet 1556636) MCX simulation speed: 23.87 photon/ms total simulated energy: 10000.00 absorbed: 69.90376%

Even though it seems correct, when plotting the photons detected, they lie in that triangle faces that we saw in the figure of the case 1.2, not in the detector:

• CASE 2.1: No detected photons, no absorbed energy:

retrieving flux ... saved 200000 trajectory positions, total: 200000 detected 0 photons, total: 0 transfer complete: 1770 ms saving data to file ... saving data complete : 1772 ms simulated 100000 photons (100000) with 1 devices (ray-tet 400000) MCX simulation speed: 235.85 photon/ms total simulated energy: 100000.00 absorbed: 0.00000%

• CASE 2.2: Hangs again, same as case 1.2, when reducing photons to 10000 the simulation randomly hangs or finishes.

retrieving flux ... saved 422819 trajectory positions, total: 422819 detected 4017 photons, total: 4017 transfer complete: 1804 ms saving data to file ... saving data complete : 1809 ms simulated 10000 photons (10000) with 1 devices (ray-tet 1618513) MCX simulation speed: 22.83 photon/ms total simulated energy: 10000.00 absorbed: 70.44712%

However, in this case, when it finishes, the results look correct:

Simulations are run from cli with: 

mmc -f config.json -V 1 -k 0 -b 1 -M P -D SEM -d 1

And the following config.json is used (changing position of source and detector accordingly with the cases):

{ "Mesh": { "MeshID": "sample_cube", "InitElem": -1 }, "Session": { "Photons": 100000, "Seed": 42, "ID": "sample_cube", "DoPartialPath": true, "DoSaveExit": true, "DoDCS": false }, "Forward": { "T0": 0.0, "T1": 5e-09, "Dt": 5e-09 }, "Optode": { "Source": { "Type": "planar", "Pos": [5, 5, 35], "Dir": [0, 0, -1], "Param1": [10.0, 0.0, 0.0, 0.0], "Param2": [0.0, 10.0, 0.0, 0.0] }, "Detector": [ { "Pos": [10, 10, 25], "R": 20 } ] } }

Specifications:

iso2mesh toolbox version: 1.9.5-Rev

2 computers tested:

• ubuntu 24 x64 RTX4060 GPU (8GB VRAM)
• ubuntu 24 x64 QuadroM1200 (4GB VRAM)

MMC compiled from commit 31d6a06.

I would appreciate if you could take a look at this, since I can't seem to be able to design widefield

experiments at the moment.

Thank you very much,

please let me know if I can help with more information.

Alejandro

[alejandro martinez]

Hi, 

Hangs seem to only happen with the -c opencl/cuda, when i try cpu with -c sse there are no hangs.

Thanks

[Qianqian Fang]

hi Alejandro

sorry for the late reply.

this is a known limitation of the current implementation of the "automated mesh-tessellation workflow" (based on convhull + meshrefine) for widefield MMC.

mmcaddsrc and mmcadddet are designed to handle only non-spatially overlapping (in relation to the convex hull of the extended mesh) src and detectors. the added source and detector can not have a direct line-of-sight to each other - for example, a planar src on the top side of a cube with the widefield detector at the bottom. However, if the mesh tessellation of one of them touches the convex hull of the previously extended mesh with a src or detector, then, the previously tagged elements (0 for interconnecting elements, -1 for src-containing elements, and -2 for det-containing elements) are removed. 

this does not mean that spatially overlapping widefield src/det can not work with mesh re-tessellation. it just mean mmcaddsrc/mmcadddet can't do this automatically. You can, however, manually label the respective elements - those extended over the original mesh based on the desired label I mentioned above. mmc simulation will still respect these labels and perform launch/detection properly.

in other words, this is a mesh generator's limitation, not a problem for mmc per se.

on a related note, my student had recently published RT-MMC in this JBO paper

https://www.spiedigitallibrary.org/journals/journal-of-biomedical-optics/volume-31/issue-3/035002/Accelerating-mesh-based-Monte-Carlo-simulations-using-contemporary-graphics-ray/10.1117/1.JBO.31.3.035002.full

compared to tetrahedral-mesh based mmc, rt-mmc uses triangular mesh and ray-tracing hardware for photon migration. it has the capability of launching photons from anywhere in the domain without needing mesh retessellation. see Section 4.4.

the RT-MMC code has already been merged to the main mmc repo. to compile it, you need to install CUDA and a small NVIDIA OptiX ray-tracing library (header-only), and run "make optix" or "make optixmex". 

the OptiX path does not have all the features as the SSE/OpenCL/CUDA paths that were already in mmc for years. While RT-MMC supports arbitrary source launch, it currently does not have the logic for widefield detection. In this summer, after releasing the v2026 of the mcx/mmc suite, I will try to port some of the key features from OpenCL/CUDA implementations to OptiX kernels.

let me know if you want to give it a try.

Qianqian

--
You received this message because you are subscribed to the Google Groups "mcx-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mcx-users+...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/mcx-users/026e065a-7c85-4bb5-8790-1f609577e627n%40googlegroups.com.

[alejandro martinez]

Hi Qianqian, 

Thank you for all the valuable information! Now I fully understand the limitation of the automatic remeshing. I will try to manually label the mesh to place the detector and light source correctly. 

Regarding this, when I remesh to include an external light source the elements that are label as -1 are just the elem that contains the full source. But when I remesh to add only an external detector, not only the detector plane elements are -2 but also some of the scaffolding elements in between (see in the following image). Is this the intended behavior? From the widefield paper (Generalized mesh-based Monte Carlo...) Fig 1 (e,f) I can kind of see the same behavior, but not sure. If so, how to do that manually? Is it even needed or can I just set as -2 the elements containing the detector plane?

This new path with RT-MMC looks promising. For my experiments I need to capture photon output (pos, dir, ppath..) because I'm using a virtual camera in an arbitrary position in space to check the captured information.

would I be able to output the photons at the mesh surface (or any other way to capture them)? If so, please instruct me on how to do so  as I didn't find any docs yet.

Happy to test the new features!

Thanks again for your time, 

Alejandro