Skip to content

divide by actual mus in otRFmus/otWP/otWPTOF in SVMC replay, fix chain rule for mus/musp Jacobians in adjoint mode#263

Open
HirviP wants to merge 1 commit into
fangq:masterfrom
HirviP:master
Open

divide by actual mus in otRFmus/otWP/otWPTOF in SVMC replay, fix chain rule for mus/musp Jacobians in adjoint mode#263
HirviP wants to merge 1 commit into
fangq:masterfrom
HirviP:master

Conversation

@HirviP
Copy link
Copy Markdown
Contributor

@HirviP HirviP commented Apr 9, 2026

1. Update the code to divide "b"/"wptof" (otWPTOF) counts also with the current scattering coefficient in SVMC mode during ray-tracing in replay. Set weight to zero if mus=0, since scattering events shouldn't happen in theory (and I got one very high value on boundary, which might result from this). Fix the scattering coefficient to be the true physical scattering coefficient for the division, if cfg.unitinmm is not 1 and mus has been scaled previously (cfg->prop[i].mus *= cfg->unitinmm).

NOTE: Currently source code compilation fails since cfg->unitinmm is undefined on line 2255 in mcx_core.cu, where mus-division is performed. I did not find where the unitinmm is stored – gcfg does not have the field – can you help me? Otherwise compilation is successful.

The below figure shows that the voxel-based (A correct, B modified) and SVMC (C; 0.5 shift not fixed yet) "scattering Jacobians" only match if we multiply the voxel-wise scattering coefficients with cfg.unitinmm=0.5 (mus/2, B) before scaling (post-processing in MCX; same optical parameters used in both simulations):

Debug_SVMC_Replay_WithShift_Det3_PH


2. Fix the inner derivative sign when converting derivatives with respect to the diffusion coefficient (D) to derivatives with respect to the scattering coefficient (mus or musp) using the chain rule in the adjoint mode, since:

dA/dmus = dA/dD * dD/dmus = J_D * ( - 1/(3*(1-g)*mus^2))

In my simulations with demo_mcx_adjoint_jacobian.m, this fix matched the signs of the adjoint and replay scattering Jacobians (dXdmus and dYdmus; should match up to a constant scaling factor depending on measurement type), as shown in the figure below (modified source code used):

Adjoint_dXYdmus_sign_PH

Check List

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

  • You have only modified the minimum number of lines that are necessary for the update; you should never add/remove white spaces to other lines that are irrelevant to the desired change.
  • You have run make pretty (requires astyle in the command line) under the src/ folder and formatted your C/C++/CUDA source codes before every commit; similarly, you should run python3 -m black *.py (pip install black first) to reformat all modified Python codes, or run mh_style --fix . (pip install miss-hit first) at the top-folder to format all MATLAB scripts.
  • Add sufficient in-code comments following the doxygen C format
  • In addition to source code changes, you should also update the documentation (README.md, mcx_utils.c and/or mcxlab.m) if any command line flag was added or changed.

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

  • Add a test in the mcx/test/testmcx.sh script, following existing examples, to test the newly added feature; or add a MATLAB script under mcxlab/examples to gives examples of the desired outputs
  • MCX's simulation speed is currently limited by the number of GPU registers. In your change, please consider minimizing the use of registers by reusing existing ones. CUDA compilers may not be able to optimize register counts, thus require manual optimization.
  • Please create a github Issue first with detailed descriptions of the problem, testing script and proposed fixes, and link the issue with this pull-request

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

…h the current scattering coefficient in SVMC mode. Fix the scattering coefficient to be the true physical scattering coefficient, if cfg.unitinmm is not 1 and mus has been scaled (cfg->prop[i].mus *= cfg->unitinmm). NOTE: Currently compilation fails since cfg->unitinmm is undefined on line 2255! Where is it stored – gcfg does not have the field? Otherwise compilation is successful. In addition, fix the inner derivative sign when converting derivatives with respect to the diffusion coefficient (D) to derivatives with respect to the scattering coefficient (mus or musp) using the chain rule in the adjoint mode. In my simulations, this fix matches the signs of the adjoint and replay scattering Jacobians.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant