/
ToolGetMeOutletsBoundaries.py
155 lines (131 loc) · 6.01 KB
/
ToolGetMeOutletsBoundaries.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
#!/usr/bin/env python
## Christophe.Chnafa@gmail.com
import argparse
import numpy
import vtk
import src.ImportData as ImportData
from src.DisplayData import DisplayModel, VtkText
from src.NetworkBoundaryConditions import FlowSplitting
def Program(fileNameCenterline, fileNameModel, outputFlowrates,
displayModel, localRadii, verboseprint):
print
print "--- Input files:"
print "Input centerline file name: ", fileNameCenterline.rsplit('/', 1)[-1]
if not(fileNameModel == ''):
print "Input model file name: ", fileNameModel.rsplit('/', 1)[-1]
print
# Check if the radii used for the flow splitting are computed locally.
localRadii = float(localRadii)
if localRadii > 0.0:
PowerLawUsesLocalRadii = True
else:
PowerLawUsesLocalRadii = False
# Load the centerline vtk data from the file 'fileNameCenterline'.
centerline = ImportData.loadFile(fileNameCenterline)
# Set the corresponding 0D network.
network = ImportData.Network()
ImportData.SetNetworkStructure(centerline, network, verboseprint,
isConnectivityNeeded=True, isLocalRadiiNeeded=PowerLawUsesLocalRadii,
localRadii=localRadii)
# Compute the outlet boundary conditions.
flowSplitting = FlowSplitting()
flowSplitting.ComputeAlphas(network, verboseprint,
PowerLawUsesLocalRadii=PowerLawUsesLocalRadii)
flowSplitting.ComputeBetas(network, verboseprint)
flowSplitting.CheckTotalFlowRate(network, verboseprint)
# Outlets coords and outlets %.
points = vtk.vtkPoints()
scalar = vtk.vtkDoubleArray()
scalarBis = vtk.vtkDoubleArray()
scalar.SetNumberOfComponents(1)
for element in network.elements:
if not(element.IsAnOutlet()):
continue
points.InsertNextPoint(element.GetOutPointsx1()[0])
scalar.InsertNextValue(100.0*element.GetBeta())
polydata = vtk.vtkPolyData()
polydata.SetPoints(points)
polydata.GetPointData().SetScalars(scalar)
if outputFlowrates:
npoints = polydata.GetNumberOfPoints()
print '{:^12} {:^12} {:^12} {:^12}'.format('X', 'Y', 'Z', '% outflow')
for i in range(0, npoints):
x = [0.0, 0.0, 0.0]
polydata.GetPoints().GetPoint(i, x)
print '{:^12.4f} {:^12.4f} {:^12.4f} {:^12.2f}'. \
format(x[0], x[1], x[2], polydata.GetPointData(). \
GetScalars().GetTuple(i)[0])
print
if displayModel:
labelMapper = vtk.vtkLabeledDataMapper()
if vtk.VTK_MAJOR_VERSION <= 5:
labelMapper.SetInputConnection(polydata.GetProducerPort())
else:
labelMapper.SetInputData(polydata)
labelMapper.SetLabelModeToLabelScalars()
labelProperties = labelMapper.GetLabelTextProperty()
labelProperties.SetFontFamilyToArial()
#tprop.SetColor(0.0,0.0,0.0)
labels = vtk.vtkActor2D()
labels.SetMapper(labelMapper)
labelMapper.SetLabelFormat("%2.2f %%")
# GUI text.
text = ''
if not(fileNameModel == ''):
text = ("Model file name: "
+ repr(fileNameModel.rsplit('/', 1)[-1]) + "\n")
else:
text = ("Centerline file name: "
+ repr(fileNameCenterline.rsplit('/', 1)[-1]) + "\n")
text = (text
+ "Q to exit.")
guiText = VtkText(text)
# Create the renderer
renderer = vtk.vtkRenderer()
renderer.AddActor(guiText.text)
renderer.AddActor(labels)
# Read 3D model if necessary.
if not(fileNameModel == ''):
model = ImportData.loadFile(fileNameModel)
opacity = 0.3
renderer.AddActor(DisplayModel().polyDataToActor(model, opacity))
else:
renderer.AddActor(DisplayModel().polyDataToActor(centerline, 1.0))
renderer.SetBackground(.2, .3, .4)
# Set the lights of the renderer
DisplayModel().setLight(renderer)
# Create the RenderWindow and RenderWindowInteractor
windowTitle = "Percents of the inflow - AneuTool version 0.0.1."
DisplayModel().renderWindow(renderer, windowTitle)
if __name__ == "__main__":
'''Command-line arguments.'''
parser = argparse.ArgumentParser(
description = "GetMeProbePoints: get probe points along the centerline.")
parser.add_argument('-v', '--verbosity', action = "store_true", dest='verbosity',
default = False, help = "Activates the verbose mode.")
parser.add_argument('-i', '--inputCenterline', type = str, required = True, dest = 'fileNameCenterline',
help = "Input file containing the centerlines data in a vtk compliant format.")
parser.add_argument('-iModel', '--inputModel', type = str, required = False, default = '', dest='fileNameModel',
help = "Input file containing the 3D model. It is for vizualization purpose only and it is not required.")
parser.add_argument('-nw', '--notWritePoints', required = False, default = True,
dest = 'writePoints', action = "store_false",
help = "Set this argument to true if you want the script to display as an output the outlet flowrates.")
parser.add_argument('-nd', '--notDisplayModel', required = False, default = True,
dest='displayModel', action = "store_false",
help = "Set this argument to true if you want to display the outlets flow of the model.")
parser.add_argument('-localRadii', '--localRadii', required = False, default = 0, type=int,
dest='localRadii',
help = "Instead of averaging a radius along the branches, a local radius can be computed.")
args = parser.parse_args()
if args.verbosity:
print(">")
print("> --- VERBOSE MODE ACTIVATED ---")
def verboseprint(*args):
for arg in args:
print arg,
print
else:
verboseprint = lambda *a: None
# Start the script.
Program(args.fileNameCenterline, args.fileNameModel, args.writePoints,
args.displayModel, args.localRadii, verboseprint)