Detected Photons Flux

[Sravan Kumar Bhagavatula]

Hello,

I was wondering if it was possible to create a flux volume using only detected photons.

[Qianqian Fang]

hi Sravan

you asked this question on the right day!

I just finish implementing a new feature in MMC (not
in MCX yet) called "replay". What this feature does is
to store the detected photon seeds and rerun them
to produce the sensitivity/flux/fluence maps using
only the detected photons. This project is in
collaboration with Mark Niedre at Northeastern Univ.
The seeding approach was based on a paper from
Angelo Sassaroli (Sassaroli, Optics Letters 2011).

To use this feature, you need to first run the full
simulation using the baseline optical properties and
save the detected photon partial path and seeds.
Then you rerun mmc with "-E output.mch" where
output.mch is the detected photon file from the
first simulation. This will "replay" the detected photons
(you can use -P det_id to replay only a single detector)
and output either flux (-O F), or Jacobian (-O J/-O T)
in a second simulation.

My modified code is already in the git repository.
You can browse it online at

http://sourceforge.net/p/mcx/git/commit_browser
and
http://sourceforge.net/p/mcx/git/ci/master/tree/

You can download it with the following command.

git clone git://git.code.sf.net/p/mcx/git mcx-git

I am planning to port this feature to MCX, but given
the complexity, I expect this can take a few months
at least.

let me know if you want to give MMC a try.


Qianqian

> --
> 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 post to this group, send email to mcx-...@googlegroups.com.
> Visit this group at http://groups.google.com/group/mcx-users.
> For more options, visit https://groups.google.com/groups/opt_out.
>
>



The information in this e-mail is intended only for the person to whom it is
addressed. If you believe this e-mail was sent to you in error and the e-mail
contains patient information, please contact the Partners Compliance HelpLine at
http://www.partners.org/complianceline . If the e-mail was sent to you in error
but does not contain patient information, please contact the sender and properly
dispose of the e-mail.

[Kiefer Forseth]

Hi Sravan and Qianqian,

I actually had exactly the same question/request, this sounds great!

Best,

Kiefer

[Sravan Kumar Bhagavatula]

Thanks! I will have a look at MMC and work with it.

On Wednesday, July 17, 2013 11:57:15 AM UTC-5, Qianqian Fang wrote:

[Kiefer Forseth]

Hi Qianqian,

I've been working with MMC using the following two commands (from the onecube example)

./mmc -f onecube.inp -n 10000 -s data -d 1

./mmc -f onecube.inp -E data -O F -s detdata

I think the first one generates all of the simulations and the latter replays only those that terminated at detectors.

I know that the -P flag lets me determine a specific detector. My question is: do I use it on the first or second call? Ideally, it would be the second. This would allow me to run a single full simulation and then extract the subset of simulations ending at each detector efficiently.

Thanks!

On Wednesday, July 17, 2013 11:57:15 AM UTC-5, Qianqian Fang wrote:

[Qianqian Fang]

On 07/25/2013 07:42 AM, Kiefer Forseth wrote:

Hi Qianqian,

I've been working with MMC using the following two commands (from the onecube example)

./mmc -f onecube.inp -n 10000 -s data -d 1

./mmc -f onecube.inp -E data -O F -s detdata

I think the first one generates all of the simulations and the latter replays only those that terminated at detectors.

I know that the -P flag lets me determine a specific detector. My question is: do I use it on the first or second call? Ideally, it would be the second. This would allow me to run a single full simulation and then extract the subset of simulations ending at each detector efficiently.

here is my example showing how to use this feature, I use
the example/validation here, as the purpose of onecube
is for testing debug info, the detector location are off
the mesh actually :)

in the first step, you run this in the validation folder:

 ../../src/bin/mmc -f cube.inp -s cube -n 1000000 -b 0 -D TP -d 1 -q 1

the -q 1 enables saving seeds for each detected photon

you will now see a cube.mch file in the folder.
You can rename it to a different file to avoid overwriting.

     mv cube.mch cube_bk.mch

in the second step, you run

 ../../src/bin/mmc -f cube.inp -s cube -n 1000000 -b 0 -D TP -d 1 -E cube_bk.mch -O F -U 0 -P 1

this will load the recorded detected photons and "replay"
all of those captured by the 1st detector (-P 1) and output
the "raw" fluence (-O F) generated by these photons alone.
As you see, I disabled normalization (-U 0) because the
normalization will not be accurate in the replay mode.
The "-n 1000000" input are automatically ignored as the
photon number is determined by the .mch file in the
replay mode.

if you use "-O J", MMC will output the sensitivity map, i.e.
a "banana" pattern between the source and the selected
detector (-P 1). The output is in the cube.dat file.

 ../../src/bin/mmc -f cube.inp -s cube -n 1000000 -b 0 -D TP -d 1 -E cube_bk.mch -O J -P 1

You can replay photons from other detectors with more mmc calls
and create the full Jacobian. The "-O T" also produces
something similar to the Jacobian, but an approximated version
after the Taylor expansion.

Keep in mind that the Jacobian is defined on each node, and
is inclusive of the volumetric and time-gate effects (I am still debating
about the time gate). This is different from the default fluence/flux
output where these effects are "normalized" with "-U 1".
The only normalization in the -O J/T outputs is the number
of photons.

just as a sanity check, you can run

 OMP_NUM_THREADS=1 ../../src/bin/mmc -f cube.inp -s cube -n 1000000 -b 0 -D TP -d 1 -q 1 -E cube_bk.mch

as the second step, this will dump a cube.mch that is
bit-to-bit identical to cube_bk.mch except the header.

as I said, this is an experimental feature, let me know if
you have any feedback.

Qianqian

[Kiefer Forseth]

Thank you so much for your thorough reply! I'll work with it for awhile and let you know how it goes.

Best,

Kiefer

[Arun Nemani]

I'm running the following replay demo file and still getting an error for no seeds data. 

MMCLAB ERROR 999 in unit mmclab.cpp:418

Error: the 'seed' field can not be empty

cfg.nphoton=1e7;

cfg.seed=27182818;

[cfg.node,cfg.elem]=genT5mesh(0:2:60,0:2:60,0:2:60);

cfg.elemprop=ones(size(cfg.elem,1),1);

cfg.srcpos=[30 30 0];

cfg.srcdir=[0 0 1];

cfg.tstart=0;

cfg.tend=5e-9;

cfg.tstep=1e-10;

cfg.prop=[0 0 1 1;0.005 1.0 0.01 1.0];

cfg.debuglevel='TP';

cfg.isreflect=0;

cfg.outputtype='J';

cfg.detpos=[30. 20. 0. 1.

   30. 40. 0. 1.

   20. 30. 0. 1.

   40. 30. 0. 1.];

%cfg.issaveexit=1;

newcfg=mmclab(cfg,'prep');

[cube detp ncfg seeds]=mmclab(newcfg);

% now replay the detected photons

newcfg.seed=seeds.data;

newcfg.isnormalized=0;

[cube2 detp2 ncfg2 seeds2]=mmclab(newcfg);

% the two detected photon arrays should be the same. however, because

% the program uses multi-threading, the orders may be different

[isreplayed, mapidx]=ismember(detp.data',detp2.data','rows');

if(all(isreplayed))

   disp('replay is successful :-)');

else

   disp('replay failed :-(');

end

qmeshcut(cfg.elem(:,1:4),cfg.node,sum(cube2.data,2),'y=30.2');

view([0 1 0]);

Your help is greatly appreciated!

-Arun

[Qianqian Fang]

On 04/23/2014 05:37 PM, Arun Nemani wrote:

I'm running the following replay demo file and still getting an error for no seeds data. 

MMCLAB ERROR 999 in unit mmclab.cpp:418

Error: the 'seed' field can not be empty

cfg.nphoton=1e7;

cfg.seed=27182818;

[cfg.node,cfg.elem]=genT5mesh(0:2:60,0:2:60,0:2:60);

cfg.elemprop=ones(size(cfg.elem,1),1);

cfg.srcpos=[30 30 0];

cfg.srcdir=[0 0 1];

cfg.tstart=0;

cfg.tend=5e-9;

cfg.tstep=1e-10;

cfg.prop=[0 0 1 1;0.005 1.0 0.01 1.0];

cfg.debuglevel='TP';

cfg.isreflect=0;

cfg.outputtype='J';

cfg.detpos=[30. 20. 0. 1.

   30. 40. 0. 1.

   20. 30. 0. 1.

   40. 30. 0. 1.];

hi Arun

as I responded in the private email, your detector locations
are all within the mesh. because the photon detection calculations
are only applied to the photons escaping the mesh surface, your
detectors are too far away from these escape sites, thus, no
photons can be captured. as a result, your seeds.data below
is empty, and mmclab can not process an empty seed input.

please position your detectors on the mesh surface and try again.

Qianqian

You received this message because you are subscribed to the Google Groups "mmc-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mmc-users+...@googlegroups.com.
To post to this group, send email to mmc-...@googlegroups.com.
Visit this group at http://groups.google.com/group/mmc-users.
For more options, visit https://groups.google.com/d/optout.

[Klaas Rackebrandt]

Hi,

 

I folllowed your instructions and the simulation ran well.

 

I'd like to plot the results from the mch file to geht the banana shaped curved path oh the photons in Matlab.

 

Is there a way to do so?

 

with kind regards,

 

Klaas

[Qianqian Fang]

On 05/12/2014 08:16 AM, Klaas Rackebrandt wrote:

Hi,

 

I folllowed your instructions and the simulation ran well.

 

I'd like to plot the results from the mch file to geht the banana shaped curved path oh the photons in Matlab.

 

Is there a way to do so?

hi Klaas

if you ran the "replay" step with "-O J" in your command line,
you should be able to load the Jacobian profiles from the mc2
file using mcx/utils/loadmc2.m. The Jacobian is time-gated,
similar to the fluence output (4D array with the last index being
the time gate index).

because this feature is currently under development, the
absolute quantity of the output may subject to a scaling factor,
which I have not finalized yet.

Qianqian

You received this message because you are subscribed to the Google Groups "mmc-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mmc-users+...@googlegroups.com.
To post to this group, send email to mmc-...@googlegroups.com.
Visit this group at http://groups.google.com/group/mmc-users.
For more options, visit https://groups.google.com/d/optout.

[Klaas Rackebrandt]

Thanks for your fast answer.

 

Klaas

[Klaas Rackebrandt]

Hi,

 

I´m getting the same error running the example code.... Has anyone solved the problem how to set the DetPos right to get results in the seeds data?

 

I tried to run the sample code without any chanings and dont get why the detectors are not positioned at the escaping surface with the given values.

 

Cheers

 

Klaas