
*******DISCUSSION THUS FAR*********
Original Question:
Hi Everyone,
I am trying to figure out how much my photon reflectance (and momentum transfer) drops if my medium is only about 1cm thick. When I run this using your online application, I specify the air layer from Infinity to 0mm, then the 1st
and 2nd layers I make the same optical properties and from 0mm5mm, 5mm10mm, respectively. The last layer is air from 10mm to Infinity. Optical properties set to 0.0019mua 1musp, which by my calculation, would give 1/e attenuation at 1cm thickness
using the typical 1/(n*3*mua*musp). When I do this with the Monte Carlo I hardly see any change in reflectance from the case where the tissue is infinitely thick. To test this limitation, I even made the thickness of the tissue layers only 0.1mm thick,
and the reflectance still only dropped by about 15% when I did this vs infinitely thick. Do you guys know what might be going on here? Does the transmitted photons get counted in the detector? Thanks for your help,
Tyler
Response from Carole
If you add TDiffuseDetectorInput to your list of DetectorInputs
in the input xml file, it will tally everything that exits the
bottom of the tissue. After you execute the MC, a TDiffuse.xml
will get created that will show how many photons contributed
to this tally (TallyCount) and the average weight (Mean).
My Response
I tried this, but when I load an XML using the VTS website and rerun the simulation, the saved settings XML doesn't reflect any of the changes. It's like it is not even loading it. Do I have to run the simulation separate from
the website? Perhaps I did not format the TDiffuseDetectorInput correctly when entering it in the XML?




Hi Tyler,
I was trying to reproduce the issue you were having and I noticed that the GUI did lock up after I modified the XML file to add the TDiffuse detector input. When I went to debug it seems to have been resolved. I will push a new version of the GUI out to
the web this week and hopefully this will resolve your issue.
Thanks,
Lisa




ok thanks. Also, is there a tally list in the .xml? Is there a way I can ask it to tally momentum transfer?
On Wed, Feb 15, 2012 at 12:49 PM, lmalenfant <notifications@codeplex.com> wrote:
From: lmalenfant
Hi Tyler,
I was trying to reproduce the issue you were having and I noticed that the GUI did lock up after I modified the XML file to add the TDiffuse detector input. When I went to debug it seems to have been resolved. I will push a new version of the GUI out to
the web this week and hopefully this will resolve your issue.
Thanks,
Lisa




Tyler  Momentum Transfer is implemented but not exposed...if you ran it now, you'd get zeros. Carole wanted to follow up with you to determine what you'd actually need, and we'd get this out. Can you describe here what features you'd like (reflectance vs
fluence vs radiance, solution domain  fx vs rho vs time)?




Sure  what would be useful for me is the histogram of total momentum transfer as a function of rho and/or frequency in each layer. The per photon tally of momentum transfer and rho/fx in each layer and weight would also do. Let me know if you guys make
any progress in this area. Thanks!
Tyler
On Wed, Feb 15, 2012 at 1:12 PM, dcuccia <notifications@codeplex.com> wrote:
From: dcuccia
Tyler  Momentum Transfer is implemented but not exposed...if you ran it now, you'd get zeros. Carole wanted to follow up with you to determine what you'd actually need, and we'd get this out. Can you describe here what features you'd like (reflectance vs
fluence vs radiance, solution domain  fx vs rho vs time)?



Feb 16, 2012 at 11:44 PM
Edited Feb 16, 2012 at 11:56 PM

The new version of the GUI is now available on the website, let me know if you are still having issues with this version and I can take a look at your XML file.




Concerning a new Momentum Transfer detector: I plan to add that this detector would be a nice feature to the Monte Carlo code to the codeplex Issues page and it will get tracked/discussed there.



Feb 22, 2012 at 11:36 PM
Edited Feb 29, 2012 at 8:34 PM

Just to close loop, here is link to Issues page:
http://virtualphotonics.codeplex.com/workitem/33




Tyler  we had an interesting discussion around this topic today and we've come up with a game plan. Out of scientific curiosity, we were curious why specifically a mean and second moment aren't sufficient for the computation, and a true histogram in MT
is required (vs just adding more bins to capture a finer representation of MT mean and second moment).




Hi David,
I haven't thought too much about using the moments of the distribution. I think the main thing is that the quantity of interest in the end is not the weighted average MT but the weighted average decorrelation exponential (=exp(1/3*k^2*MT*msdisplacement)).
Thus I need the full probability distribution of MT to calculate this. If I wanted to use moments I think we would have to tally the decorrelation exponential itself and use the moments of
this quantity. This would be OK but then I think we would have to assume some mean square displacement during the monte carlo  not ideal since this is the quantity to be fit and would be nice to have this outside the monte carlo. I may be way off here
though, let me know if I am not understanding you right.
Tyler
On Wed, Feb 22, 2012 at 5:31 PM, dcuccia <notifications@codeplex.com> wrote:
From: dcuccia
Tyler  we had an interesting discussion around this topic today and we've come up with a game plan. Out of scientific curiosity, we were curious why specifically a mean and second moment aren't sufficient for the computation, and a true histogram in MT
is required (vs just adding more bins to capture a finer representation of MT mean and second moment).




That makes sense  I always forget that expectation values of nonlinear functions aren't the same as nonlinear operations on expected values. In the more familiar cases to us (say reflectance or fluence), finer spacetime bins result in "flatter" distributions,
where the distribution is pretty well represented by a mean. In your case, the histogram is probably highly nonlinear. (Aside: Is there a good reference paper that shows these distributions as a function of radial distance from the source?)
Re: implementation:
 Do you typically set bounds of [0, 1] for your MT histograms, or do you narrow it down (and/or make the bounds spatiallydependent)?
 How many bins do you typically use?




Hi David,
The histogram is indeed nonlinear, and also varies quite a bit depending on radial distance from source. It doesn't look like any of the typical distributions (i.e. Gaussian, Chi2, etc) that I know of. Thus, I don't think we could accurately reconstruct
the distribution with only a couple moments. Incidentally, since you are asking my paper shows this histogram at a couple of source detector separations!
http://www.opticsinfobase.org/abstract.cfm?URI=josaa28102108 . If you are looking for references outside BLI, I know that David Boas has a couple of these plots in his thesis.
I see similar distributions to him, with a peak at dimensionless MT = 2, which would be hard backscattering. Typically my histograms range from [0300] (note we histogram the total momentum transfer after it's scattered a bunch so it will usually be above
1). I use 3000 bins with a spacing of 0.1. For radial bins, I use about 300 with a spacing of 0.05mm. I messed around with it quite a bit and found this to accurately represent the distribution, but there is probably some room for optimization still.
Tyler
On Thu, Feb 23, 2012 at 12:00 PM, dcuccia <notifications@codeplex.com> wrote:
From: dcuccia
That makes sense  I always forget that expectation values of nonlinear functions aren't the same as nonlinear operations on expected values. In the more familiar cases to us (say reflectance or fluence), finer spacetime bins result in "flatter"
distributions, where the distribution is pretty well represented by a mean. In your case, the histogram is probably highly nonlinear. (Aside: Is there a good reference paper that shows these distributions as a function of radial distance from the source?)
Re: implementation:
 Do you typically set bounds of [0, 1] for your MT histograms, or do you narrow it down (and/or make the bounds spatiallydependent)?
 How many bins do you typically use?




Hi Tyler,
I've implemented a momentum transfer detector in our code and would like to test it out with some data that you typically use. If possible, could you please let me know:
 tissue optical properties (e.g. mua, mus, g, n)
 tissue top and bottom layer thicknesses
 rho bins (e.g. start, stop, delta)
 MT bins (start, stop, delta)
 how many photons you run
Thanks!
Carole




Hi Tyler,
Since the momentum transfer results agrees with yours, I have pushed the latest Monte Carlo code to codeplex. It is version 1.0.4 Beta and can be downloaded from http://virtualphotonics.codeplex.com/releases/view/85890
A couple of notes about this new capability:
 The detector results are now divided by N, the number of photons launched. While debugging I gave you results that were not normalized by N, but with this release this normalization is included so that the second moment results make
sense.
 A sample input file specifying this detector is bundled with download, "infile_two_layer_ReflectedMTOfRhoAndSubRegionHist.xml".
 After running the simulation using the command "mc.exe infile=infile_two_layer_ReflectedMTOfRhoAndSubRegionHist.xml", you can view the results using bundled matlab files. By editing the field "datanames" in load_results_script.m to 'two_layer_ReflectedMTOfRhoAndSubRegionHist"and
running it, you should be able to view the MT results. Note that the sample infile specifies a 4layer system, with air above and below a twolayer slab, so 4 plots will come up showing MT and a 5th plot showing R(rho).
Also, in SimulationOptions in the infile you can set "track statistics" to true. This will output a 'statistics.xml' that will let you know how many photons exited out the bottom of the tissue. I think you had a question about this earlier.
Please let me know how you're doing and if you have any questions/comments/requests.
Thanks,
Carole




Hi Carole,
Testing this out a bit now. Can you tell me how the histograms are created on a layer by layer basis? Do you tally the total MT in layer 1 and layer 2 then put the weight of the photon in each histogram in the appropriate bin? Shouldn't all photon weight
be accounted for in each histogram in this case? In the ones I am looking at, the total weight of each histogram changes depending on layer thickness, etc. What happens if the MT for a layer is 0?
On Wed, Apr 11, 2012 at 5:42 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
Hi Tyler,
Since the momentum transfer results agrees with yours, I have pushed the latest Monte Carlo code to codeplex. It is version 1.0.4 Beta and can be downloaded from
http://virtualphotonics.codeplex.com/releases/view/85890
A couple of notes about this new capability:
 The detector results are now divided by N, the number of photons launched. While debugging I gave you results that were not normalized by N, but with this release this normalization is included so that the second moment results make sense.
 A sample input file specifying this detector is bundled with download, "infile_two_layer_ReflectedMTOfRhoAndSubRegionHist.xml".
 After running the simulation using the command "mc.exe infile=infile_two_layer_ReflectedMTOfRhoAndSubRegionHist.xml", you can view the results using bundled matlab files. By editing the field "datanames" in load_results_script.m to 'two_layer_ReflectedMTOfRhoAndSubRegionHist"and
running it, you should be able to view the MT results. Note that the sample infile specifies a 4layer system, with air above and below a twolayer slab, so 4 plots will come up showing MT and a 5th plot showing R(rho).
Also, in SimulationOptions in the infile you can set "track statistics" to true. This will output a 'statistics.xml' that will let you know how many photons exited out the bottom of the tissue. I think you had a question about this earlier.
Please let me know how you're doing and if you have any questions/comments/requests.
Thanks,
Carole




OK I believe the algorithm is just getting the MT from each layer, moving to the appropriate bin, and putting the weight of the photon there. The problem is that we lose some information this way because I can't tell which photons hit both and how much
each photon deposits in the first/second layer, respectively. Is there any way I can access the individual photon histories? I would need to histogram the total MT transfer, and histogram the fraction of that total MT transfer deposited in the first and second
layer. I could achieve this with access to the photon histories.
Tyler
On Thu, Apr 26, 2012 at 1:03 PM, Tyler Bywaters Rice <tybyrice@gmail.com> wrote:
Hi Carole,
Testing this out a bit now. Can you tell me how the histograms are created on a layer by layer basis? Do you tally the total MT in layer 1 and layer 2 then put the weight of the photon in each histogram in the appropriate bin? Shouldn't all photon weight
be accounted for in each histogram in this case? In the ones I am looking at, the total weight of each histogram changes depending on layer thickness, etc. What happens if the MT for a layer is 0?
On Wed, Apr 11, 2012 at 5:42 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
Hi Tyler,
Since the momentum transfer results agrees with yours, I have pushed the latest Monte Carlo code to codeplex. It is version 1.0.4 Beta and can be downloaded from
http://virtualphotonics.codeplex.com/releases/view/85890
A couple of notes about this new capability:
 The detector results are now divided by N, the number of photons launched. While debugging I gave you results that were not normalized by N, but with this release this normalization is included so that the second moment results make sense.
 A sample input file specifying this detector is bundled with download, "infile_two_layer_ReflectedMTOfRhoAndSubRegionHist.xml".
 After running the simulation using the command "mc.exe infile=infile_two_layer_ReflectedMTOfRhoAndSubRegionHist.xml", you can view the results using bundled matlab files. By editing the field "datanames" in load_results_script.m to 'two_layer_ReflectedMTOfRhoAndSubRegionHist"and
running it, you should be able to view the MT results. Note that the sample infile specifies a 4layer system, with air above and below a twolayer slab, so 4 plots will come up showing MT and a 5th plot showing R(rho).
Also, in SimulationOptions in the infile you can set "track statistics" to true. This will output a 'statistics.xml' that will let you know how many photons exited out the bottom of the tissue. I think you had a question about this earlier.
Please let me know how you're doing and if you have any questions/comments/requests.
Thanks,
Carole




Hi Tyler,
I was just about to respond to your previous post. Lets consider one rho bin in the
one tissue layer system. For each photon, the cumulative MT is determined for the entire tissue (one region), if it is >0, the appropriate MT bin is found, and the exiting photon weight is deposited into this MT bin. So here
you can imagine a histogram plot with MT along xaxis and photon weight along yaxis for one rho. If you integrate this plot over MT you should get the reflected photon weight out this rho bin.
Now lets look at the two layer tissue system. The way the code works is similar. For each photon, the cumulative MT
for each layer is determined, if it is >0, the appropriate MT bin is found, and the exiting photon weight is deposited in these two region bins. So here again you get two histogram plots with MT along xaxis and photon weight along
yaxis. But when you integrate these plots over MT you do not get the reflected photon weight out this rho bin because many weights are double counted if there is MT in both layers.
So from your most recent post, it appears you would like to have the exiting photon weight appropriately spread across the layers such that the integration over MT, over both layers gives you back the reflected photon weight. If this is what you want,
lets meet to discuss this.
Carole




Hello,
I am not sure how the photon weight spreading would work, but what I believe I need is the overall MT histogram (for all layers) and then the fraction of contribution from bottom and top as a function of radius. Would this be possible? Could also do it
post analysis with access to individual photon histories. I would be happy to meet any time this week,
Tyler
On Fri, Apr 27, 2012 at 3:50 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
Hi Tyler,
I was just about to respond to your previous post. Lets consider one rho bin in the
one tissue layer system. For each photon, the cumulative MT is determined for the entire tissue (one region), if it is >0, the appropriate MT bin is found, and the exiting photon weight is deposited into this MT bin. So here you can imagine
a histogram plot with MT along xaxis and photon weight along yaxis for one rho. If you integrate this plot over MT you should get the reflected photon weight out this rho bin.
Now lets look at the two layer tissue system. The way the code works is similar. For each photon, the cumulative MT
for each layer is determined, if it is >0, the appropriate MT bin is found, and the exiting photon weight is deposited in these two region bins. So here again you get two histogram plots with MT along xaxis and photon weight along yaxis.
But when you integrate these plots over MT you do not get the reflected photon weight out this rho bin because many weights are double counted if there is MT in both layers.
So from your most recent post, it appears you would like to have the exiting photon weight appropriately spread across the layers such that the integration over MT, over both layers gives you back the reflected photon weight. If this is what you want, lets
meet to discuss this.
Carole




Not sure if you guys have met yet, but I'd be interested in listening in...perhaps Friday? I'm at the OSA BIOMED until Wed night. DC




Friday 910 or 1112 works with me. David, I'd also like to talk to you about the detector for acerussi.




Hi Carole, let's meet at 11then !
On Wed, May 2, 2012 at 1:53 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
Friday 910 or 1112 works with me. David, I'd also like to talk to you about the detector for acerussi.




Hi Carole,
So sorry running late ill be in about 20minutes. Can you still meet?
On Thu, May 3, 2012 at 4:47 PM, Tyler Bywaters Rice <tybyrice@gmail.com> wrote:
Hi Carole, let's meet at 11then !
On Wed, May 2, 2012 at 1:53 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
Friday 910 or 1112 works with me. David, I'd also like to talk to you about the detector for acerussi.




Michael met with us a couple of weeks ago and then again this morning. At the end of our meeting today, Michael, Adam and I came up with a proposed plan for a revised ReflectedMTOfRhoAndSubregionHistDetector. Michael described that the desired
information is a break down of the MT histogram of (# photons) vs (MT) defining what fraction of each histogram element came from each tissue subregion (layer). Now say for a particular rho bin, there are 100 photons that contribute to MT bin [3031].
This would be one histogram element and for each of those 100 photons, the cumulative MT over all its collisions prior to exiting the rho bin summed to a value in [3031]. The desire is to know what fraction of this cumulative MT came from each layer.
The proposed algorithm is the following:
1) Tally the photon weight in an array R(rho,MT). This would get normalized by the number of photons launched and the rho ring area. Summing over the MT bins should give back R(rho).
2) Create an additional array FractionalMT(rho,MT,subregion). When a photon is tallied to the array in 1), determine the fraction of the cumulative MT that occurred in each tissue subregion and tally that fraction in this array. To create a mean
value here, the sum of the fractions in each subregion would need to be divided by the number of photons contributing, i.e. 100. This would provide a break down of the plot in 1) showing for each rho,MT bin, the fraction contribution from each layer.
I welcome any comments or suggestions. For those who were part of our recent discussions, please let me know if I captured our consensus correctly. For those who were not there, please let me know if the information in 2) provides meaningful
information or not.
Thanks!
Carole




Hi Michael,
I coded up the algorithm described in above post dated July 20. I ran 1e5 photons, 100 rho bins in [010]mm, 51 MT bins in [0500]. The tissue is a 4layer tissue air, layer1, layer2, air. Layer 1 is 1mm thick, layer 2 goes from 1mm to
100mm. The matlab plots of the results are shown
here. The first row shows R(rho) and R(rho,MT). The second row shows the averaged fractional MC in each subregion for a given rho (rho midpoint is shown in title of plot).
Let me know what you think.
Carole




Hi Carole,
Not sure I have the answer here, but the fractional diagrams are currently pretty difficult to understand (pretty bad statistics...too many MT bins?). Maybe consider overlaid 1D plots, or bar charts that add up to total fraction, like this: http://www.stattutorials.com/Images/Stackedbar350.jpg
David
can you possibly plot the fractional MT as 1D histograms? fr




Hi David,
I appreciate your comment. I'll see if I can improve it.
Thanks, Carole




How about this? Or multiplied by the reflected photon weight is
here.
Carole




Hi Carole,
I think we are almost there, but what we are really after is the distribution of fractional MT. For example, in the first plot, we see in the second MT bin that about 60% of the MT comes from layer 1, with 40% therefore coming from layer 2. But surely
this is the average split over many photons. What we need is the full distribution. E.g. for this second MT bin, how many photons had a 90%/10% split? How many had 80%/20%? etc. etc. We know this will average out to a 60%/40% split, but because of the nonlinearity
of our function, we can't just use the average, we need the whole distribution
On Tue, Jul 31, 2012 at 12:17 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
How about this? Or multiplied by the reflected photon weight is
here.
Carole




They're both very helpful and each provide unique info. The first is perfect. Can you try graphing the second (weighted) set as two 1D line plots in "semilogy"?
On Tue, Jul 31, 2012 at 12:17 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
How about this? Or multiplied by the reflected photon weight is
here.
Carole




Hi, Thanks for your responses.
Tyler, glad to hear we're converging. I understand what you need. I'll see what I can do.
David, I tried taking the log of the yaxis values for the stacked bar plot and the bars are shown downward since on a log scale, if you know what I mean. It might work and be easier to read compared to 2 separate line plots. I can show you what
that looks like next iteration.
Carole




I used 10 fractional bins (00.1],(0.10.2]...(0.91] and counted how many photons contributed to that bin per subregion. The resulting plots for Layer 1 and Layer 2 are shown
here. In my algorithm, I only add a count of 1 to a fractional bin if there is > 0.0 MT in that subregion. This means that in the first column (Layer 1 is on top, Layer 2 is on bottom)
for the smallest rho, if a photon accumulates 0 MT in the second layer, then a 1 gets added to bin (0.91] for Layer 1, but nothing gets added for Layer 2. The plots are a bit difficult to read for some bars, but the data is there in arrays in matlab
and can be excised for your purposes. I was able to zoom in and see the smaller bars.
David, I tried the log scale on these stacked bar plots, but couldn't make them look correct, possibly because the counts span 1 (0 on log).
Tyler, Michael, let me know what you think. Possibly the exiting photon weight should be tallied rather than 1, not sure.
Carole




Hi Carole, I'd actually suggested two 1D line plots  bars don't "visually add" in log space




Hi Carole, this is great but I think definitely the weight should go there instead of a 1. Otherwise it looks good. Can we change the bin size for the fraction stuff in the xml?
On Wed, Aug 1, 2012 at 11:58 AM, dcuccia <notifications@codeplex.com> wrote:
From: dcuccia
Hi Carole, I'd actually suggested two 1D line plots  bars don't "visually add" in log space




Great feedback.
Tyler, I'll change to using the weight and I'll add the capability for the user to specify the number of fractional bins in the input xml. I'm glad you think the rest of algorithm is okay. I'll post results again for you to view prior to pushing
to codeplex.
David, yes, I understood the 1D line plots. Since the algorithm seems to be settling, I will present some visual options to view the results including that one. Thanks for your suggestions.
Carole




The results using the exiting photon weight are shown
here. The results using semilogy are shown
here. In the second set, I "stacked" the data. In other words, I plotted the fractional value for the first bin, then added in the results of the second bin, etc. This then allows you to see the "total". Not sure it works for you,
I can plot without stacking if that would be of more interest.
Carole




Hi Carole,
I am hoping to renew effort on this. I think these results are exactly what we are looking for. Can you point me to a version of this code I could run to test things out? Thanks,
Tyler
On Wed, Aug 1, 2012 at 7:34 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
The results using the exiting photon weight are shown
here. The results using semilogy are shown
here. In the second set, I "stacked" the data. In other words, I plotted the fractional value for the first bin, then added in the results of the second bin, etc. This then allows you to see the "total". Not sure it works for you, I
can plot without stacking if that would be of more interest.
Carole




Hi Tyler,
Great! I'll put up a new download of Monte Carlo Command Line on our codeplex site tomorrow. It will have the updated code that you can test out.
Carole




I just updated the Monte Carlo Command Line Download on the codeplex site. It is version 1.0.8 Beta. There is a sample infile for the MT detector called infile_two_layer_ReflectedMTOfRhoAndSubregionHist.xml. Let me know if you have any
questions or problems with the new release. Thanks, Carole




Hi Carole,
Everything is working well, and I got the sweeps to work. I have one question though, the fractionalMT has data for layers 14. I don't understand what we are recording here for layer 2 and 3. I would expect all photons that get detected to have a fractional
MT value and this is what gets recorded. What does fractional MT in layer 2 or layer 3 mean? Just that they reached this layer?
On Wed, Oct 10, 2012 at 4:27 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
I just updated the Monte Carlo Command Line Download on the codeplex site. It is version 1.0.8 Beta. There is a sample infile for the MT detector called infile_two_layer_ReflectedMTOfRhoAndSubregionHist.xml. Let me know if you have any questions or problems
with the new release. Thanks, Carole




Hi Tyler,
First a quick question. When you say you have layers 14, did you specify a 2layer tissue with air above and below? If so, then layers 1 and 4 (air) will not have any fractionalMT recorded and should be 0.
The algorithm determines two arrays, Mean[ir,imt] and FractionalMT[ir,imt,isr,ifrac]. The algorithm in the code is the following:
for each photon exiting surface of tissue:
determine rho bin index it exited, ir (if beyond last rho bin, it is put in last bin)
for each collision photon experienced:
determine which tissue subregion (layer) index collision occurred in, isr
determine MT = 1  cosineBetweenTrajectories;
add MT into array subregionMT[isr] += MT
end
totalMT = sum of subregionMT over isr;
if totalMT > 0
determine which MT bin totalMT is in, imt
add photon exiting weight to Mean[ir,imt] += photon weight
for all subregions isr
if (subregionMT[isr] > 0
determine which fractional bin the fraction (subregionMT[isr]/totalMT) is in, ifrac
add photon exiting weight to FractionalMT[ir,imt,isr,ifrac] += photon weight
endif
endfor
endif
endfor
after all photons have been processed, normalize both arrays by the area of the surface rho bin ring and the number of photons launched.
So the FractionMT array for a particular rho and MT bin, has for layers 2,3 the averaged fractional MT accumulated in that layer. I tally the photon exiting weight so that results parallel what is tallied in the Mean array.
Let me know if this algorithm makes sense or not.
Carole




OK I think I've got it. However, if there are 4 layers air  layer 2  layer 3  air, then I would expect that the "fractional MT" in layer 2 to be one minus the fractional MT in layer 3. E.g. if layer 2 has 10 photons with 90% MT from this layer,
then by default layer 3 should have 10 photons with 10% MT from this layer, since these are the only possibilities. With this logic, the total photon weight in the FractionalMT variable in the results for layer 2 and 3 should be the same, and equal to the
total weight in the mean variable. However, I find this is not the case, and I think it has something to do with not including the photons with 0 fraction. I think we should remove the if subregion MT>0 line, and instead have a 0% bin which tracks it. I
might be able to work it out with what is currently in there, will keep you posted.
On Thu, Nov 1, 2012 at 3:07 PM, hayakawa <notifications@codeplex.com> wrote:
From: hayakawa
Hi Tyler,
First a quick question. When you say you have layers 14, did you specify a 2layer tissue with air above and below? If so, then layers 1 and 4 (air) will not have any fractionalMT recorded and should be 0.
The algorithm determines two arrays, Mean[ir,imt] and FractionalMT[ir,imt,isr,ifrac]. The algorithm in the code is the following:
for each photon exiting surface of tissue:
determine rho bin index it exited, ir (if beyond last rho bin, it is put in last bin)
for each collision photon experienced:
determine which tissue subregion (layer) index collision occurred in, isr
determine MT = 1  cosineBetweenTrajectories;
add MT into array subregionMT[isr] += MT
end
totalMT = sum of subregionMT over isr;
if totalMT > 0
determine which MT bin totalMT is in, imt
add photon exiting weight to Mean[ir,imt] += photon weight
for all subregions isr
if (subregionMT[isr] > 0
determine which fractional bin the fraction (subregionMT[isr]/totalMT) is in, ifrac
add photon exiting weight to FractionalMT[ir,imt,isr,ifrac] += photon weight
endif
endfor
endif
endfor
after all photons have been processed, normalize both arrays by the area of the surface rho bin ring and the number of photons launched.
So the FractionMT array for a particular rho and MT bin, has for layers 2,3 the averaged fractional MT accumulated in that layer. I tally the photon exiting weight so that results parallel what is tallied in the Mean array.
Let me know if this algorithm makes sense or not.
Carole

