fangq / mcx

Monte Carlo eXtreme (MCX) - GPU-accelerated photon transport simulator
http://mcx.space
Other
132 stars 73 forks source link

For MEDIA_ASGN_F2H where mua/mus are nan ensure that the provided g/n… #224

Closed lkeegan closed 2 months ago

lkeegan commented 4 months ago

… values are nontheless read

Check List

Before you submit your pull-request, please verify and check all below items

If your commits included in this PR contain changes that did not follow the above guidelines, you are strongly recommended to create a clean patch using git rebase and git cherry-pick to prevent in-compliant history from appearing in the upstream code.

Moreover, you are highly recommended to

Please copy/paste the corresponding Issue's URL after the below dash

fangq commented 4 months ago

@lkeegan, can you test my above fix in https://github.com/fangq/mcx/commit/80b57941cea33e10ab369723846108606d1cdba6 ?

this should allow mcx to read g/n regardless mua/mus=nan or not.

lkeegan commented 4 months ago

With 80b5794 I see very different behaviour in the MEDIA_ASGN_F2H case with boundary padding.

I think this is because

https://github.com/fangq/mcx/blob/80b57941cea33e10ab369723846108606d1cdba6/src/mcx_utils.c#L3526-L3532

https://github.com/fangq/mcx/blob/80b57941cea33e10ab369723846108606d1cdba6/src/mcx_utils.c#L3534-L3539

So I think the g/n values are now being read, but the boundary mua/mus values are now nan instead of 0

fangq commented 3 months ago

sorry for taking a while to get back to this - let me know if the above commit fixes the issue.

lkeegan commented 3 months ago

Thanks! With this commit the behaviour is the same as #224 - the ref total ref=-inf issue is fixed, but I still don't get any reflectance with b=1:

https://gist.github.com/lkeegan/ee77309c2b48883278fd9716c95d1c34

And there are also many error messages of the form:

ERROR: should never happen! mediaid=0 idx1d=88F isreflect=1 gcfg->doreflect=1 n1=1.299805 n2=1.299805 isdet=0 flipdir[3]=2 p=(31.562723 18.700874 1.000000)[31 18 0]

(see https://github.com/fangq/mcx/issues/225, https://groups.google.com/g/mcx-users/c/bToGluYYdao/m/5XGk9jeHAAAJ)

fangq commented 3 months ago

ok, let me do more test on this.

kdreher commented 3 months ago

I was wondering if there are any updates on this.

fangq commented 2 months ago

@lkeegan and @kdreher, sorry for my late response.

I could not find the files vol_bg.bin and vol_bg.json in your jupyter notebook, so I wrote a simple test case based on a built-in demo

clear cfg;
cfg.nphoton = 1e7;
cfg.vol = uint8(ones(60, 60, 60));
cfg.srcpos = [30 30 1];
cfg.srcdir = [0 0 1];
cfg.gpuid = 1;
% cfg.gpuid='11'; % use two GPUs together
cfg.autopilot = 1;
cfg.prop = [0 0 1 1; 0.005 1 0 1.37];
cfg.tstart = 0;
cfg.tend = 5e-9;
cfg.tstep = 5e-9;

mua1 = 0.003;
mua2 = 0.1;
mua = single(reshape(repmat([mua1:(mua2 - mua1) / (60 - 1):mua2]', 60, 60), 60, 60, 60));

mus1 = 1.0;
mus2 = 5.0;
mus = single(reshape(repmat([mus1:(mus2 - mus1) / (60 - 1):mus2]', 60, 60), 60, 60, 60));

cfg.vol = reshape([mua(:)'; mus(:)'], [2 60 60 60]);
cfg.vol(3, :, :, :) = 0;
cfg.vol(4, :, :, :) = 1.37;

cfg.vol(1,:,:,1)=nan;
cfg.issaveref=1;
cfg.isreflect=1;

flux = mcxlab(cfg);
mcxplotvol(log10(flux.dref))

with the above script, I was able to get non-empty dref, see screenshot below - can you reproduce this in your environments (mcx binary+python)?

image

fangq commented 2 months ago

closing for now as the issue can no longer be reproduced.

feel free to reopen if the problem persists.

lkeegan commented 2 months ago

Thanks, I'll try to reproduce your example script with python+mcx.

To reproduce the original reported issue, here is a zipfile with vol_bg.bin and vol_bg.json:

vol_bin_and_json.zip

And the command line args for mcx:

mcx -f vol_bg.json -O F -a 1 -F jnii -Z 2 -b 1 --saveref 1 -s b0_bg

Running the above results in no reflectance data and many error messages of the form:

ERROR: should never happen! mediaid=0 idx1d=1981 isreflect=1 gcfg->doreflect=1 n1=1.299805 n2=1.299805 isdet=0 flipdir[3]=2 p=(49.131416 54.812988 1.000000)[49 54 0]

If I set -b 0 instead then it runs without error messages and produces the expected reflectance data.

Many thanks for your help!