VirtualPhotonics / VTS

Virtual Tissue Simulator
https://virtualphotonics.org
Other
34 stars 9 forks source link

MC detectors with NA #18

Closed hayakawa16 closed 1 year ago

hayakawa16 commented 5 years ago

Goal: to enable defining: (a) detector NA capabilities in our current detectors, and (b) a fiber detector with a index of refraction, n, and numerical aperture, NA, adjacent to the tissue. That is, change the tissue-air interface in the location of the detector to determine reflection and refraction based on the n of the detector. This is the first time a detector definition would modify the photon transport in the tissue.

Proposed solution: after reviewing the code classes and trying not to mangle classes, I came up with a plan. Tissue typically has no information about detectors, however, to modify transport, tissue needs to know that a detector is there instead of air. So if a fiber detector is defined, tissue would have a infinitesimally thin TissueRegion at the tissue-air interface under the detector that has different OPs.

I wrote out specs for the code. Please feel free to challenge any one of them.

Specs: 1) Add a ReflectanceFiberDetector and TransmittanceFiberDetector that has parameters: a) center b) radius c) refractive index d) NA These detectors would be the only detectors that modify the photon transport along the tissue-air surface. The rest of the resident detectors would assume a consistent tissue-air interface.

2) For the resident detectors enable a NA parameter but assume detection is always in the air. This would mean users would specify NA of detector relative to n=1.0. I came to this conclusion because it wouldn't make sense for say ROfRho (or all of the rest of the resident detectors) to be adjacent to the tissue since then the entire tissue surface would interface with the n of the detector and that didn't make sense to me.

3) Currently Photon class CrossRegionOrReflect method asks the tissue to provide the current n and the neighbor region index (index within the array of regions). This method then determines the n of the neighbor region. To do this, tissue checks the layers in MultiLayerTissue or the inclusion in SingleInclusionTissue to determine the neighbor index.

Proposal: If ReflectanceFiberDetector of TransmittanceFiberDetector are defined, also modify TissueInput to define a infinitesimally thin DetectorRegion that has center, radius, n and NA. Then when CrossRegionOrReflect asks tissue for the neighbor index, tissue will know the index of the DetectorRegion and Photon can then know the n of that region. No change to Photon class would be needed.

The DetectorRegion would be handled internal to the code when an infile specifies a ReflectanceFiberDetector. It would have dimensions based on the user input for this detector.

dcuccia commented 5 years ago

I think I understand the need(s). "Coupling" the detector and tissue in an elegant way sounds tricky.

What is the limitation keeping us from wrapping both index and NA into just the detector definition (defining acceptance criteria - angle and position)? Would be best if the tissue didn't have to worry about how the photons we're being detected (on the flip side, totally fine of the detector "knows" what tissue region it's in. Imagine the case of one homogeneous tissue with 100 different detectors placed in various locations to look at fluence/flux for dosimetry. Seems silly to pollute the elegance of the tissue with the detail of all the light detectors.

I believe we've done "2D surface" detectors before, which could model a fiber - what can we learn from what worked/didn't work there?

Other alternative is that we actually model the physical fiber as part of the tissue structure, separate from the detector aspect, if we need to model the fiber-tissue interaction (reflects photons outside of a certain angle, within that angle infinite absorption).

Whatever we come up with should work inside and outside of the tissue (submerged fiber in intralipid, e.g.)

On Fri, Apr 12, 2019, 11:45 AM Carole Hayakawa notifications@github.com wrote:

Goal: to enable defining: (a) detector NA capabilities in our current detectors, and (b) a fiber detector with a index of refraction, n, and numerical aperture, NA, adjacent to the tissue. That is, change the tissue-air interface in the location of the detector to determine reflection and refraction based on the n of the detector. This is the first time a detector definition would modify the photon transport in the tissue.

Proposed solution: after reviewing the code classes and trying not to mangle classes, I came up with a plan. Tissue typically has no information about detectors, however, to modify transport, tissue needs to know that a detector is there instead of air. So if a fiber detector is defined, tissue would have a infinitesimally thin TissueRegion at the tissue-air interface under the detector that has different OPs.

I wrote out specs for the code. Please feel free to challenge any one of them.

Specs:

1.

Add a ReflectanceFiberDetector and TransmittanceFiberDetector that has parameters: a) center b) radius c) refractive index d) NA These detectors would be the only detectors that modify the photon transport along the tissue-air surface. The rest of the resident detectors would assume a consistent tissue-air interface. 2.

For the resident detectors enable a NA parameter but assume detection is always in the air. This would mean users would specify NA of detector relative to n=1.0. I came to this conclusion because it wouldn't make sense for say ROfRho (or all of the rest of the resident detectors) to be adjacent to the tissue since then the entire tissue surface would interface with the n of the detector and that didn't make sense to me. 3.

Currently Photon class CrossRegionOrReflect method asks the tissue to provide the current n and the neighbor region index (index within the array of regions). This method then determines the n of the neighbor region. To do this, tissue checks the layers in MultiLayerTissue or the inclusion in SingleInclusionTissue to determine the neighbor index.

Proposal: If ReflectanceFiberDetector of TransmittanceFiberDetector are defined, also modify TissueInput to define a infinitesimally thin DetectorRegion that has center, radius, n and NA. Then when CrossRegionOrReflect asks tissue for the neighbor index, tissue will know the index of the DetectorRegion and Photon can then know the n of that region. No change to Photon class would be needed.

The DetectorRegion would be handled internal to the code when an infile specifies a ReflectanceFiberDetector. It would have dimensions based on the user input for this detector.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/VirtualPhotonics/VTS/issues/18, or mute the thread https://github.com/notifications/unsubscribe-auth/AAdRgVX5raJtUZrBLZFjyG0ehbU0Axhnks5vgNRlgaJpZM4cstqf .

hayakawa16 commented 5 years ago

David, thank you very much for your feedback. Since you were part of the original design, I knew you'd catch consistencies and also know the best code design.

The paradigm shift that I'm wrestling with is that all of our detectors are currently "virtual". They are like surface or volume sieves tallying photon trajectory activity but not changing that trajectory. So, for instance, the RadianceOfRhoAtZDetector was defined to tally radiance as a function of rho at a certain depth Z within the tissue. This was actually designed for dosimetry estimates not for a submerged fiber. It is an internal surface detector that tallies as photons go by the surface. There is no physical detector defined. I don't think anyone has used this detector so no lessons learned so far.

To define detectors that mimick reality closer, then the detectors will need to be part of the tissue structure and take a more "real" role in that they will change the photon trajectory. This would be your "other alternative". So a submerged fiber in intralipid would have a vertical cylinder tissue region representing the fiber and photons would interact with this region.

This is somewhat consistent with my proposal for the fiber at the surface except here there is a virtual boundary at the tissue-air interface, so I thought that the fiber would lie between the tissue and that virtual boundary and be infinitesimally thin. I thought of having an air layer above the tissue that has a vertical cylinder representing the fiber and transport would continue into that layer. I didn't like this idea because the tally is occurring solely at the fiber tissue interface and no other.

You gave me good food for thought. I'll ponder some more.

dcuccia commented 5 years ago

"all of our detectors are currently "virtual"...tallying photon trajectory activity but not changing that trajectory."

I know what you mean, but technically I disagree: if a reflectance detector tallies (is within acceptance angle of detector, survives the reflection off the air-tissue surface due to index mismatch), the photon is killed. That's a pretty significant change in trajectory (both angle and weight are modified). If a fiber tallies (is within acceptance angle of detector, survives reflection/index mismatch), I assume we all agree it should be killed. Whether or not you choose to augment the actual volume with more complex/more realistic inclusions would depend on what kind of simulation/experiment you're trying to accomplish.


David Cuccia, PhD

Chief Executive Officer

Chief Technical Officer

Office: 949-825-5070 | Mobile: 949-697-1982

Modulim https://modulim.com/ | 2400 Barranca Parkway | Irvine | CA 92606

On Tue, Apr 16, 2019 at 11:54 AM Carole Hayakawa notifications@github.com wrote:

David, thank you very much for your feedback. Since you were part of the original design, I knew you'd catch consistencies and also know the best code design.

The paradigm shift that I'm wrestling with is that all of our detectors are currently "virtual". They are like surface or volume sieves tallying photon trajectory activity but not changing that trajectory. So, for instance, the RadianceOfRhoAtZDetector was defined to tally radiance as a function of rho at a certain depth Z within the tissue. This was actually designed for dosimetry estimates not for a submerged fiber. It is an internal surface detector that tallies as photons go by the surface. There is no physical detector defined. I don't think anyone has used this detector so no lessons learned so far.

To define detectors that mimick reality closer, then the detectors will need to be part of the tissue structure and take a more "real" role in that they will change the photon trajectory. This would be your "other alternative". So a submerged fiber in intralipid would have a vertical cylinder tissue region representing the fiber and photons would interact with this region.

This is somewhat consistent with my proposal for the fiber at the surface except here there is a virtual boundary at the tissue-air interface, so I thought that the fiber would lie between the tissue and that virtual boundary and be infinitesimally thin. I thought of having an air layer above the tissue that has a vertical cylinder representing the fiber and transport would continue into that layer. I didn't like this idea because the tally is occurring solely at the fiber tissue interface and no other.

You gave me good food for thought. I'll ponder some more.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/VirtualPhotonics/VTS/issues/18#issuecomment-483799865, or mute the thread https://github.com/notifications/unsubscribe-auth/AAdRgY68WaBbv3DUE0-6lWXpEbrWw-uiks5vhhxagaJpZM4cstqf .

hayakawa16 commented 5 years ago

I think we agree about current reflectance tallies. But just to recap. We run a photon trajectory and all reflections/refractions at the tissue surface are based on comparing the n=1.0 of air with n of tissue. This progresses until the photon refracts out of the tissue and is killed (so a direction change but no weight change - the refraction is a pseudo-collision). Then the photon trajectory that occurred is "replayed" for any detector definitions and for ROfRho say, the exiting location is checked to see if it is within Rho (and the refracted direction is checked against NA on the branch code). So the detector plays no part in the photon trajectory.

If we were to have a fiber adjacent to the surface, then the photon trajectory would check reflections/refractions based on n of the detector if the photon is at the detector location. So if I were to run the photon trajectory described above, the fiber n (which is most likely >1.0) could potentially change what happens to the photon when the photon is at the surface at the detector. For example, if the fiber were not there the photon could possibly reflect at that location, however, if the fiber is there and say had an n that matched the tissue, it would pass through to the detector and be killed.

Are we with the same understanding or am I missing a subtlety of your comment?

hayakawa16 commented 5 years ago

I have discussed this over with others and have come to the following plan: For fiber detectors, the user can define a CylinderTissueRegion as part of the tissue definition so that the photon transport will interact with this region. We already have a CylinderTissueRegion class with properties Center, Radius, Height and OPs. For a fiber that is at the surface, the Height specification should be set to 0. I will code up RayIntersectBoundary in CylinderTissueRegion (it is currently not implemented) and handle the usual interaction with a cylinder. I will add in a case for when the height is 0 so that the photon interacts with a circle at the surface. I think we have an assumption that any layer inclusion lies entirely within the layer. I'll see if I run into this problem coding this up. With this plan, the user can define whether the detector is there virtually or realistically. It will also accommodate defining a realistic source for example, for the one fiber source we have, LambertialSurfaceEmittingCylindricalFiberSourceInput. I plan to proceed with this code on the branch. I will post my status as I go.

hayakawa16 commented 5 years ago

I have progressed with the code development. I am now working on the VirtualBoundary related to the detector. I was going to use the RadianceVirtualBoundary for the submerged fiber detector, CylindricalFiberDetector, because they are both internal surface detectors. The one detector we have associated with this VB is RadianceOfRhoAtZ which was designed with dosimetry in mind. The VB has code particular to this detector within it (i.e. checking if Uz>0 downward flow) which makes it different from all the rest of the VBs. If I use this same VB for CylindricalFiberDetector, I would need to add code particular to this detector.

As I'm writing this, I think the best thing to do is add a VirtualBoundary for the submerged fiber and not try to merge it with RadianceVirtualBoundary. I plan to proceed with this plan. I think we should rename RadianceVB to DosimetryVB though or something like that because it doesn't handle all internal surface radiance detectors, it is only capturing downward radiance.

hayakawa16 commented 5 years ago

I have kind of gone full circle with this effort. First, I tried to code up a literal CylindricalFiberDetector that could be at the surface or submerged. This required an equivalent tissue addition of a CylinderTissueRegion that represented the fiber. I ran into problems with this code. Technically, all inclusions are assumed to be contained within the layer. This assumption and trying to get the surface fiber to work caused me to ditch this code.

So I restarted and reread David's last post. In terms of how our software is architected, I think his idea has merit and so I coded it up. I added a reflection detector SurfaceFiberDetector that attached to the DiffuseReflectance Virtual Boundary. This coded up very easily because of the structure we have in place. However, all detectors that hang off the ReflectanceVB are tallying because they escaped the surface after passing Fresnel comparison of tissue n with air n. I added a check in the new detector class to compare it's n with that of the tissue and called Fresnel to see if I should tally, but all of the prior pseudo-collisions that bounced off the surface due to the tissue/air n comparison could have tallied to the detector if they were under it.

This is what I did to remedy this. I made SurfaceFiberDetector a IHistoryDetector which attaches to a GenericVolume VB. The key difference with a IHistoryDetector is that it goes back through all of the photon collisions and tallies accordingly. This allowed me to retrace and see if any prior pseudo-collisions occurred at the surface and check if under detector and if so, Fresnel compare the tissue n with the detector n. If transmitted to the detector, I stop processing that photon.

I have written unit tests and stepped through the code. It behaves correctly. I plan to build an executable to test.

hayakawa16 commented 5 years ago

The method I posted in the prior post did not work. After returning to my old C code and coding up what I thought was needed, I got results that agreed with theory. So to code up something similar in the C# code, I created a MultiLayerWithSurfaceFiberTissue class. This class is similar to MultiLayerTissue except that a circular Tissue region, SurfaceFiberTissueRegion, can be defined that exists at the surface of the tissue. Then I modified GetRegionIndex and GetNeightRegionIndex to check if the photon is at the surface within the circular region. The key interaction with Photon class is in CrossRegionOrReflect asking the tissue: int neighborIndex = _tissue.GetNeighborRegionIndex(this); double nNext = _tissue.Regions[neighborIndex].RegionOP.N; I needed to have the N of the fiber to be nNext when crossing the surface under the fiber. This gave me results that were more in line with theory. I'll see if Owen finds that it works for his test cases or not.

dcuccia commented 1 year ago

How did this go, ultimately?

hayakawa16 commented 1 year ago

The ability to define a detector NA for both reflectance and transmittance tallies has been done for a while. I also added a surface fiber detector with user-specified radius, NA, and final tissue region index (so can have adjacent to tissue or air gap). I also created a surface fiber tissue region so that if adjacent, I could handle the refractive mismatch between the fiber and the tissue. I used this code to try to recreate results that a user was interested in, however never obtained a complete match. I reached out to the authors of the results we wanted to match to check that they took such care with the MC simulation under the fiber and they thought the coder most likely did.

hayakawa16 commented 1 year ago

I feel this issue was closed when branch branch ch-141202-detector-NA-changes was merged into master.