@@ -73,110 +73,66 @@ void CLD::genCLD() {
7373 binaryThresholding (FDoG, result, this ->tau );
7474}
7575
76-
77-
7876/* *
79- * Private Functions
77+ * Flow-based DoG filtering
8078 */
8179void CLD::flowDoG (Mat & src, Mat & dst, const double sigma_m) {
8280 vector<double > gau_m;
8381 MakeGaussianVector (sigma_m, gau_m);
8482
85- vector<double > vt (2 , 0 );
86- double x, y, d_x, d_y;
87- double weight1, w_sum1, sum1;
88-
89- int i_x, i_y;
90- int x1, y1;
91- double val;
92-
93- int image_x = src.rows ;
94- int image_y = src.cols ;
95-
96- int half_l = gau_m.size () - 1 ;
97-
98- int flow_DOG_sign = 0 ;
99-
100- for (int i = 0 ; i < image_x; i++) {
101- for (int j = 0 ; j < image_y; j++) {
102- sum1 = w_sum1 = weight1 = 0.0 ;
103-
104- val = src.at <float >(i, j);
105- weight1 = gau_m[0 ];
106- sum1 = val * weight1;
107- w_sum1 += weight1;
108-
109- d_x = (double )i;
110- d_y = (double )j;
111- i_x = i;
112- i_y = j;
113- for (int k = 0 ; k < half_l; k++) {
114- vt[0 ] = etf.flowField .at <Vec3f>(i_x, i_y)[0 ];
115- vt[1 ] = etf.flowField .at <Vec3f>(i_x, i_y)[1 ];
116-
117- if (vt[0 ] == 0.0 && vt[1 ] == 0.0 ) break ;
118-
119- x = d_x;
120- y = d_y;
121-
122- if (x >(double )image_x - 1 || x < 0.0 || y >(double )image_y - 1 || y < 0.0 ) break ;
123-
124- x1 = min (max ((int )round (x), 0 ), image_x - 1 );
125- y1 = min (max ((int )round (y), 0 ), image_y - 1 );
126-
127- val = src.at <float >(x1, y1);
128-
129- weight1 = gau_m[k];
130-
131- sum1 += val * weight1;
132- w_sum1 += weight1;
133-
134- d_x += vt[0 ] * STEPSIZE;
135- d_y += vt[1 ] * STEPSIZE;
136-
137- i_x = round (d_x);
138- i_y = round (d_y);
139-
140- if (d_x < 0 || d_x > image_x - 1 || d_y < 0 || d_y > image_y - 1 ) break ;
83+ const int img_w = src.cols ;
84+ const int img_h = src.rows ;
85+ const int kernel_half = gau_m.size () - 1 ;
86+
87+ for (int y = 0 ; y < img_h; y++) {
88+ for (int x = 0 ; x < img_w; x++) {
89+ double gau_m_acc = -gau_m[0 ] * src.at <float >(y, x);
90+ double gau_m_weight_acc = -gau_m[0 ];
91+
92+ // Intergral alone ETF
93+ Point2f pos (x, y);
94+ for (int step = 0 ; step < kernel_half; step++) {
95+ Vec3f tmp = etf.flowField .at <Vec3f>((int )round (pos.y ), (int )round (pos.x ));
96+ Point2f direction = Point2f (tmp[1 ], tmp[0 ]);
97+
98+ if (direction.x == 0 && direction.y == 0 ) break ;
99+ if (pos.x > (double )img_w - 1 || pos.x < 0.0 || pos.y >(double )img_h - 1 || pos.y < 0.0 ) break ;
100+
101+ float value = src.at <float >((int )round (pos.y ), (int )round (pos.x ));
102+ float weight = gau_m[step];
103+
104+ gau_m_acc += value*weight;
105+ gau_m_weight_acc += weight;
106+
107+ // move alone ETF direction
108+ pos += direction;
109+
110+ if ((int )round (pos.x ) < 0 || (int )round (pos.x ) > img_w - 1 || (int )round (pos.y ) < 0 || (int )round (pos.y ) > img_h - 1 ) break ;
141111 }
142112
143- d_x = (double )i;
144- d_y = (double )j;
145- i_x = i;
146- i_y = j;
147- for (int k = 0 ; k < half_l; k++) {
148- vt[0 ] = -etf.flowField .at <Vec3f>(i_x, i_y)[0 ];
149- vt[1 ] = -etf.flowField .at <Vec3f>(i_x, i_y)[1 ];
150- if (vt[0 ] == 0.0 && vt[1 ] == 0.0 ) break ;
151-
152- x = d_x;
153- y = d_y;
113+ // Intergral alone inverse ETF
114+ pos = Point2f (x, y);
115+ for (int step = 0 ; step < kernel_half; step++) {
116+ Vec3f tmp = -etf.flowField .at <Vec3f>((int )round (pos.y ), (int )round (pos.x ));
117+ Point2f direction = Point2f (tmp[1 ], tmp[0 ]);
154118
155- if (x >(double )image_x - 1 || x < 0.0 || y >(double )image_y - 1 || y < 0.0 ) break ;
119+ if (direction.x == 0 && direction.y == 0 ) break ;
120+ if (pos.x >(double )img_w - 1 || pos.x < 0.0 || pos.y >(double )img_h - 1 || pos.y < 0.0 ) break ;
156121
157- x1 = min ( max (( int )round (x ), 0 ), image_x - 1 );
158- y1 = min ( max (( int ) round (y), 0 ), image_y - 1 ) ;
122+ float value = src. at < float >(( int )round (pos. y ), ( int ) round (pos. x ) );
123+ float weight = gau_m[step] ;
159124
160- val = src.at <float >(x1, y1);
125+ gau_m_acc += value*weight;
126+ gau_m_weight_acc += weight;
161127
162- weight1 = gau_m[k];
128+ // move alone ETF direction
129+ pos += direction;
163130
164- sum1 += val * weight1;
165- w_sum1 += weight1;
166-
167- d_x += vt[0 ] * STEPSIZE;
168- d_y += vt[1 ] * STEPSIZE;
169-
170- i_x = round (d_x);
171- i_y = round (d_y);
172-
173- if (d_x < 0 || d_x > image_x - 1 || d_y < 0 || d_y > image_y - 1 ) break ;
131+ if ((int )round (pos.x ) < 0 || (int )round (pos.x ) > img_w - 1 || (int )round (pos.y ) < 0 || (int )round (pos.y ) > img_h - 1 ) break ;
174132 }
175133
176- sum1 /= w_sum1;
177-
178-
179- dst.at <float >(i, j) = (sum1 > 0 ) ? 1.0 : 1.0 + tanh (sum1);
134+
135+ dst.at <float >(y, x) = (gau_m_acc / gau_m_weight_acc) > 0 ? 1.0 : 1 + tanh (gau_m_acc / gau_m_weight_acc);
180136 }
181137 }
182138
0 commit comments