1+ # This is an example of python code using VTS to Compute fluence for a two-layer
2+ # medium as a function of radial extent and depth at a given set of optical properties
3+ #
4+ # Import PythonNet
5+ from pythonnet import load
6+ load ('coreclr' )
7+ import clr
8+ # Import the Operating System so we can access the files for the VTS library
9+ import os
10+ file = '../libraries/Vts.dll'
11+ clr .AddReference (os .path .abspath (file ))
12+ import numpy as np
13+ import plotly .graph_objects as go
14+ import plotly .express as px
15+ from Vts import *
16+ from Vts .Common import *
17+ from Vts .Extensions import *
18+ from Vts .Modeling .Optimizers import *
19+ from Vts .Modeling .ForwardSolvers import *
20+ from Vts .SpectralMapping import *
21+ from Vts .Factories import *
22+ from Vts .MonteCarlo import *
23+ from Vts .MonteCarlo .Sources import *
24+ from Vts .MonteCarlo .Tissues import *
25+ from Vts .MonteCarlo .Detectors import *
26+ from Vts .MonteCarlo .Factories import *
27+ from Vts .MonteCarlo .PhotonData import *
28+ from Vts .MonteCarlo .PostProcessing import *
29+ from System import Array , Double , Object , Func , Math
30+ clr .AddReference ("System.Core" )
31+ from System .Linq import Enumerable
32+
33+ solver = TwoLayerSDAForwardSolver ()
34+ solver .SourceConfiguration = SourceConfiguration .Distributed
35+ topLayerThickness = 5
36+ opRegions = Array .CreateInstance (IOpticalPropertyRegion , 2 )
37+ opRegions [0 ] = LayerOpticalPropertyRegion (DoubleRange (0 , topLayerThickness , 2 ), OpticalProperties (0.1 , 1 , 0.8 , 1.4 ))
38+ opRegions [1 ] = LayerOpticalPropertyRegion (DoubleRange (topLayerThickness , Double .PositiveInfinity , 2 ), OpticalProperties (0.01 , 1 , 0.8 , 1.4 ))
39+ # Create the DoubleRange instance
40+ rhos_range = DoubleRange (0.1 , 19.9 , 100 ) # range of s-d separations in mm
41+
42+ # Convert to .NET array
43+ rho_delta = rhos_range .GetDelta ()
44+ print (rho_delta )
45+ # Start at 0.1, increment by 0.2, 100 elements
46+ rhos = 0.1 + rho_delta * np .arange (100 )
47+ print (rhos )
48+
49+ zs_range = DoubleRange (0.1 , 19.9 , 100 ) # range of depths in mm
50+
51+ # Convert to .NET array
52+ zs_delta = zs_range .GetDelta ()
53+ print (zs_delta )
54+ # Start at 0.1, increment by delta, 100 elements
55+ zs = 0.1 + zs_delta * np .arange (100 )
56+ print (zs )
57+
58+ arrayOP = Array [Array [IOpticalPropertyRegion ]]([opRegions ])
59+ # predict the tissue's fluence(rho, z) for the given optical properties
60+ fluenceOfRhoAndZ = solver .FluenceOfRhoAndZ (arrayOP , rhos , zs );
61+ #print("*********************************** FLUENCE *****************************************")
62+ #print(list(fluenceOfRhoAndZ))
63+
64+ allRhos = np .concatenate ((- rhos [::- 1 ], rhos ))
65+ print ("*********************************** RHOS *****************************************" )
66+ print (allRhos .tolist ())
67+
68+ # log transform
69+ log_fluence = [Math .Log (f ) for f in fluenceOfRhoAndZ ]
70+
71+ size = len (zs )
72+ # split into rows
73+ fluenceRowsToPlot = np .array ([log_fluence [i :i + size ] for i in range (0 , len (log_fluence ), size )])
74+ #print("*********************************** FLUENCE ROWS TO PLOT *****************************************")
75+ #print(fluenceRowsToPlot)
76+
77+ # reverse and concatenate
78+ allFluenceRowsToPlot = np .concatenate ((fluenceRowsToPlot [::- 1 ], fluenceRowsToPlot ))
79+ #print("*********************************** ALL FLUENCE ROWS TO PLOT *****************************************")
80+ #print(allFluenceRowsToPlot.tolist())
81+
82+ def heatmap (values , x , y , x_label = "" , y_label = "" , title = "" ):
83+ """Create a heatmap chart."""
84+ # values should be a 2D array-like (list of lists or 2D numpy array)
85+ fig = go .Figure (data = go .Heatmap (
86+ z = values ,
87+ x = x ,
88+ y = y ,
89+ transpose = True ,
90+ colorscale = 'Hot' ,
91+ colorbar = dict (title = title )
92+ ))
93+ fig .update_layout (
94+ title = title ,
95+ xaxis_title = x_label ,
96+ yaxis_title = y_label ,
97+ yaxis_autorange = 'reversed'
98+ )
99+ return fig
100+
101+ fluenceChart = heatmap (allFluenceRowsToPlot .tolist (), allRhos .tolist (), list (zs ), "ρ [mm]" , "z [mm]" , "log(Φ(ρ, z) [mm-3])" )
102+ fluenceChart .show (renderer = "browser" )
0 commit comments