-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmidvein_finder.py
195 lines (179 loc) · 9.9 KB
/
midvein_finder.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
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
import numpy as np
import networkx as nx
import utils
from skimage import measure
def limited_pairs_longest_shortest_path_length(graph, nodes_list):
max_length = 0
max_length_path = None
for i in range(len(nodes_list)):
for j in range(i+1, len(nodes_list)):
shortest_path = nx.dijkstra_path(graph, tuple(nodes_list[i]), tuple(nodes_list[j]))
shortest_path_length = nx.dijkstra_path_length(graph, tuple(nodes_list[i]), tuple(nodes_list[j]))
if shortest_path_length > max_length:
max_length = shortest_path_length
max_length_path = shortest_path
return max_length_path, max_length
def mask2contour(mask, approx_factor=0.0, longest_contour_only=False):
"""Find the outlines of the mask. If the approx_factor is not 0.0, it will return approximated
with the step length by approx_factor*total_length
"""
contours = measure.find_contours(~mask.astype(bool), 0)
if approx_factor != 0.0:
approx_contours = [measure.approximate_polygon(contour, .05 * utils.contour_length(contour)) for contour in contours]
if longest_contour_only:
return contours[np.argmax(list(map(len, contours)))]
else:
return contours
class SorghumLeafMeasure:
# TODO add a method to get edge_point_list
# TODO add a option to down sample the image
def __init__(self, mask, xyz_map, max_neibor_pixel=5, downsample=False):
self.mask = mask
self.xyz_map = xyz_map
self.max_neibor_pixel = max_neibor_pixel
self.leaf_edge_approx_factor = 0.05
self.leaf_graph = self.create_graph()
self._leaf_mask_for_edge_detect = np.pad(self.mask, 1, 'constant', constant_values=(0, 0))
self.leaf_edge = mask2contour(self._leaf_mask_for_edge_detect, approx_factor=0.0, longest_contour_only=True).astype(int) - 1
self.leaf_edge_approx = measure.approximate_polygon(self.leaf_edge, self.leaf_edge_approx_factor * utils.contour_length(self.leaf_edge))
self.leaf_len = None
# leaf_len_path is the midvein
self.leaf_len_path = None
self.leaf_width = None
self.leaf_width_path = None
def create_graph(self):
leaf_graph = nx.Graph()
contain_points_list = np.nonzero(self.mask)
contain_points_list = np.stack(contain_points_list, axis=1)
# add points to graph
# TODO could be optimized by add_nodes()
for point in contain_points_list:
leaf_graph.add_node(tuple(point))
# add edges to graph
for node in leaf_graph:
node_x = node[0]
node_y = node[1]
coord_0 = np.array([self.xyz_map[node_x, node_y, 0],
self.xyz_map[node_x, node_y, 1],
self.xyz_map[node_x, node_y, 2]])
neighbor_list = []
# TODO may be optimized by add only one side neighbor
for m in range(-self.max_neibor_pixel, self.max_neibor_pixel):
for n in range(-self.max_neibor_pixel, self.max_neibor_pixel):
neighbor_list.append([node[0]+m, node[1]+n])
neighbor_list.remove([node[0], node[1]])
for neighbor in neighbor_list:
if tuple(neighbor) not in leaf_graph:
continue
x = neighbor[0]
y = neighbor[1]
coord_1 = np.array([self.xyz_map[x, y, 0],
self.xyz_map[x, y, 1],
self.xyz_map[x, y, 2]])
weight = np.linalg.norm(coord_0 - coord_1)
leaf_graph.add_edge(node, tuple(neighbor), weight=weight)
return leaf_graph
def calc_leaf_length(self):
self.leaf_len_path, self.leaf_len = limited_pairs_longest_shortest_path_length(self.leaf_graph, self.leaf_edge_approx)
return self.leaf_len
def split_leaf_edge(self):
"""Using the two endpoints of the leaf to seperate the leaf edge as two part"""
endpoint_0 = self.leaf_len_path[0]
endpoint_1 = self.leaf_len_path[-1]
endpoint_0_idx = np.where((self.leaf_edge == endpoint_0).all(axis=1))[0][0]
endpoint_1_idx = np.where((self.leaf_edge == endpoint_1).all(axis=1))[0][0]
side_0_start = min(endpoint_0_idx, endpoint_1_idx)
side_0_end = max(endpoint_0_idx, endpoint_1_idx)
self.side_0 = self.leaf_edge[side_0_start:side_0_end]
self.side_1 = np.concatenate([self.leaf_edge[side_0_end:-1], self.leaf_edge[:side_0_start]])
return self.side_0, self.side_1
def calc_leaf_width(self, split_n=5):
# get points of n-section
split_length = self.leaf_len / split_n
next_point_position = split_length
prev_travelled = 0.0
n_section_points_list = []
travelled = 0.0
target_length = split_length
for i in range(1, len(self.leaf_len_path)):
# |----+|--------+-----+---
vector_length = self.leaf_graph[self.leaf_len_path[i-1]][self.leaf_len_path[i]]['weight']
travelled += vector_length
if travelled >= target_length:
n_section_points_list.append(self.leaf_len_path[i])
while travelled >= target_length:
target_length += split_length
# for i in range(1, len(self.leaf_len_path) - 1):
# vector_to_next = np.array(self.leaf_len_path[i]) - np.array(self.leaf_len_path[i-1])
# points_distance = np.linalg.norm(vector_to_next)
# current_travelled = prev_travelled + points_distance
# while current_travelled > next_point_position:
# length_ratio = (next_point_position - prev_travelled) / points_distance
# n_section_point = np.array(self.leaf_len_path[i-1]) + length_ratio * (vector_to_next)
# n_section_points_list.append(n_section_point.round().astype(int))
# next_point_position += split_length
# if len(n_section_points_list) >= split_n:
# break
# prev_travelled = current_travelled
self.leaf_width_mid_points = np.array(n_section_points_list)
# Throw outside n_section point
# n_section_points_list = list(filter (lambda x: tuple(x) in self.leaf_graph.node and ~any(np.equal(self.leaf_edge, x).all(1)), n_section_points_list))
aux_edges = []
for edge_point in self.leaf_edge:
aux_edges.append(('aux_point', tuple(edge_point)))
###############
aux_edges_side_0 = []
aux_edges_side_1 = []
side_0, side_1 = self.split_leaf_edge()
for edge_point in side_0:
aux_edges_side_0.append(('aux_point', tuple(edge_point)))
for edge_point in side_1:
aux_edges_side_1.append(('aux_point', tuple(edge_point)))
n_section_leaf_width_list = []
n_section_leaf_width_path_list = []
for n_section_point in n_section_points_list:
# find the nearest point from the point of n-section
"""
The nearest point can be find by two method:
1. norm
2. add an point that connected to all the edge points with same weight
then get the shortest path between that point and n_section_point
"""
temp_graph = self.leaf_graph.copy()
temp_graph.add_node('aux_point')
temp_graph.add_edges_from(aux_edges_side_0)
aux_shortest_path = nx.dijkstra_path(temp_graph, tuple(n_section_point), 'aux_point')
nearest_point_on_edge_0 = np.array(aux_shortest_path[-2])
temp_graph.remove_edges_from(aux_edges_side_0)
temp_graph.add_edges_from(aux_edges_side_1)
aux_shortest_path = nx.dijkstra_path(temp_graph, tuple(n_section_point), 'aux_point')
temp_graph.remove_edges_from(aux_edges_side_1)
nearest_point_on_edge_1 = np.array(aux_shortest_path[-2])
leaf_width_path = nx.dijkstra_path(self.leaf_graph,
tuple(nearest_point_on_edge_0),
tuple(nearest_point_on_edge_1))
leaf_width = nx.dijkstra_path_length(self.leaf_graph,
tuple(nearest_point_on_edge_0),
tuple(nearest_point_on_edge_1))
# remove vectors at same side
# point_cosine = np.dot(nearest_point_on_edge - n_section_point, (self.leaf_edge - n_section_point).T)/(np.linalg.norm(self.leaf_edge)* np.linalg.norm(n_section_point))
# self._point_dot_prod = point_cosine
# same_side_points = list(map(tuple, self.leaf_edge[np.where(point_cosine>0)]))
# self._same_side_points = np.array(same_side_points)
# temp_graph.remove_nodes_from(same_side_points)
# find the nearest point for another sid
# aux_other_side_shortest_path = nx.dijkstra_path(temp_graph, tuple(n_section_point), 'aux_point')
# nearest_point_on_edge_other_side = np.array(aux_other_side_shortest_path[-2])
# add two up
# leaf_width_path = nx.dijkstra_path(self.leaf_graph,
# tuple(nearest_point_on_edge),
# tuple(nearest_point_on_edge_other_side))
# leaf_width = nx.dijkstra_path_length(self.leaf_graph,
# tuple(nearest_point_on_edge),
# tuple(nearest_point_on_edge_other_side))
n_section_leaf_width_list.append(leaf_width)
n_section_leaf_width_path_list.append(leaf_width_path)
max_idx = np.argmax(n_section_leaf_width_list)
self.leaf_width = n_section_leaf_width_list[max_idx]
self.leaf_width_path = n_section_leaf_width_path_list[max_idx]
return self.leaf_width