Skip to content

Commit 87380f8

Browse files
committed
ENH: Add Python test
1 parent 7c1a826 commit 87380f8

File tree

2 files changed

+109
-0
lines changed

2 files changed

+109
-0
lines changed

wrapping/test/CMakeLists.txt

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
set(TESTING_OUTPUT_PATH "${CMAKE_CURRENT_BINARY_DIR}/Output")
2+
make_directory(${TESTING_OUTPUT_PATH})
3+
set(DATA_DIR ${CMAKE_CURRENT_SOURCE_DIR}/../../test/Input)
4+
5+
itk_python_add_test(
6+
NAME MinimalPathExtraction_Real_DSA_01_02_RegularStepGradientDescent_Python
7+
COMMAND DigitalSubtractionAngiographyVesselPath.py
8+
${TESTING_OUTPUT_PATH}/Real-DSA-01-02-RegularStepGradientDescentPython.png
9+
${DATA_DIR}/Real-DSA-01-Speed-02.mhd
10+
${DATA_DIR}/Real-DSA-01.path
11+
3.0
12+
4000
13+
0.25
14+
0.5
15+
)
Lines changed: 94 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,94 @@
1+
#!/usr/bin/env python
2+
3+
import sys
4+
from pathlib import Path
5+
6+
import itk
7+
import numpy as np
8+
9+
if len(sys.argv) < 8:
10+
print('Usage: ' + sys.argv[0] + ' OutputFilename SpeedFilename PathFilename TerminationValue NumberOfIterations StepLengthFactor StepLengthRelax')
11+
sys.exit(1)
12+
output_filename = sys.argv[1]
13+
speed_filename = sys.argv[2]
14+
path_filename = sys.argv[3]
15+
termination_value = float(sys.argv[4])
16+
number_of_iterations = int(sys.argv[5])
17+
step_length_factor = float(sys.argv[6])
18+
step_length_relax = float(sys.argv[7])
19+
20+
Dimension = 2
21+
22+
PointType = itk.Point[itk.D, Dimension]
23+
PathInformationType = itk.SpeedFunctionPathInformation[PointType]
24+
25+
path_information = PathInformationType.New()
26+
with open(path_filename, 'r') as fp:
27+
for line in fp:
28+
line = line.replace("Path: ", "")
29+
line = line.replace("[", "").strip()
30+
points = line.split(']')[:-1]
31+
for index, point in enumerate(points):
32+
point_float = PointType()
33+
point_float[0] = float(point.split(',')[0])
34+
point_float[1] = float(point.split(',')[1])
35+
print(point_float)
36+
if index == 0:
37+
path_information.SetStartPoint(point_float)
38+
elif index == len(points)-1:
39+
path_information.SetEndPoint(point_float)
40+
else:
41+
path_information.AddWayPoint(point_float)
42+
43+
# path_information2 = PathInformationType.New()
44+
# start = PointType()
45+
# start[0] = 202.
46+
# start[1] = 369.
47+
# path_information2.SetStartPoint(start)
48+
# way1 = PointType()
49+
# way1[0] = 201
50+
# way1[1] = 335
51+
# path_information2.AddWayPoint(way1)
52+
# end = PointType()
53+
# end[0] = 165
54+
# end[1] = 326
55+
# path_information2.SetEndPoint(end)
56+
57+
speed_image = itk.imread(speed_filename, itk.F)
58+
59+
interpolator = itk.LinearInterpolateImageFunction.New(speed_image)
60+
interpolator.SetInputImage(speed_image)
61+
62+
cost_function = itk.SingleImageCostFunction[type(speed_image)].New(interpolator=interpolator)
63+
cost_function.SetInterpolator(interpolator)
64+
65+
spacing = list(speed_image.GetSpacing())
66+
min_spacing = min(spacing)
67+
68+
optimizer = itk.RegularStepGradientDescentOptimizer.New(
69+
number_of_iterations=number_of_iterations,
70+
maximum_step_length = 1.0*step_length_factor*min_spacing,
71+
minimum_step_length = 0.5*step_length_factor*min_spacing,
72+
relaxation_factor=step_length_relax)
73+
74+
path_filter = itk.SpeedFunctionToPathFilter.New(speed_image,
75+
cost_function=cost_function,
76+
optimizer=optimizer,
77+
termination_value=termination_value)
78+
path_filter.SetInput(speed_image)
79+
path_filter.AddPathInformation(path_information)
80+
path_filter.Update()
81+
82+
number_of_paths = path_filter.GetNumberOfOutputs()
83+
print('Number of paths: ' + str(number_of_paths))
84+
for ii in range(number_of_paths):
85+
path = path_filter.GetOutput(ii)
86+
number_of_vertices = path.GetVertexList().Size()
87+
print('Number of Path {0} vertices: {1}'.format(ii, number_of_vertices))
88+
assert(number_of_vertices != 0)
89+
90+
path_filter.GetNumberOfOutputs()
91+
vl = path.GetVertexList()
92+
path_filter.GetOptimizer()
93+
arrival = path_filter.GetCurrentArrivalFunction()
94+
print(np.where(itk.array_from_image(arrival) < 100.0))

0 commit comments

Comments
 (0)