Monte Carlo Simulation Penetration Depth

[Liban Hussein]

Hello,

I've been meaning to reach out regarding a Monte Carlo simulation using the Virtual Photonics framework. From what I've seen, the simulations for the MCCL simulator depend on several input parameters, including source detector separation.

Assuming I just want to make a basic simulation that compares the penetration depths of two laser lights like 450 nm and 520 nm (blue and green, respectively) based on their optical power densities, would the Virtual Photonics tool support this? As mentioned, the simulation would be very basic, just showcasing the penetration depth in a 3D multi-layered tissue for two varying wavelengths.

Looking forward to hearing from you soon,

Liban

[Carole Hayakawa]

Hi Liban,

Thanks for reaching out!  Yes, I think the Monte Carlo code could work for you.  Please check out these links:

1) For download of our code:

https://github.com/VirtualPhotonics/Vts.MonteCarlo/releases

and select version per your operating system

2) Documenation:

https://github.com/VirtualPhotonics/Vts.MonteCarlo/wiki

on the right hand panel there are several links to help users.  Please check out first "Getting Started" to learn how to download and install .net 6.0. 

We have several detectors that determine "MaxDepth", for example ROfRhoAndMaxDepthDetector.  An example of how to define this detector is in our example infile "infile_one_layer_all_detectors.txt".  This detector determines each photon's maximum depth of penetration prior to detection and then tallies it's weight to the appropriate exiting rho bin and max depth bin.  If this sounds like what you are interested in, please give the code a try.  If you have questions about how to set up a simulation or how to post-process these depth-dependent results, please post again.

Best,

Carole

[Liban Hussein]

Hello Carole,

Thanks greatly for following up! I've been able to get the system up and running :)

My only question is, given that there are so many parameters that can be tweaked, is there any way to simplify the "infile_one_layer_all_detectors.txt" file to do the following:
1. Generate multi-layered tissue phantom with user-defined optical properties

2. Vary the input wavelength while keeping the rest of these user-defined properties constant.

It would probably be some variation of the image I have attached below.

For reference, I've fabricated a tissue phantom using the following technique for a laser speckle contrast imaging application (https://escholarship.org/uc/item/94q4r29v). If the simulation parameters allow me to add custom parameters like wavelength, refractive index, layer thickness, etc, then I'm convinced Virtual Photonics is the right tool for my applications!

Thanks greatly,
Liban

[Carole Hayakawa]

Hi Liban,

Great to hear you got something running!  The infile_one_layer_all_detectors.txt was developed to show all users how to specify each detector.  You can cut out all other detectors that are not of interest (just keep track of matching curly braces).  To define a multi-layer tissue, you can add in layers (see example infile infile_two_layer_RORho.txt from running "mc geninfiles" to see how a 2 layer tissue looks like), just make sure that the first tissue layer (2nd layer in file, 1st is air above) "Stop" and the next layer "Start" are equal to define contiguous tissue layers.  Also note there is air as last layer, however you can make last tissue layer (penultimate layer) "Stop" be anything (100mm is equivalent to infinity depending on optical properties).

Now to vary the input wavelength optical properties, I have done this in the past by finalizing my infile and setting the mua, mus, g, n, and mus' to a variable, like a1,s1,g1,n1,p1 for those optical properties in layer one, and so on for all layers in the infile.  Then I wrote a MATLAB tool (I'm on linux) to read spectra from an xlsx file and substitute accordingly.  If you are on linux too, I'd be happy to share.

We also have a command line parameter sweep option (see https://github.com/VirtualPhotonics/Vts.MonteCarlo/wiki/MCCL-Examples) however currently the increments specified are uniform and might not work for your application.

Let me know if you have any more questions.

Best,

Carole

[Liban Hussein]

Hi Carole,

Thank you so much for the swift help!!

I’m using a macOS- would it still be possible to do the xlsx file read?

I can spend some time this weekend working to define the tissue layers and properties. If any questions arise, I’ll keep you posted.

Thanks again and have a great weekend,

Liban 

[Carole Hayakawa]

Hi Liban,

The general process that I use is to have a Matlab script read the optical properties from an xlsx file, copy an infile template to an actual named infile, and then loop through the layers and substitute a1,s1,etc. with  the optical properties using the attached sub_ops.sh.  It runs from a bash shell and uses the command sed.  It might run directly for you on a mac.  If not, I googled it and their might be an equivalent command to use on mac os.

If this sounds like it might work for you, I will create some example files that illustrate this concept.

Best, Carole

[Liban Hussein]

Hello Carole,

Thanks for clarifying this. Yes, I'd love to see if there are example files to illustrate this concept.

At the moment, I've defined two air layers (top and bottom) and three layers mimicking the epidermis, dermis, and hypodermis with varying optical properties. I've attached the corresponding code.

My question for this is, I'm seeing a plot (see attached file) for the infile_two_layer_RORho.txt instead of a graphic of the tissue layers after running the Matlab script- is this supposed to be the case for the "infile_two_layer_RORho.txt", or does the separate file for each wavelength need to be added to see this? I guess my question is, is there a specific output file where I can see the tissue layers defined based on what I have in "infile_two_layer_RORho.txt" after running the Matlab script?

Best,

Liban

[Carole Hayakawa]

Hi Liban,

I see the problem.  In your infile, you need to specify non-zero "Musp"  values for the tissue layers.  The "Mus" and the "Musp" are somewhat redundant because one can be determined from the other if "G" is known, however "Musp" rules, that is, if there is a disagreement between the "Musp" and calculated (1-G)*Mus, then the Mus is calculated from the Musp.  Please update the infile with the appropriate "Musp" values and you should see a difference.

Also, I notice from your screen shot that you are plotting negative rho on the x-axis.  Rho can never be <0 because it relies on a cylindrical geometry.  So like your fluence(rho,z) plot above, just plot rho>=0.

Let me know how it goes.

Best, Carole

[Liban Hussein]

Hi Carole,

This helped and after running the "one_layer_ROfRho_FluenceOfRhoAndZ", I can see the light propagation through the tissue layers. I had a couple last questions:

1. Is there a way to see segmentation between different layers? Figure 4 in this paper for example shows segmentation between layers, although I don't know if the code itself does this or if they added this manually.

2. Is there a way to show a blood vessel in the "one_layer_ROfRho_FluenceOfRhoAndZ" simulations? Something similar to Figure 1A or 1C in this link.

After I verify these, I can work to test each Monte Carlo for each wavelength.

Thanks again,
Liban

[Carole Hayakawa]

Hi Liban,

Concerning your questions:

1) I'm not sure what you mean here.  Do you mean the black dashed lines at the tissue layer interfaces?  If so, you would add those using Matlab commands.  If you mean how the fluence changes at the layer interface, you should see that in the plot you generated.

2) You can add a cylinder running along the y-axis to mimic a blood vessel.  We have a SingleInfiniteCylinderTissueInput specification, the example infile infile_infinite_cylinder_AOfXAndYAndZ.txt shows how to specify this type of tissue.  One caveat is that the cylinder needs to be entirely contained within a layer (no straddling layers).

I have attached a zip file that includes files to use a template infile and run gen_infiles.m using sub_ops.sh to generate infiles for various wavelengths read from a xlsx spreadsheet.

Best, Carole

[Liban Hussein]

Thanks for sharing this Carole and apologies for the delay in my response.

I've been able to resolve some issues and not others.

1. I'm having trouble editing the multi-wavelength simulation to include the "infile_infinite_cylinder_AOfXAndYAndZ.txt"- the simulation doesn't run. Could the file I attached be checked to see if the formatting is correct?

2. Is it possible to change the "infile_infinite_cylinder_AOfXAndYAndZ.txt" cylinder to appear like a blood vessel in terms of shape and color?

Best,

Liban

[Liban Hussein]

Hi Carole,

Apologies for the delay in my response- I've submitted a question a few days back that goes into the specifics of the issue I'm currently facing.

Any help is greatly appreciated,

Liban Hussein

Ph.D. Student | GEM Fellow | ARCS Fellow | NSF GRFP Recipient

Department of Electrical and Computer Engineering

University of Washington-Seattle

Contact: lib...@uw.edu | LinkedIn

--
You received this message because you are subscribed to a topic in the Google Groups "Virtual Photonics" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/virtual-photonics/va8bUcKqT8g/unsubscribe.
To unsubscribe from this group and all its topics, send an email to virtual-photon...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/virtual-photonics/68e805ff-98d0-4ab9-a3ff-a6a9faba874fn%40googlegroups.com.

[Carole Hayakawa]

Hi Liban,

Sorry for the delay in responding, I have been out of work for a bit.  I have edited your infile and now it runs.  I just changed a couple of lines but they are crucial.  Let me know if you have any questions regarding my edits.

Best,

Carole

[Liban Hussein]

Thanks for following up Carole!

In the original folder sent that had the files attached, was there supposed to be a sub_ops.sh file included? When running the program in the correct directory, I see the following error (see attached below).

Any clarification is greatly appreciated.

Best,

Liban

[Carole Hayakawa]

HI Liban,

Sorry for the omission.  Here it is.  You'll see from the first line in the script that it runs in bash shell on linux.  It *could* possibly run on Mac.  Let me know how it goes.

Best, Carole

[Liban Hussein]

Hi Carole,

Thanks for following up.

I see the following error even after changing the permissions on my Mac (see below):

Could you share a screenshot of what the simulation looks like on your end? If I can see the propagation of the light through the tissue medium for varying wavelengths with the cylinder resembling the blood vessel, then after resolving this issue ideally this case should be resolved.

Thanks so much,

Liban

[Carole Hayakawa]

Hi Liban,

I forgot a critical step.  You need to give "execute" privileges to that file.  In linux the command is:

chmod 755 sub_ops.sh

Or (sudo chmod 755 sub_ops.sh) if you don't own file.  Once you do this, on linux if you type:

ls -l sub_ops.sh

you should see
-rwxr-xr-x 1 hayakawa hayakawa 35 Oct 24 10:46 sub_ops.sh*

The code at the start is -=not a directory,r=read,w=write,x=execute in 3s for (from left to right) owner,group,world.

I hope you can do this on Mac.  

Sorry for the omission.  Let me know how it goes.

Best, Carole

[Liban Hussein]

Thanks Carole,

I did the chmod before submitting the question, but realize the error showed up as a result of the formatting for sub_ops.sh being unix based and not being recognized by Mac OS. That error is now resolved and I appreciate the help.

I had some other general questions:

1. When I run the gen_infile.m program, it generates the infile_lam850.txt file. However, I'd like to know how to interface this infile with the load_results_script.m file- there seem to be if else statements for generating the actual plots.

2. I also want to clarify the following- with the infile_lam850.txt being run and plotted, will this display both layers, the light propagation through the layers and include the cylinder shape resembling the blood vessel? I want to use this to compare between varying wavelengths and optical properties and their penetration depths through a blood vessel.

3. If possible, could you show what the result from the infile_lam850.txt may look like for reference? That way I can tweak and have a better understanding of the results I should expect.

Thanks greatly,

Liban

[Carole Hayakawa]

Hi Liban,

I'm glad you got the sub_ops.sh file to work.  Here are my answers to your questions:

1) Once you generate the infile_lam850.txt you need to run it.  Read about running a simulation here:

https://github.com/VirtualPhotonics/Vts.MonteCarlo/wiki/MCCL-Getting-Started-on-Mac

After running you should have a folder with the results created named whatever you substituted for "out1"_combined.  Following the directions on the above link and modifying "datanames" to this name you should be able to run load_results_script.m and see Fluence(x,y,z) plot.  Editing load_results_script.m in the section that starts:

if isfield(results{di}, 'FluenceOfXAndYAndZ') && show.FluenceOfXAndYAndZ

should allow you to view plots of interest.

2) I looked closer at your infile and notice something that won't work.  The cylinder need to be entirely contained within a layer (I should have a validation check on this and currently don't -> will add in the future).  So you have a 2-layer system with a top layer of thickness 1mm and bottom layer of thickness 99mm.  The cylinder is specified with Center=(X=0,Y=0,Z=1) and Radius=1.  So the cylinder straddles the top and bottom layer.  Is it possible to put the cylinder entirely within the top or bottom layer?

3) I don't have a similar result that I've created, but you should see the "fluence" within the 3-region system.  You have a point source, perpendicular to the tissue and so the fluence with be greatest near the source and diminish as you move away.  Of course, the optical properties of the layers and cylinder will affect how fast that decay will be.

Let me know if you have any questions about this.

Carole

[Liban Hussein]

Hello Carole,

Thanks for the continued help. I've been able to get the code to run, but am still having trouble adding a second cylinder 3 mm deep to mimic a blood vessel. What I'm seeing is the following (see screenshot).

How can the infile_850.txt file attached be updated to add a cylinder 3 mm deep mimicking a blood vessel? I've tried to use the "infile_infinite_cylinder_AOfXAndYAndZ.txt" to add the cylinder but have been coming across issues with the compatibility.

Any clarification is greatly appreciated,

Liban

[Carole Hayakawa]

Hi Liban,

I see the problem.  You have specified a tissue type of ["TissueType": "SingleInfiniteCylinder"'] however have specified 2 infinite cylinders.  We have a MultiConcentricInfiniteCylinderTissue but that places the one cylinder inside the other.  I think you'd like to have a MultiInfiniteCylinderTissue which has two infinite cylinders (not intersecting) within a single layer of tissue.  We don't have that tissue code yet.  If that is what you are interested in, I can work on that.

Happy Holidays!

Carole

[Liban Hussein]

Hi Carole,

Thanks for clarifying this. Yes, I'd like to see if it would be possible to place a cylinder, but not within the same single layer of tissue.

The current code has the light source which resembles a cylinder.

I now want to see if it can be included in the code to add a cylinder that resembles a blood vessel at some chosen depth to be confined within a layer.

Looking forward to hearing from you,

Liban

[Carole Hayakawa]

Hi Liban,

Just to make sure of what you'd like, you'd like the ability to define 2 cylinders, 1 within one layer of tissue and another within another layer of tissue.  Is that correct?

I'm not sure I understand your statement: "The current code has the light source which resembles in cylinder."  In your infile_lam850.txt you have specified a directional point source starting in air.

Best,

Carole

[Liban Hussein]

Hello Carole,

Yes, that is correct. I currently have a directional point source starting in air.

I'd now like to define a single cylinder within a layer of tissue 3 mm deep that resembles a blood vessel.

I apologize for the confusion, and appreciate the continued support.

All the best,

Liban

[Carole Hayakawa]

Hi Liban,

Then you should be able to use the SingleInfiniteCylinderTissue.  In your infile_lam850.txt, remove one of the sections that define a "TissueRegionType": "InfiniteCylinder", currently you have 2 (make sure matching curly braces are removed).  Also, the optical properties of either of these InifiniteCylinders are not blood-like (highly absorbing).  In order to see the change to fluence due to the blood vessel you will need appropriate optical properties for the blood vessel at lambda=850nm. 

Let me know if you have questions or problems setting up your infile.

Best,

Carole

[Liban Hussein]

Thanks Carole!

I'll go ahead and try this and will keep you posted on progress.

Hoping you had a restful holiday break,

Liban

[Carole Hayakawa]

You're welcome Liban!  

Just in case this is of interest.  We used this SingleInfiniteCylinderTissue code to determine the simulation results in this paper.  The wavelength range is different from what you are interested in but it gives you an idea of the kinds of plots that we generated from the simulation.

Best, Carole

[Liban Hussein]

Thanks so much for sharing this Carole- these results look great!

I'm trying to update the infile_850.txt to run as the "FluenceOfRhoAndZ". I've attached it below but see some errors. I've also attached a visual graphic of the phantom I'm trying to simulate with the Monte Carlo- the diameter of the cylinder mimicking the blood vessel is also mentioned.

Error:

The given key 'SingleInfiniteCylinder' was not present in the dictionary.

Is there anyway I can use the "FluenceOfRhoAndZ" to show what I have attached in this message? From the paper, the cylinder looks like a front angle view, and I'd like to see if from the fluence image I can have the image confined in 1 layer to resemble a blood vessel.

Thanks so much for your continued help,

Liban Hussein

[Carole Hayakawa]

Hi Liban,

Concerning your diagram and your infile 850_final.txt that you attach, I have a few questions:

1) The source in the diagram (gold triangle) appears to be a converging circular source (like a cone).  Is that what you would like to model?  We have that, we just need to know if it is a cone (circle at surface), what is the circle radius and angle of convergence.  The infile 850_final.txt defines a DirectionalPointSource and this would be a collimated directed straight into the tissue at location (x,y,z)=(0,0,0).  Let me know which you would like to model.

2) The blood cylinder in the diagram (red cylinder) looks to have axis along the z-axis (vertical), and is finite (has circles at ends).  We currently cannot model that.  Possibly I'm misreading your diagram. Our  SingleInfiniteCylinderTissue has the infinite cylinder with axis along the y-axis (into the page).  The paper shows a slice of the fluence(x,y,z) at y=0 (axis we use is x-axis in plane of page horizontally, y-axis in/out of the page, positive z-axis into tissue vertically).  The infile 850_final.txt defines a MultiLayer tissue with no cylinder inclusion.  The infile you posted earlier infile_lam850.txt is closer to what you want I think.

Also, you state: From the paper, the cylinder looks like a front angle view, and I'd like to see if from the fluence image I can have the image confined in 1 layer to resemble a blood vessel.  I'm not sure what you mean here.  If you define FluenceOfXAndYAndZ in the infile, you can slice this 3D matrix any way you want.

Please give me more details about questions in 1) and 2).  If the blood vessel has axis along the y-axis,  I can create an infile for you that would model what is in the diagram.

Best,

Carole

[Liban Hussein]

Hello Carole,

Thanks greatly for following up. To clarify your questions:

1. This was a mistake on my part for the source- I chose an arbitrary shape for the light source. This should be the normal shaped light source from the code (collimated light directed at tissue location (x,y,z)=(0,0,0).

2. I've attached an updated diagram for the blood vessel where the vessel is modeled into the page. The goal is to show the fluence for different wavelengths and different penetration depths, and if the blood vessel is from a view where it goes into the page similar to the paper, this is sufficient.

Please let me know if this clarifies things or if there's anything else I can provide!

Thanks again for the continued support,

Liban H.

[Carole Hayakawa]

Hi Liban,

I understand diagram better now.  I am working on infile for you.  A couple of questions.  

1) Is there any reason you have separated Layer 1 from Layer 2?  Since they have the same optical properties, this is not needed. 

2) From your diagram it looks like the cylinder center axis is at z=1.5mm?

3) You show rings for fluence, however if the blood cylinder is going into the page, then this breaks cylindrical symmetry and we need to define a Cartesian fluence.  And if you'd like to zero in on region around blood vessel, selection of the fluence (x,y,z) indices will allow you to do this.  Do you have optical properties that you'd like to define for blood?

Best,

Carole

[Liban Hussein]

Hi Carole,

Thanks for your swift response!

1. I'm thinking about developing a future iteration of a phantom I have where the optical properties will differ between the two layers- for the time being I have two layers with the same optical properties, but if it makes it easier for demonstrative purposes to have a single layer 2 mm thick and a block layer below that is thicker, then this is sufficient.

2. Ah, my apologies, you are correct. The top of the blood vessel would be right below the second layer (and subsequently at the top of the third layer). Thanks for the catch.

3. Thanks for clarifying this- the goal would be to show the vessel and the cartesian fluence depth for varying wavelengths on the blood cell. Once verified, one more blood vessel can be placed, but this isn't a priority. As for optical properties:
mu_a--> 0.02

mu_s (reduced scattering)--> 0.05

mu_s (scattering coefficient)--> 5

g--> 0.98

N--> 1.42

I will adjust these accordingly for varying wavelengths after testing with an initial wavelength.

All the best,

Liban

[Carole Hayakawa]

Hi Liban,

I have a draft of an infile for you (see attached).  This specifies a point source.  Three tissue layers (thicknesses 1mm, 1mm, 25mm) currently with the same optical properties.  An infinite cylinder in second tissue layer centered at 1.5mm with radius 0.49mm (I reduced this so that there is no intersection with the layer planes) with optical properties you give above but with g=0.99 (which with mu_s' (reduced scattering)=0.05/mm will give mu_s(scattering coeff)=5/mm).  I had to give the blood region the same N=1.4 as the surrounding layer (that mismatch capability might be coming later).  The infile specifies an absorbed energy detector in Cartesian coordinates.  I typically use about a 0.2mm resolution in the x,y and z grid.  I also made the y bin be a single bin.

I just ran 100 photons.  This creates a folder named "lam850_update231108".  I edited load_results_script.m and put this string in "datanames", i.e.

datanames = { 'lam850_update231108' };

and ran using Octave (should also work using MATLAB).  Attached is a plot of the absorbed energy.

Give this a try and make sure it works on your side.  If you have further questions, don't hesitate to ask.

Best,

Carole

[Liban Hussein]

Hello Carole,

Thanks so much for providing this and apologies for the delay in my response.

I've been able to run the code, with a small change- I ran the code with both the 1.5 mm centering for the blood vessel but feel like the point source covered the blood vessel up.

I then ran a 10,000 photon run shifting the blood vessel to be located at 9 mm. However (see attached image), I am only seeing the point source and not the blood vessel.

Can this be clarified?

Best,

Liban

[Carole Hayakawa]

Hi Liban,

I'm on holiday for the next two weeks so am keeping it brief.  I suggest that you run a simulation without the blood vessel and everything else the same.  This would mean defining a MultiLayerTissue instead of the SingleInfiniteCylinderTissue.  Look for an example of that tissue portion of the infile in infile_one_layer_all_detectors and edit with your optical properties.  Then plotting a ratio (or difference) between the absorbed energy plots at the same y=0 slice should show the effects of the blood vessel and remove the strong effects of the source above the cylinder.

Keep me posted.

Best,

Carole

[Carole Hayakawa]

Hi Liban,

This was bothering me so I debugged it.  Turns out the infile has some subtle errors that need to be adjusted (see new infile).  In "TissueInput", the "InfiniteCylinderRegions" needs to be "InfiniteCylinderRegion" (no ending "s") because the tissue type is SingleInfiniteCylinder (not MultiConcentricInfiniteCylinder).  And the square brackets need to be removed indicating more than one InfiniteCylinderRegion.  Since these were in error, the InfiniteCylinderRegion properties were in error (e.g. Center and Radius).  I ran this infile and plotted the results (see attached jpg) and now you can see the cylinder.

I'm sorry about the confusion and leading you astray with my prior infile post.

Let me know if you have any further questions.

Carole

[Liban Hussein]

Carole,

Thank you so much for your help and for clarifying this, and no worries at all! It means a lot that you've taken time to help me resolve this, and I'm very grateful.

I'll go ahead and play around with this further, and wish you a restful and well deserved break.

All the best,

Liban