bug traj output pmcx.mcxlab(), traj.id 0 for all data

65 views
Skip to first unread message

G Buist

unread,
Nov 23, 2023, 9:54:54 AM11/23/23
to mcx-users
Dear professor Fang,

When testing the output of the pmxc.mcxlab(cfg) (version 0.2.6) function for the trajectory output mode ( cfg['debuglevel'] = 'M') the id array contained only zeros, code below.  When porting this function from matlab to pmcx I suspect a bug was introduced.
I also tried to construct an id array from the 'data' array in both the pmcx.run() and pmcx.mcxlab() cases which only worked for the pmcx.run() method.

Kind Regards,
Gijs

###### code ########

#%% []
import numpy as np
import pmcx
import matplotlib.pyplot as plt
from struct import unpack

#%%[]

print(pmcx.__version__)

cfg={}
cfg['nphoton']=1e1
cfg['vol']=np.ones([60, 60, 60], dtype='uint8')
cfg['vol'][20:40, 30:40, 20:30]=2
cfg['tstart']=0
cfg['tend']=5e-9
cfg['tstep']=5e-9
cfg['srcpos']=[30,30,0]
cfg['srcdir']=[0,0,1]
cfg['prop']=[[0, 0, 1, 1], [0.005, 1, 0.01, 1.37], [0.1, 10, 0.9, 1]]
cfg['detpos']=[[30,20,0,1]]        # to detect photons, one must first define detectors
cfg['debuglevel'] = 'M'

res = pmcx.mcxlab(cfg)

print(res.keys())
print(res['traj'].keys())
print(np.all(res['traj']['id'] == 0))


### output is c struct that needs to be converted to integers in the photon id case (other cases are floats which numpy automatically assumes I guess)
structAr = np.asarray(res['traj']['data'][0], order = 'C')
lenstr = str(len(structAr))
traj_id =  np.array( unpack( lenstr + 'I', structAr))
print(np.all(traj_id == 0))

plt.figure()
plt.plot(res['traj']['id'])
plt.show()


res2 = pmcx.run(cfg)

print(res.keys())
print(res2['traj'].shape)


### output is c struct that needs to be converted to integers in the photon id case (other cases are floats which numpy automatically assumes I guess)
structAr = np.asarray(res2['traj'][0], order = 'C')
lenstr = str(len(structAr))
traj_id =  np.array( unpack( lenstr + 'I', structAr))
print(np.all(traj_id == 0))

plt.figure()
plt.plot(traj_id)
plt.show()

######### Output#####################

0.2.6
nphoton: 10
tstart: 0
tstep: 5e-09
tend: 5e-09
srcpos: [30, 30, 0, 1]
srcdir: [0, 0, 1, 0]
dict_keys(['traj', 'flux', 'stat'])
dict_keys(['pos', 'id', 'data'])
True
True
nphoton: 10
tstart: 0
tstep: 5e-09
tend: 5e-09
srcpos: [30, 30, 0, 1]
srcdir: [0, 0, 1, 0]
dict_keys(['traj', 'flux', 'stat'])
(6, 1720)
False

###############################################################################
#                      Monte Carlo eXtreme (MCX) -- CUDA                      #
#          Copyright (c) 2009-2023 Qianqian Fang <q.fang at neu.edu>          #
#                https://mcx.space/  &  https://neurojson.org/                #
#                                                                             #
# Computational Optics & Translational Imaging (COTI) Lab- http://fanglab.org #
#   Department of Bioengineering, Northeastern University, Boston, MA, USA    #
###############################################################################
#    The MCX Project is funded by the NIH/NIGMS under grant R01-GM114365      #
###############################################################################
#  Open-source codes and reusable scientific data are essential for research, #
# MCX proudly developed human-readable JSON-based data formats for easy reuse,#
#  Please consider using JSON (https://neurojson.org/) for your research data #
###############################################################################
$Rev::996580$ v2023  $Date::2023-09-23 00:39:39 -04$ by $Author::Qianqian Fang$
###############################################################################
- code name: [Fermi MCX] compiled by nvcc [10.2] for CUDA-arch [350] on [Sep 23 2023]
- compiled with: RNG [xorshift128+] with Seed Length [4]

GPU=1 (NVIDIA GeForce RTX 3080 Ti) threadph=0 extra=10 np=10 nthread=327680 maxgate=1 repetition=1
initializing streams ... init complete : 1 ms
requesting 1536 bytes of shared memory
launching MCX simulation for time window [0.00e+00ns 5.00e+00ns] ...
simulation run# 1 ...
kernel complete:   101 ms
retrieving fields ... saved 1720 trajectory positions, total: 1720 detected 0 photons, total: 0 transfer complete: 109 ms
normalizing raw data ... source 1, normalization factor alpha=20000000.000000
data normalization complete : 113 ms
simulated 10 photons (10) with 327680 threads (repeat x1)
MCX simulation speed: 2.00 photon/ms
total simulated energy: 10.00 absorbed: 39.98276%
(loss due to initial specular reflection is excluded in the total)
###############################################################################
#                      Monte Carlo eXtreme (MCX) -- CUDA                      #
#          Copyright (c) 2009-2023 Qianqian Fang <q.fang at neu.edu>          #
#                https://mcx.space/  &  https://neurojson.org/                #
#                                                                             #
# Computational Optics & Translational Imaging (COTI) Lab- http://fanglab.org #
#   Department of Bioengineering, Northeastern University, Boston, MA, USA    #
###############################################################################
#    The MCX Project is funded by the NIH/NIGMS under grant R01-GM114365      #
###############################################################################
#  Open-source codes and reusable scientific data are essential for research, #
# MCX proudly developed human-readable JSON-based data formats for easy reuse,#
#  Please consider using JSON (https://neurojson.org/) for your research data #
###############################################################################
$Rev::996580$ v2023  $Date::2023-09-23 00:39:39 -04$ by $Author::Qianqian Fang$
###############################################################################
- code name: [Fermi MCX] compiled by nvcc [10.2] for CUDA-arch [350] on [Sep 23 2023]
- compiled with: RNG [xorshift128+] with Seed Length [4]

GPU=1 (NVIDIA GeForce RTX 3080 Ti) threadph=0 extra=10 np=10 nthread=327680 maxgate=1 repetition=1
initializing streams ... init complete : 1 ms
requesting 1536 bytes of shared memory
launching MCX simulation for time window [0.00e+00ns 5.00e+00ns] ...
simulation run# 1 ...
kernel complete:   101 ms
retrieving fields ... saved 1720 trajectory positions, total: 1720 detected 0 photons, total: 0 transfer complete: 110 ms
normalizing raw data ... source 1, normalization factor alpha=20000000.000000
data normalization complete : 114 ms
simulated 10 photons (10) with 327680 threads (repeat x1)
MCX simulation speed: 2.00 photon/ms
total simulated energy: 10.00 absorbed: 39.98276%
(loss due to initial specular reflection is excluded in the total)



Qianqian Fang

unread,
Nov 28, 2023, 3:21:11 AM11/28/23
to mcx-...@googlegroups.com

thanks Gijs for reporting this, and the debug result for why this happens.

you were correct that pmcx.mcxlab() fails to "typecast" the traj.id buffer from float32 to uint32, like we did in mcxlab (https://github.com/fangq/mcx/blob/v2023/mcxlab/mcxlab.m#L549C1-L550C1)

I just created a new issue on github (https://github.com/fangq/mcx/issues/199) and was able to fix this with the below patch

https://github.com/fangq/mcx/commit/fae5b7273ca535385be3020ec4ec32e084dbdb39

the sample code seems to work well after the patch. The pmcx version has been bumped to 0.2.7 and is currently being compiled - once completed, should appear on pypi shortly, please do an upgrade and let me know if you see any additional issues.

--
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 on the web visit https://groups.google.com/d/msgid/mcx-users/c1f3f613-96a2-4f63-be17-cf4a214184f5n%40googlegroups.com.

Buist, G. (Gijs)

unread,
Nov 28, 2023, 3:37:45 PM11/28/23
to mcx-...@googlegroups.com

 

I upgraded pmcx to 0.2.7 and the sample code now seems to work and seems to produce correct photon id’s.

 

Thank you for the fix.

 

Kind regards,
Gijs Buist

--
You received this message because you are subscribed to a topic in the Google Groups "mcx-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/mcx-users/O79Yk4G8_xM/unsubscribe.
To unsubscribe from this group and all its topics, send an email to mcx-users+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/mcx-users/a18156f4-750f-4deb-ba7b-1d97de5025b0%40northeastern.edu.

Reply all
Reply to author
Forward
0 new messages