TPSF with MCX

98 views
Skip to first unread message

cate.ame...@gmail.com

unread,
Nov 9, 2020, 6:00:29 AM11/9/20
to mcx-users
Dear professor Fang,

I have used MMCLAB many times to simulate photon propagation inside head meshes, without any problem.
Now I am using MCXLAB, and I found a strange problem when I try to create the time point spread function (TPSF). I attach you the code I am using:

w690 = mcxdetweight(detp,ncfg.prop,ncfg.unitinmm);    
            t690 = mcxdettime(detp,ncfg.prop,ncfg.unitinmm);     

            nTG690 = round((cfg.tend-cfg.tstart)/cfg.tstep);
            tpsf690 = zeros(nTG690,1);

            ntg690 = ceil((t690-cfg.tstart)/cfg.tstep);

            for i=1:length(w690)
                
                tpsf690(ntg690(i)) = tpsf690(ntg690(i))+w690(i);
               
            end

The error I get, is that tha variable ntg690(i) sometimes has a number higher than nTG690! Meaning that the arrival time of the photons is higher than the ending time of the simulation. Is that possible? 
Could you please help me in understanding what is the error I am doing?

Thank you in advance,

Caterina Amendola

Fang, Qianqian

unread,
Nov 9, 2020, 1:45:17 PM11/9/20
to mcx-...@googlegroups.com, cate.ame...@gmail.com
On 11/9/20 6:00 AM, cate.ame...@gmail.com wrote:
Dear professor Fang,

I have used MMCLAB many times to simulate photon propagation inside head meshes, without any problem.
Now I am using MCXLAB, and I found a strange problem when I try to create the time point spread function (TPSF). I attach you the code I am using:



hi Caterina

we already have a function for this, can you try that?

https://github.com/fangq/mcx/blob/master/utils/mcxdettpsf.m


w690 = mcxdetweight(detp,ncfg.prop,ncfg.unitinmm);    
            t690 = mcxdettime(detp,ncfg.prop,ncfg.unitinmm);     

            nTG690 = round((cfg.tend-cfg.tstart)/cfg.tstep);
            tpsf690 = zeros(nTG690,1);

            ntg690 = ceil((t690-cfg.tstart)/cfg.tstep);

            for i=1:length(w690)
                
                tpsf690(ntg690(i)) = tpsf690(ntg690(i))+w690(i);
               
            end

The error I get, is that tha variable ntg690(i) sometimes has a number higher than nTG690! Meaning that the arrival time of the photons is higher than the ending time of the simulation. Is that possible? 
Could you please help me in understanding what is the error I am doing?

can you print out the value that are larger than nTG690? please let me know if they are slightly over the maximum time gate (within numerical precision round-off errors) or significantly longer?


Qianqian


Thank you in advance,

Caterina Amendola

--
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/3ba4e4f1-e66f-4bf6-89d2-ed3d2dc4e663n%40googlegroups.com.


Caterina Amendola

unread,
Nov 9, 2020, 2:13:09 PM11/9/20
to Fang, Qianqian, mcx-...@googlegroups.com
Dear Qianqian,
Thank you so much! With this function I have no problem anymore.
Just one question, my unit is not 1mm, but 0.86mm. Do I have to modify the code in mcxdettpsf  Inserting my cfg.unitinmm? So from:
replayweight=mcxdetweight(detp,prop,1);
replaytime=mcxdettime(detp,prop,1);
to this:
replayweight=mcxdetweightdetp,prop,cfg.unitinmm);
replaytime=mcxdettime(detp,prop,cfg.unitinmm);

More over, the code seems to me very similar to the one I attached. But in mcxdettpsf you put equal to the maximum time (nTG) all the photons that arrive after tend, why?

The values I obtained were not significantly higher than tend (a few ps after tend).

Thank you very much,

Caterina


Virus-free. www.avast.com

Fang, Qianqian

unread,
Nov 9, 2020, 2:20:31 PM11/9/20
to mcx-...@googlegroups.com, Caterina Amendola
On 11/9/20 2:12 PM, Caterina Amendola wrote:
Dear Qianqian,
Thank you so much! With this function I have no problem anymore.
Just one question, my unit is not 1mm, but 0.86mm. Do I have to modify the code in mcxdettpsf  Inserting my cfg.unitinmm? So from:


hi Caterina,

mcxdetweight and mcxdettime take care of the voxel size scaling internally, see

https://github.com/fangq/mcx/blob/master/utils/mcxdettime.m#L29
https://github.com/fangq/mcx/blob/master/utils/mcxdettime.m#L29

you don't need to manually set it using the 3rd input.


Qianqian


Reply all
Reply to author
Forward
0 new messages