Simulation problem

[Long Chen]

Dear Dr. Fang,

I've been using MCXLAB recently and encountered a strange problem. I did a simple two-dimensional simulation, and the organizational structure is shown in Figure 1 below. The background with a rectangle in the middle, using a Gaussian light source, with parameters like radius and position set according to Figure 2, and optical parameters set according to Figure 3. The simulation results are shown in Figure 4. However, when I kept other parameters fixed and only changed the optical parameters of the background as shown in Figures 5 and 6, the results obtained were shown in Figures 7 and 8. A strange phenomenon occurs where the simulated values increase as they get closer to the square. Normally, this value should gradually decrease. I haven't encountered this situation while using mcmatlab and mcxyz. Is this phenomenon normal, or am I setting the parameters incorrectly? And, what is the relationship between the values simulated by MCXLAB and the light flux we usually use? Could you please clarify this, professor? Thank you very much!

Best regards,
cl

Figure 1

Figure 2

Figure 3

Figure 4

Figure 5

Figure 6

Figure 7

Figure 8

[Qianqian Fang]

hi Long,

instead of posting the screenshot of your code, can you please directly post your script as text?

the screenshot not only takes ~200x larger space (197 kB vs 1 kB) to store and transmit, it also won't allow me to copy/paste to matlab to reproduce your issue. I tried using google lens to ORC the code from your screenshot, but it produced a lot of errors. So sharing your code in the form of image highly inefficient.

we have done thorough comparison with mcxyz using its skinvessel benchmark, and mcx matches the result excellently - here is the script - if there is any major difference, I am pretty sure it was because you set the simulation condition differently.

https://github.com/fangq/mcx/blob/master/mcxlab/examples/demo_mcxyz_skinvessel.m

--
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/114025f8-de7b-4a77-934c-b2b37c56bf23n%40googlegroups.com.

[Long Chen]

Dear Dr. Fang,

ok, thank you for your reply, the codes are as follows:

% ==========================================================================
% A sample 2D MCX simulation
%
% Author: Qianqian Fang <q.fang at neu.edu>
% ==========================================================================

% you must define a 3D array with one singleton dimension (with length 1).
% unfortunately if you define z as singleton, matlab will make the array 2D
% instead of 3D, so we have to purmute it to make the 1st dimension
% singleton.9+

% only clear cfg to avoid accidentally clearing other useful data
clear cfg;

cfg.nphoton = 1e6;

cfg.vol = permute(uint8(ones(400, 400)), [3, 1, 2]); % from 2d to 3d
cfg.vol(1, 100:300, 100:300) = 2;
% cfg.vol(1, 80:90, 10:20) = 0;
% cfg.vol(:,:,21:25) = 2;
% cfg.vol(:,:,26:35) = 3;
% cfg.vol(:,:,36:end) = 4;
cfg.issrcfrom0 = 1;
% cfg.srctype = 'pencil';
cfg.srctype = 'gaussian';
cfg.srcparam1 = [40 0 0 0];  
cfg.srcparam2 = [0 0 0 0];  


cfg.srcpos = [0,200, 0];   % src position must be located in the 2D plane
cfg.srcdir = [0 0 1];    % src dir must align in the plan (y-z in this case)

cfg.gpuid = 1;        
cfg.autopilot = 1;
cfg.gscatter = 1000;


cfg.prop = [0 0 1 1            % medium 0: the environment
    0.0004 0.0009 0.9 1.37   % water
    0.01  1.0 0.9 1.37      % absorber
    0.08 40.9 0.9 1.37];   %
cfg.tstart = 0;
cfg.tend = 1e-7;
cfg.tstep = 1e-7;
flux = mcxlab(cfg);

cw = squeeze(flux.data);

cw = cw / max(cw(:));          
figure;
imagesc(cw);
colormap('jet');
axis equal;

 
figure;
imagesc(log10(abs(cw)));
colormap('jet');
axis equal;

I originally thought the result was light fluence, but after thinking about it, this result feels like it is the energy accumulated over time that would lead to such an outcome. So what is the relationship between this result and the light fluence we commonly use? Can the accumulation over a shorter period of time (10^-9) be considered light fluence? Could you please help me with this question? Thank you!

Best regards,
cl

[Qianqian Fang]

hi Long,

the 2D mode in mcx is special - it is a theoretical 2D space photon transport, which does not exist in real world. I do not believe any other MC solvers has this type of output.

in the 2D mode, photon is limited in propagate in a 2D plane, without the ability to have a scattering angle escaping from the plane. This is entirely different from the 2D solution produced by MCML - which performs 3D photon transport, but accumulates solutions in a 2D radial plane due to the rotational symmetry of the space.

all other MC solvers that I know of perform 3D photon transport. They do not have the option to limit photon's propagation entirely inside a 2D space.

again, this only has theoretical value, and you shall not compare it with other 3D solver solutions or use it to analyze physical world measurements.

if you are interested in comparing with mcxyz, you should start with the skinvessel script that I linked previously.

Qianqian

To view this discussion visit https://groups.google.com/d/msgid/mcx-users/c8504fe2-eaf0-4b66-b8b1-c48815931846n%40googlegroups.com.

[Long Chen]

Dear Dr. Fang,
Thank you for your reply, I'm not trying to compare the quality of MCXLAB with other simulation software, but in the process of using it, I think MCXLAB has many similarities and even better results than mcxyz and others. I went to experiment with the code you gave, and the effect is indeed similar to MCXYZ, or even better, as shown in Figure 1 and 2, you can see that the value of the result is attenuated with the depth of light. But going back to my first question, I want to simulate a situation where the surroundings are water, which means that the optical parameters will be much lower than the values in the 'skinvessel' example, and if the mua and mus are set to 0.0001mm-1, there will be a case where the light decreases first with the increase in depth, and then increases suddenly when it is close to the object, which is the case in both 3D and 2D. Since I can't simulate this situation in mcxyz for the time being, I took MCmatalb (the matlab version of mcxyz) and MCXALB to do a three-dimensional comparison experiment on this situation, and took the xy section containing the light source to show, the organizational structure and parameters are shown in Figure 3 and 4, and the results are shown in Figure 5 and 6. After I reduced cfg.tstart, cfg.tend, cfg.tstep to 10^-9, it gradually became a case of decay with depth, which made me realize that the result might not be the luminous flux I thought it was, but probably the flunence rate that accumulated over time. Is that right? I may not fully understand the results of the MCXLAB simulation. Can I ask you about the relationship between the simulation results (flux.data) and the luminous flux, and how to convert them? Could you please help me with this question? Thank you!

Best regards,

cl

figure 1

figure 2

figure 3

figure 4

figure 5  ( MCmatlab )

figure 6  (MCXLAB)

[Long Chen]

Luminous flux similar to the formula below

[Qianqian Fang]

hi Long,

mcx output these data types (see this section of the README: https://github.com/fangq/mcx?tab=readme-ov-file#output-files)

1. fluence-rate (unit: 1/mm^2/s): photon area density across an infinitesimal sphere around any given spatial position per second

   this is the default output for mcx (cfg.outputtype='flux') - for historical reasons, it was called "flux" but it is actually fluence-rate

2. fluence (unit: 1/mm^2): photon area density across an infinitesimal sphere around any given spatial position

  fluence (cfg.outputtype='fluence') is basically time-integration of fluence-rate within any given time-gate (cfg.tstep), or the entire time window if there is only a single time gate (cfg.tend=cfg.tstep)

  if the single time-gate is long enough to capture all the changes in the TPSF, then this is the CW fluence, otherwise, it is a time-domain/time-resolved/time-gated fluence

3. per-voxel energy-deposition (unit: count, or J) - it is the multiplication of fluence Phi(r) (1/mm^2) with the local absorption coefficient (mua(r), 1/mm), integrated for every voxel, i.e. E[i,j,k] = Phi[i,j,k] * mua[i,j,k] * Vvox, where Vvox (mm^3) is the volume of a voxel.

  this is the accumulation of the energy-loss for all photons passing by any given voxel. if you sum all energy-deposition across all voxels, plus the diffuse-reflected/transmitted energy, it should be 1 (if normalization is on)

energy-deposition allows discontinuity due to mua(r), but fluence/fluence-rate are expected to be spatially differentiable.

let me know if this is what you are looking for

Qianqian

To view this discussion visit https://groups.google.com/d/msgid/mcx-users/97e2db79-1953-4bd5-8704-ae96ab8d1e54n%40googlegroups.com.

[Long Chen]

Dear Dr. Fang,
Thank you for your reply. Yes, what I'm looking for should be the second case, I'll try again, thank you very much!

Best wishes

cl