1
+ ##################################################################################
2
+ # Author: Ricardo Sanchez Matilla
3
+
4
+ # Created Date: 2020/02/13
5
+ # Modified Date: 2020/02/28
6
+
7
+ # Centre for Intelligent Sensing, Queen Mary University of London, UK
8
+ #
9
+ ##################################################################################
10
+ # License
11
+ # This work is licensed under the Creative Commons Attribution-NonCommercial 4.0
12
+ # International License. To view a copy of this license, visit
13
+ # http://creativecommons.org/licenses/by-nc/4.0/ or send a letter to
14
+ # Creative Commons, PO Box 1866, Mountain View, CA 94042, USA.
15
+ ##################################################################################
16
+ #
17
+ # OpenCV
18
+ import cv2
19
+ #
20
+ # Numpy
21
+ import numpy as np
22
+ import numpy .ma as ma
23
+ from numpy import linalg as LA
24
+ from numpy .linalg import inv
25
+ #
26
+ import pickle
27
+ #
28
+ # ScyPy
29
+ from scipy .spatial .distance import cdist
30
+ from scipy .signal import find_peaks
31
+ #
32
+ import math
33
+ #
34
+ import copy
35
+ #
36
+ #
37
+ def getCentroid (mask ):
38
+ # Get the largest contour
39
+ contours , _ = cv2 .findContours (mask , cv2 .RETR_TREE , cv2 .CHAIN_APPROX_SIMPLE )
40
+ contour_sizes = [(cv2 .contourArea (contour ), contour ) for contour in contours ]
41
+ largest_contour = max (contour_sizes , key = lambda x : x [0 ])[1 ]
42
+
43
+ # Get centroid of the largest contour
44
+ M = cv2 .moments (largest_contour )
45
+
46
+ try :
47
+ centroid = np .array ((M ['m10' ]/ M ['m00' ], M ['m01' ]/ M ['m00' ]))
48
+ return centroid , largest_contour .squeeze ()
49
+ except :
50
+ print ('Centroid not found' )
51
+ return None , None
52
+
53
+ #
54
+ def triangulate (c1 , c2 , point1 , point2 , undistort = True ):
55
+
56
+ if (point1 .dtype != 'float64' ):
57
+ point1 = point1 .astype (np .float64 )
58
+
59
+ if (point2 .dtype != 'float64' ):
60
+ point2 = point2 .astype (np .float64 )
61
+
62
+ point3d = cv2 .triangulatePoints (c1 .extrinsic ['rgb' ]['projMatrix' ], c2 .extrinsic ['rgb' ]['projMatrix' ], point1 .reshape (2 ,1 ), point2 .reshape (2 ,1 )).transpose ()
63
+ for point in point3d :
64
+ point /= point [- 1 ]
65
+ return point3d .reshape (- 1 )
66
+
67
+ #
68
+ def get3D (c1 , c2 , mask1 , mask2 , glass , _img1 = None , _img2 = None , drawCentroid = False , drawDimensions = False ):
69
+
70
+ img1 = copy .deepcopy (_img1 )
71
+ img2 = copy .deepcopy (_img2 )
72
+
73
+ centr1 = getCentroid (c1 , mask1 )
74
+ centr2 = getCentroid (c2 , mask2 )
75
+ if centr1 is not None and centr2 is not None :
76
+ glass .centroid = triangulate (c1 , c2 , centr1 , centr2 )[:- 1 ].reshape (- 1 ,3 )
77
+
78
+ # Draw centroid
79
+ if drawCentroid :
80
+
81
+ # Draw 2D centroid of tracking mask
82
+ #cv2.circle(img1, tuple(centr1.astype(int)), 10, (0,128,0), -1)
83
+ #cv2.circle(img2, tuple(centr2.astype(int)), 10, (0,128,0), -1)
84
+
85
+ # Draw 3D centroid projected to image
86
+ point1 , _ = cv2 .projectPoints (glass .centroid , c1 .extrinsic ['rgb' ]['rvec' ], c1 .extrinsic ['rgb' ]['tvec' ], c1 .intrinsic ['rgb' ], c1 .distCoeffs )
87
+ point2 , _ = cv2 .projectPoints (glass .centroid , c2 .extrinsic ['rgb' ]['rvec' ], c2 .extrinsic ['rgb' ]['tvec' ], c2 .intrinsic ['rgb' ], c2 .distCoeffs )
88
+
89
+ point1 = point1 .squeeze ().astype (int )
90
+ point2 = point2 .squeeze ().astype (int )
91
+
92
+ cv2 .circle (img1 , tuple (point1 ), 6 , (128 ,0 ,0 ), - 1 )
93
+ cv2 .circle (img2 , tuple (point2 ), 6 , (128 ,0 ,0 ), - 1 )
94
+
95
+ # Draw height and width lines
96
+ if drawDimensions :
97
+
98
+ # Get top/bottom points
99
+ top = copy .deepcopy (glass .centroid )
100
+ bottom = copy .deepcopy (glass .centroid )
101
+ top [0 ,2 ] += glass .h / 2
102
+ bottom [0 ,2 ] -= glass .h / 2
103
+ topC1 , _ = cv2 .projectPoints (top , c1 .extrinsic ['rgb' ]['rvec' ], c1 .extrinsic ['rgb' ]['tvec' ], c1 .intrinsic ['rgb' ], c1 .distCoeffs )
104
+ bottomC1 , _ = cv2 .projectPoints (bottom , c1 .extrinsic ['rgb' ]['rvec' ], c1 .extrinsic ['rgb' ]['tvec' ], c1 .intrinsic ['rgb' ], c1 .distCoeffs )
105
+ topC2 , _ = cv2 .projectPoints (top , c2 .extrinsic ['rgb' ]['rvec' ], c1 .extrinsic ['rgb' ]['tvec' ], c2 .intrinsic ['rgb' ], c2 .distCoeffs )
106
+ bottomC2 , _ = cv2 .projectPoints (bottom , c2 .extrinsic ['rgb' ]['rvec' ], c2 .extrinsic ['rgb' ]['tvec' ], c2 .intrinsic ['rgb' ], c2 .distCoeffs )
107
+ topC1 = topC1 .squeeze ().astype (int )
108
+ bottomC1 = bottomC1 .squeeze ().astype (int )
109
+ topC2 = topC2 .squeeze ().astype (int )
110
+ bottomC2 = bottomC2 .squeeze ().astype (int )
111
+
112
+ # Get rigth/left points
113
+ right = copy .deepcopy (glass .centroid )
114
+ left = copy .deepcopy (glass .centroid )
115
+ right [0 ,0 ] += glass .w / 2
116
+ left [0 ,0 ] -= glass .w / 2
117
+ rightC1 , _ = cv2 .projectPoints (right , c1 .extrinsic ['rgb' ]['rvec' ], c1 .extrinsic ['rgb' ]['tvec' ], c1 .intrinsic ['rgb' ], c1 .distCoeffs )
118
+ leftC1 , _ = cv2 .projectPoints (left , c1 .extrinsic ['rgb' ]['rvec' ], c1 .extrinsic ['rgb' ]['tvec' ], c1 .intrinsic ['rgb' ], c1 .distCoeffs )
119
+ rightC2 , _ = cv2 .projectPoints (right , c2 .extrinsic ['rgb' ]['rvec' ], c2 .extrinsic ['rgb' ]['tvec' ], c2 .intrinsic ['rgb' ], c2 .distCoeffs )
120
+ leftC2 , _ = cv2 .projectPoints (left , c2 .extrinsic ['rgb' ]['rvec' ], c2 .extrinsic ['rgb' ]['tvec' ], c2 .intrinsic ['rgb' ], c2 .distCoeffs )
121
+ rightC1 = rightC1 .squeeze ().astype (int )
122
+ leftC1 = leftC1 .squeeze ().astype (int )
123
+ rightC2 = rightC2 .squeeze ().astype (int )
124
+ leftC2 = leftC2 .squeeze ().astype (int )
125
+
126
+ cv2 .line (img1 , tuple (topC1 ), tuple (bottomC1 ), (128 ,0 ,0 ), 2 )
127
+ cv2 .line (img1 , tuple (rightC1 ), tuple (leftC1 ), (128 ,0 ,0 ), 2 )
128
+ cv2 .line (img2 , tuple (topC2 ), tuple (bottomC2 ), (128 ,0 ,0 ), 2 )
129
+ cv2 .line (img2 , tuple (rightC2 ), tuple (leftC2 ), (128 ,0 ,0 ), 2 )
130
+
131
+
132
+ return glass , img1 , img2
133
+
134
+ #
135
+ def pointsOnContour (p2d_c1 , p2d_c2 , _contour1 , _contour2 , _c1 , _c2 , draw = True ):
136
+ contour1 = copy .deepcopy (_contour1 )
137
+ contour2 = copy .deepcopy (_contour2 )
138
+
139
+ c1 = copy .deepcopy (_c1 )
140
+ c2 = copy .deepcopy (_c2 )
141
+
142
+ # Closer re-projected point to contour [ < 5 pixels]
143
+ distances_c1 = cdist (p2d_c1 , contour1 )
144
+ distances_c2 = cdist (p2d_c2 , contour2 )
145
+
146
+ avg_dist_c1 = np .mean (distances_c1 )
147
+ avg_std_c1 = np .std (distances_c1 )
148
+
149
+ avg_dist_c2 = np .mean (distances_c2 )
150
+ avg_std_c2 = np .std (distances_c2 )
151
+ print ('C1: avgDistance {:.2f} +- {:.2f}' .format (avg_dist_c1 , avg_std_c1 ))
152
+ print ('C2: avgDistance {:.2f} +- {:.2f}' .format (avg_dist_c2 , avg_std_c2 ))
153
+
154
+ #
155
+ def getObjectDimensions (_c1 , _c2 , centroid , draw = False ):
156
+
157
+ c1 = copy .deepcopy (_c1 )
158
+ c2 = copy .deepcopy (_c2 )
159
+
160
+ # Radiuses
161
+ step = 0.001 #meters
162
+ minDiameter = 0.005 #meters
163
+ maxDiameter = 0.15 #meters
164
+ radiuses = np .linspace (maxDiameter / 2 , minDiameter / 2 , num = int ((maxDiameter - minDiameter )/ step ))
165
+
166
+ angularStep = 18 #degrees
167
+ angles = np .linspace (0. , 359. , num = int ((359. )/ angularStep ))
168
+
169
+ # Heights
170
+ step = 0.001 #meters
171
+ minHeight = - 0.1 #meters
172
+ maxHeight = 0.4 #meters
173
+ estRadius = []
174
+ converged = []
175
+
176
+ heights = np .linspace (minHeight , maxHeight , num = int ((maxHeight - minHeight )/ step ))
177
+ for height in heights :
178
+ for rad in radiuses :
179
+ seg1 = copy .deepcopy (c1 ['seg' ])
180
+ seg2 = copy .deepcopy (c2 ['seg' ])
181
+
182
+ # Sample 3D circunference
183
+ p3d = []
184
+ for angle_d in angles :
185
+ angle = math .radians (angle_d )
186
+ p3d .append (np .array ((centroid [0 ]+ (rad * math .cos (angle )), centroid [1 ]+ (rad * math .sin (angle )), height )).reshape (1 ,3 ))
187
+ p3d = np .array (p3d )
188
+
189
+ # Reproject to C1
190
+ p2d_c1 , _ = cv2 .projectPoints (p3d , c1 ['extrinsic' ]['rvec' ], c1 ['extrinsic' ]['tvec' ], c1 ['intrinsic' ], np .array ([0. ,0. ,0. ,0. ,0. ]))
191
+ p2d_c1 = p2d_c1 .squeeze ().astype (int )
192
+
193
+ # Reproject to C2
194
+ p2d_c2 , _ = cv2 .projectPoints (p3d , c2 ['extrinsic' ]['rvec' ], c2 ['extrinsic' ]['tvec' ], c2 ['intrinsic' ], np .array ([0. ,0. ,0. ,0. ,0. ]))
195
+ p2d_c2 = p2d_c2 .squeeze ().astype (int )
196
+
197
+ # Check if imaged points are in the segmentation
198
+ areIn_c1 = seg1 [p2d_c1 [:,1 ], p2d_c1 [:,0 ]]
199
+ areIn_c2 = seg2 [p2d_c2 [:,1 ], p2d_c2 [:,0 ]]
200
+
201
+ if (np .count_nonzero (areIn_c1 ) == areIn_c1 .shape [0 ]) and (np .count_nonzero (areIn_c2 ) == areIn_c2 .shape [0 ]):
202
+ estRadius .append (rad )
203
+ converged .append (True )
204
+ break
205
+ if rad == minDiameter / 2 :
206
+ estRadius .append (rad )
207
+ converged .append (False )
208
+ break
209
+
210
+
211
+ estRadius = np .array (estRadius )
212
+ converged = np .array (converged )
213
+ estHeights = heights [converged ]
214
+
215
+ width = np .max (estRadius ) * 2. * 1000
216
+ height = (estHeights [- 1 ] - estHeights [0 ]) * 1000
217
+
218
+ # Draw final dimensions
219
+ if draw :
220
+ img1 = copy .deepcopy (c1 ['rgb' ])
221
+ img2 = copy .deepcopy (c2 ['rgb' ])
222
+
223
+ for i , rad in enumerate (estRadius ):
224
+
225
+ p3d = []
226
+ for angle_d in angles :
227
+ angle = math .radians (angle_d )
228
+ p3d .append (np .array ((centroid [0 ]+ (rad * math .cos (angle )), centroid [1 ]+ (rad * math .sin (angle )), heights [i ])).reshape (1 ,3 ))
229
+ p3d = np .array (p3d )
230
+
231
+ # Reproject to C1
232
+ p2d_c1 , _ = cv2 .projectPoints (p3d , c1 ['extrinsic' ]['rvec' ], c1 ['extrinsic' ]['tvec' ], c1 ['intrinsic' ], np .array ([0. ,0. ,0. ,0. ,0. ]))
233
+ p2d_c1 = p2d_c1 .squeeze ().astype (int )
234
+
235
+ # Reproject to C2
236
+ p2d_c2 , _ = cv2 .projectPoints (p3d , c2 ['extrinsic' ]['rvec' ], c2 ['extrinsic' ]['tvec' ], c2 ['intrinsic' ], np .array ([0. ,0. ,0. ,0. ,0. ]))
237
+ p2d_c2 = p2d_c2 .squeeze ().astype (int )
238
+
239
+ # Check if imaged points are in the segmentation
240
+ areIn_c1 = seg1 [p2d_c1 [:,1 ], p2d_c1 [:,0 ]]
241
+ areIn_c2 = seg2 [p2d_c2 [:,1 ], p2d_c2 [:,0 ]]
242
+
243
+ for p , isIn in zip (p2d_c1 , areIn_c1 ):
244
+ if isIn :
245
+ cv2 .circle (img1 , (int (p [0 ]), int (p [1 ])), 2 , (0 ,255 ,0 ), - 1 )
246
+ else :
247
+ cv2 .circle (img1 , (int (p [0 ]), int (p [1 ])), 2 , (0 ,0 ,255 ), - 1 )
248
+
249
+ for p , isIn in zip (p2d_c2 , areIn_c2 ):
250
+ if isIn :
251
+ cv2 .circle (img2 , (int (p [0 ]), int (p [1 ])), 2 , (0 ,255 ,0 ), - 1 )
252
+ else :
253
+ cv2 .circle (img2 , (int (p [0 ]), int (p [1 ])), 2 , (0 ,0 ,255 ), - 1 )
254
+
255
+ return height , width , np .concatenate ((img1 ,img2 ), axis = 1 )
0 commit comments