vbl_local_minima.hxx
Go to the documentation of this file.
1 #ifndef VBL_LOCAL_MINIMA_TXX_
2 #define VBL_LOCAL_MINIMA_TXX_
3 #include <iostream>
4 #include <limits>
5 #include "vbl_local_minima.h"
6 //
7 #ifdef _MSC_VER
8 # include <vcl_msvc_warnings.h>
9 #endif
10 #include <cassert>
11 
12 template <class T>
13 bool local_minima(vbl_array_1d<T> const& in, vbl_array_1d<T>& minima, T thresh)
14 {
15  const unsigned int n = (unsigned int)(in.size());
16  assert(minima.size()==n);
17  //special cases
18  // minimum is not defined for n<3
19  if (n<3)
20  return false;
21  bool minima_found = false;
22  for (unsigned int i=0; i<n; ++i)
23  minima[i] = T(0);
24  //the general case
25  for (unsigned int c=1; c<n-1; ++c) {
26  T dm = in[c-1]-in[c], dp = in[c+1]-in[c];
27  if (dm>thresh && dp > thresh) {
28  T dmin = dm;
29  if (dp<dmin)
30  dmin = dp;
31  minima[c]=dmin;
32  minima_found = true;
33  }
34  }
35  // Check the ends of the array for minima
36  // left end
37  if ((in[1]-in[0])>thresh) {
38  minima[0]=in[1]-in[0];
39  minima_found = true;
40  }
41  // right end
42  if ((in[n-2]-in[n-1])>thresh) {
43  minima[n-1]=in[n-2]-in[n-1];
44  minima_found = true;
45  }
46 
47  return minima_found;
48 }
49 
50 template <class T>
51 bool local_minima(vbl_array_2d<T> const& in, vbl_array_2d<T>& minima, T thresh)
52 {
53  const unsigned int nr = (unsigned int)(in.rows()), nc = (unsigned int)(in.cols());
54  assert(nr==minima.rows() && nc==minima.cols());
55  //special case
56  // actually a 1-d array or null
57  if (nr<2||nc<2)
58  return false;
59  bool minima_found = false;
60  T mval = std::numeric_limits<T>::max();
61  for (unsigned r = 0; r<nr; ++r)
62  for (unsigned c = 0; c<nc; ++c)
63  minima[r][c] = T(0);
64  T ul, um, ur;
65  T lf,/*m,*/ri;
66  T ll, lm, lr;
67  T dmin;
68  //general case
69  if (nr>2&&nc>2)
70  for (unsigned r = 1; r<nr-1; ++r)
71  for (unsigned c = 1; c<nc-1; ++c) {
72  // xxxxxxxxxxxxxxxxxxxxxxx
73  // xxxx ul um ur xxxxx
74  // xxxx lf m ri xxxxx
75  // xxxx ll lm lr xxxxx
76  // xxxxxxxxxxxxxxxxxxxxxxx
77  ul = in[r-1][c-1] - in[r][c];
78  um = in[r-1][c] - in[r][c];
79  ur = in[r-1][c] - in[r][c];
80  lf = in[r][c-1] - in[r][c];
81  ri = in[r][c+1] - in[r][c];
82  ll = in[r+1][c-1] - in[r][c];
83  lm = in[r+1][c] - in[r][c];
84  lr = in[r+1][c+1] - in[r][c];
85  dmin = mval;
86  if (ul<=thresh) continue; if (ul<dmin) dmin = ul;
87  if (um<=thresh) continue; if (um<dmin) dmin = um;
88  if (ur<=thresh) continue; if (ur<dmin) dmin = ur;
89  if (lf<=thresh) continue; if (lf<dmin) dmin = lf;
90  if (ri<=thresh) continue; if (ri<dmin) dmin = ri;
91  if (ll<=thresh) continue; if (ll<dmin) dmin = ll;
92  if (lm<=thresh) continue; if (lm<dmin) dmin = lm;
93  if (lr<=thresh) continue; if (lr<dmin) dmin = lr;
94  if (dmin>thresh) {
95  minima[r][c] = dmin;
96  minima_found = true;
97  }
98  }
99 
100  // special cases at the borders
101  if (nc>2) {
102  for (unsigned c = 1; c<nc-1; ++c) {
103  // first row case
104  // xxxx lf m ri xxxxx
105  // xxxx ll lm lr xxxxx
106  // xxxxxxxxxxxxxxxxxxxxxxx
107  lf = in[0][c-1] - in[0][c];
108  ri = in[0][c+1] - in[0][c];
109  ll = in[1][c-1] - in[0][c];
110  lm = in[1][c] - in[0][c];
111  lr = in[1][c+1] - in[0][c];
112  dmin = mval;
113  if (lf<=thresh) continue; if (lf<dmin) dmin = lf;
114  if (ri<=thresh) continue; if (ri<dmin) dmin = ri;
115  if (ll<=thresh) continue; if (ll<dmin) dmin = ll;
116  if (lm<=thresh) continue; if (lm<dmin) dmin = lm;
117  if (lr<=thresh) continue; if (lr<dmin) dmin = lr;
118  if (dmin>thresh) {
119  minima[0][c] = dmin;
120  minima_found = true;
121  }
122  }
123  for (unsigned c = 1; c<nc-1; ++c) {
124  // last row case
125  // xxxxxxxxxxxxxxxxxxxxxxx
126  // xxxx ul um ur xxxxx
127  // xxxx lf m ri xxxxx
128  ul = in[nr-2][c-1] - in[nr-1][c];
129  um = in[nr-2][c] - in[nr-1][c];
130  ur = in[nr-2][c+1] - in[nr-1][c];
131  lf = in[nr-1][c-1] - in[nr-1][c];
132  ri = in[nr-1][c+1] - in[nr-1][c];
133  dmin = mval;
134  if (ul<=thresh) continue; if (ul<dmin) dmin = ul;
135  if (um<=thresh) continue; if (um<dmin) dmin = um;
136  if (ur<=thresh) continue; if (ur<dmin) dmin = ur;
137  if (lf<=thresh) continue; if (lf<dmin) dmin = lf;
138  if (ri<=thresh) continue; if (ri<dmin) dmin = ri;
139  if (dmin>thresh) {
140  minima[nr-1][c] = dmin;
141  minima_found = true;
142  }
143  }
144  }
145  if (nr>2) {
146  //first column case
147  for (unsigned r = 1; r<nr-1; ++r) {
148  // um ur xxxxx
149  // m ri xxxxx
150  // lm lr xxxxx
151  um = in[r-1][0] - in[r][0];
152  ur = in[r-1][1] - in[r][0];
153  ri = in[r][1] - in[r][0];
154  lm = in[r+1][0] - in[r][0];
155  lr = in[r+1][1] - in[r][0];
156  dmin = mval;
157  if (um<=thresh) continue; if (um<dmin) dmin = um;
158  if (ur<=thresh) continue; if (ur<dmin) dmin = ur;
159  if (ri<=thresh) continue; if (ri<dmin) dmin = ri;
160  if (lm<=thresh) continue; if (lm<dmin) dmin = lm;
161  if (lr<=thresh) continue; if (lr<dmin) dmin = lr;
162  if (dmin>thresh) {
163  minima[r][0] = dmin;
164  minima_found = true;
165  }
166  }
167  //last column case
168  for (unsigned r = 1; r<nr-1; ++r) {
169  // xxxxx ul um
170  // xxxxx lf m
171  // xxxxx ll lm
172  ul = in[r-1][nc-2] - in[r][nc-1];
173  um = in[r-1][nc-1] - in[r][nc-1];
174  lf = in[r][nc-2] - in[r][nc-1];
175  ll = in[r+1][nc-2] - in[r][nc-1];
176  lm = in[r+1][nc-1] - in[r][nc-1];
177  dmin = mval;
178  if (ul<=thresh) continue; if (ul<dmin) dmin = ul;
179  if (um<=thresh) continue; if (um<dmin) dmin = um;
180  if (lf<=thresh) continue; if (lf<dmin) dmin = lf;
181  if (ll<=thresh) continue; if (ll<dmin) dmin = ll;
182  if (lm<=thresh) continue; if (lm<dmin) dmin = lm;
183  if (dmin>thresh) {
184  minima[r][nc-1] = dmin;
185  minima_found = true;
186  }
187  }
188  }
189  // check the corners for minima
190  // upper left corner
191  // m ri xxxxx
192  // lm lr xxxxx
193  bool fail = false;
194  ri = in[0][1] - in[0][0];
195  lm = in[1][0] - in[0][0];
196  lr = in[1][1] - in[0][0];
197  dmin = mval;
198  if (ri<=thresh) fail = true; if (ri<dmin) dmin = ri;
199  if (lm<=thresh) fail = true; if (lm<dmin) dmin = lm;
200  if (lr<=thresh) fail = true; if (lr<dmin) dmin = lr;
201  if (!fail) {
202  minima[0][0] = dmin;
203  minima_found = true;
204  }
205  // upper right corner
206  // xxxxx lf m
207  // xxxxx ll lm
208  fail = false;
209  lf = in[0][nc-2] - in[0][nc-1];
210  lm = in[1][nc-1] - in[0][nc-1];
211  ll = in[1][nc-2] - in[0][nc-1];
212  dmin = mval;
213  if (lf<=thresh) fail = true; if (lf<dmin) dmin = lf;
214  if (lm<=thresh) fail = true; if (lm<dmin) dmin = lm;
215  if (ll<=thresh) fail = true; if (ll<dmin) dmin = ll;
216  if (!fail) {
217  minima[0][nc-1] = dmin;
218  minima_found = true;
219  }
220  // lower right corner
221  // xxxxx ul um
222  // xxxxx lf m
223  fail = false;
224  ul = in[nr-2][nc-2] - in[nr-1][nc-1];
225  um = in[nr-2][nc-1] - in[nr-1][nc-1];
226  lf = in[nr-1][nc-2] - in[nr-1][nc-1];
227  dmin = mval;
228  if (ul<=thresh) fail = true; if (ul<dmin) dmin = ul;
229  if (um<=thresh) fail = true; if (um<dmin) dmin = um;
230  if (lf<=thresh) fail = true; if (lf<dmin) dmin = lf;
231  if (!fail) {
232  minima[nr-1][nc-1] = dmin;
233  minima_found = true;
234  }
235  // lower left corner
236  // xxxxx um ur
237  // xxxxx m ri
238  fail = false;
239  ur = in[nr-2][1] - in[nr-1][0];
240  um = in[nr-2][0] - in[nr-1][0];
241  ri = in[nr-1][1] - in[nr-1][0];
242  dmin = mval;
243  if (ur<=thresh) fail = true; if (ur<dmin) dmin = ur;
244  if (um<=thresh) fail = true; if (um<dmin) dmin = um;
245  if (ri<=thresh) fail = true; if (ri<dmin) dmin = ri;
246  if (!fail) {
247  minima[nr-1][0] = dmin;
248  minima_found = true;
249  }
250  return minima_found;
251 }
252 
253 template <class T>
254 bool local_minima(vbl_array_3d<T> const& in, vbl_array_3d<T>& minima, T thresh)
255 {
256  const unsigned int n1=(unsigned int)(in.get_row1_count()),
257  n2=(unsigned int)(in.get_row2_count()),
258  n3=(unsigned int)(in.get_row3_count());
259  assert(n3==minima.get_row3_count() &&
260  n2==minima.get_row2_count() &&
261  n1==minima.get_row1_count() );
262  //special case
263  // actually a 2-d array or null
264  if (n3<2||n2<2||n1<2)
265  return false;
266  bool minima_found = false;
267  unsigned int x3=0, x2=0, x1=0;
268  for (x3 = 0; x3<n3; ++x3)
269  for (x2 = 0; x2<n2; ++x2)
270  for (x1 = 0; x1<n1; ++x1)
271  minima[x1][x2][x3] = T(0);
272  T v, d, mind;
273  bool fail;
274  const T mval = std::numeric_limits<T>::max();
275  // general case
276  if (n3>2||n2>2||n1>2)
277  for (x3 = 1; x3<n3-1; ++x3)
278  for (x2 = 1; x2<n2-1; ++x2)
279  for (x1 = 1; x1<n1-1; ++x1) {
280  mind = mval;
281  fail = false;
282  v = in[x1][x2][x3];
283  for (int k3 = -1; k3<=1;++k3)
284  for (int k2 = -1; k2<=1;++k2)
285  for (int k1 = -1; k1<=1;++k1)
286  if (k1!=0||k2!=0||k3!=0) {
287  d = in[x1+k1][x2+k2][x3+k3]-v;
288  if (d<=thresh) {fail = true; break;}
289  if (d<mind) mind = d;
290  }
291  if (!fail) { // local min found
292  minima[x1][x2][x3]=mind;
293  minima_found = true;
294  }
295  }
296  //face border cases six in total
297  if (n3>2&&n2>2) { // top and bottom faces of the array vary x3, x2
298  //top face
299  x1 = 0;
300  for (x3 = 1; x3<n3-1; ++x3)
301  for (x2 = 1; x2<n2-1; ++x2) {
302  mind = mval; fail = false;
303  v = in[x1][x2][x3];
304  for (int k3 = -1; k3<=1;++k3)
305  for (int k2 = -1; k2<=1;++k2)
306  for (int k1 = 0; k1<=1;++k1)
307  if (k3!=0||k2!=0||k1!=0) {
308  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
309  if (d<mind) mind = d;
310  }
311  if (!fail) { // local min found
312  minima[x1][x2][x3]=mind;
313  minima_found = true;
314  }
315  }
316  //bottom face
317  x1 = n1-1;
318  for (x3 = 1; x3<n3-1; ++x3)
319  for (x2 = 1; x2<n2-1; ++x2) {
320  mind = mval; fail = false;
321  v = in[x1][x2][x3];
322  for (int k3 = -1; k3<=1;++k3)
323  for (int k2 = -1; k2<=1;++k2)
324  for (int k1 = -1; k1<=0;++k1)
325  if (k3!=0||k2!=0||k1!=0) {
326  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
327  if (d<mind) mind = d;
328  }
329  if (!fail) { // local min found
330  minima[x1][x2][x3]=mind;
331  minima_found = true;
332  }
333  }
334  }
335  if (n3>2&&n1>2) { // front and back faces of the array vary x3, x1
336  //front face
337  x2 = 0;
338  for (x3 = 1; x3<n3-1; ++x3)
339  for (x1 = 1; x1<n1-1; ++x1) {
340  mind = mval; fail = false;
341  v = in[x1][x2][x3];
342  for (int k3 = -1; k3<=1;++k3)
343  for (int k2 = 0; k2<=1;++k2)
344  for (int k1 = -1; k1<=1;++k1)
345  if (k3!=0||k2!=0||k1!=0) {
346  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
347  if (d<mind) mind = d;
348  }
349  if (!fail) { // local min found
350  minima[x1][x2][x3]=mind;
351  minima_found = true;
352  }
353  }
354  //back face
355  x2 = n2-1;
356  for (x3 = 1; x3<n3-1; ++x3)
357  for (x1 = 1; x1<n1-1; ++x1) {
358  mind = mval; fail = false;
359  v = in[x1][x2][x3];
360  for (int k3 = -1; k3<=1;++k3)
361  for (int k2 = -1; k2<=0;++k2)
362  for (int k1 = -1; k1<=1;++k1)
363  if (k3!=0||k2!=0||k1!=0) {
364  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
365  if (d<mind) mind = d;
366  }
367  if (!fail) { // local min found
368  minima[x1][x2][x3]=mind;
369  minima_found = true;
370  }
371  }
372  }
373 
374  if (n2>2&&n1>2) { // left and right faces of the array vary x2, x1
375  //left face
376  x3 = 0;
377  for (x2 = 1; x2<n2-1; ++x2)
378  for (x1 = 1; x1<n1-1; ++x1) {
379  mind = mval; fail = false;
380  v = in[x1][x2][x3];
381  for (int k3 = 0; k3<=1;++k3)
382  for (int k2 = -1; k2<=1;++k2)
383  for (int k1 = -1; k1<=1;++k1)
384  if (k3!=0||k2!=0||k1!=0) {
385  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
386  if (d<mind) mind = d;
387  }
388  if (!fail) { // local min found
389  minima[x1][x2][x3]=mind;
390  minima_found = true;
391  }
392  }
393  //right face
394  x3 = n3-1;
395  for (x2 = 1; x2<n2-1; ++x2)
396  for (x1 = 1; x1<n1-1; ++x1) {
397  mind = mval; fail = false;
398  v = in[x1][x2][x3];
399  for (int k3 = -1; k3<=0;++k3)
400  for (int k2 = -1; k2<=1;++k2)
401  for (int k1 = -1; k1<=1;++k1)
402  if (k3!=0||k2!=0||k1!=0) {
403  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
404  if (d<mind) mind = d;
405  }
406  if (!fail) { // local min found
407  minima[x1][x2][x3]=mind;
408  minima_found = true;
409  }
410  }
411  }
412 
413  //edge border cases, 12 in total
414  if (x1>2)
415  {
416  //edges along x1
417  x3 = 0, x2 = 0; //vary x1
418  for (x1 = 1; x1<n1-1; ++x1) {
419  mind = mval; fail = false;
420  v = in[x1][x2][x3];
421  for (int k3 = 0; k3<=1;++k3)
422  for (int k2 = 0; k2<=1;++k2)
423  for (int k1 = -1; k1<=1;++k1)
424  if (k3!=0||k2!=0||k1!=0) {
425  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
426  if (d<mind) mind = d;
427  }
428  if (!fail) { // local min found
429  minima[x1][x2][x3]=mind;
430  minima_found = true;
431  }
432  }
433  x3 = n3-1; x2 = 0; //vary x1
434  for (x1 = 1; x1<n1-1; ++x1) {
435  mind = mval; fail = false;
436  v = in[x1][x2][x3];
437  for (int k3 = 0; k3<=1;++k3)
438  for (int k2 = 0; k2<=1;++k2)
439  for (int k1 = -1; k1<=1;++k1)
440  if (k3!=0||k2!=0||k1!=0) {
441  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
442  if (d<mind) mind = d;
443  }
444  if (!fail) { // local min found
445  minima[x1][x2][x3]=mind;
446  minima_found = true;
447  }
448  }
449  x3 = 0; x2 = n2-1; //vary x1
450  for (x1 = 1; x1<n1-1; ++x1) {
451  mind = mval; fail = false;
452  v = in[x1][x2][x3];
453  for (int k3 = 0; k3<=1;++k3)
454  for (int k2 = -1; k2<=0;++k2)
455  for (int k1 = -1; k1<=1;++k1)
456  if (k3!=0||k2!=0||k1!=0) {
457  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
458  if (d<mind) mind = d;
459  }
460  if (!fail) { // local min found
461  minima[x1][x2][x3]=mind;
462  minima_found = true;
463  }
464  }
465  x3 = n3-1; x2 = n2-1; //vary x1
466  for (x1 = 1; x1<n1-1; ++x1) {
467  mind = mval; fail = false;
468  v = in[x1][x2][x3];
469  for (int k3 = -1; k3<=0;++k3)
470  for (int k2 = -1; k2<=0;++k2)
471  for (int k1 = -1; k1<=1;++k1)
472  if (k3!=0||k2!=0||k1!=0) {
473  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
474  if (d<mind) mind = d;
475  }
476  if (!fail) { // local min found
477  minima[x1][x2][x3]=mind;
478  minima_found = true;
479  }
480  }
481  }
482  if (x2>2)
483  {
484  // edges along x2
485  x3 = 0, x1 = 0; //vary x2
486  for (x2 = 1; x2<n2-1; ++x2) {
487  mind = mval; fail = false;
488  v = in[x1][x2][x3];
489  for (int k3 = 0; k3<=1;++k3)
490  for (int k2 = -1; k2<=1;++k2)
491  for (int k1 = 0; k1<=1;++k1)
492  if (k3!=0||k2!=0||k1!=0) {
493  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
494  if (d<mind) mind = d;
495  }
496  if (!fail) { // local min found
497  minima[x1][x2][x3]=mind;
498  minima_found = true;
499  }
500  }
501  x3 = n3-1; x1 = 0; //vary x2
502  for (x2 = 1; x2<n2-1; ++x2) {
503  mind = mval; fail = false;
504  v = in[x1][x2][x3];
505  for (int k3 = -1; k3<=0;++k3)
506  for (int k2 = -1; k2<=1;++k2)
507  for (int k1 = 0; k1<=1;++k1)
508  if (k3!=0||k2!=0||k1!=0) {
509  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
510  if (d<mind) mind = d;
511  }
512  if (!fail) { // local min found
513  minima[x1][x2][x3]=mind;
514  minima_found = true;
515  }
516  }
517  x3 = 0; x1 = n1-1; //vary x2
518  for (x2 = 1; x2<n2-1; ++x2) {
519  mind = mval; fail = false;
520  v = in[x1][x2][x3];
521  for (int k3 = 0; k3<=1;++k3)
522  for (int k2 = -1; k2<=1;++k2)
523  for (int k1 = -1; k1<=0;++k1)
524  if (k3!=0||k2!=0||k1!=0) {
525  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
526  if (d<mind) mind = d;
527  }
528  if (!fail) { // local min found
529  minima[x1][x2][x3]=mind;
530  minima_found = true;
531  }
532  }
533  x3 = n3-1; x1 = n1-1; //vary x2
534  for (x2 = 1; x2<n2-1; ++x2) {
535  mind = mval; fail = false;
536  v = in[x1][x2][x3];
537  for (int k3 = -1; k3<=0;++k3)
538  for (int k2 = -1; k2<=1;++k2)
539  for (int k1 = -1; k1<=0;++k1)
540  if (k3!=0||k2!=0||k1!=0) {
541  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
542  if (d<mind) mind = d;
543  }
544  if (!fail) { // local min found
545  minima[x1][x2][x3]=mind;
546  minima_found = true;
547  }
548  }
549  }
550  if (x3>2)
551  {
552  // edges along x3
553  x2 = 0, x1 = 0; //vary x3
554  for (x3 = 1; x3<n3-1; ++x3) {
555  mind = mval; fail = false;
556  v = in[x1][x2][x3];
557  for (int k3 = -1; k3<=1;++k3)
558  for (int k2 = 0; k2<=1;++k2)
559  for (int k1 = 0; k1<=1;++k1)
560  if (k3!=0||k2!=0||k1!=0) {
561  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
562  if (d<mind) mind = d;
563  }
564  if (!fail) { // local min found
565  minima[x1][x2][x3]=mind;
566  minima_found = true;
567  }
568  }
569  x2 = n2-1; x1 = 0; //vary x3
570  for (x3 = 1; x3<n3-1; ++x3) {
571  mind = mval; fail = false;
572  v = in[x1][x2][x3];
573  for (int k3 = -1; k3<=1;++k3)
574  for (int k2 = -1; k2<=0;++k2)
575  for (int k1 = 0; k1<=1;++k1)
576  if (k3!=0||k2!=0||k1!=0) {
577  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
578  if (d<mind) mind = d;
579  }
580  if (!fail) { // local min found
581  minima[x1][x2][x3]=mind;
582  minima_found = true;
583  }
584  }
585  x2 = 0; x1 = n1-1; //vary x3
586  for (x3 = 1; x3<n3-1; ++x3) {
587  mind = mval; fail = false;
588  v = in[x1][x2][x3];
589  for (int k3 = -1; k3<=1;++k3)
590  for (int k2 = 0; k2<=1;++k2)
591  for (int k1 = -1; k1<=0;++k1)
592  if (k3!=0||k2!=0||k1!=0) {
593  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
594  if (d<mind) mind = d;
595  }
596  if (!fail) { // local min found
597  minima[x1][x2][x3]=mind;
598  minima_found = true;
599  }
600  }
601  x2 = n2-1; x1 = n1-1; //vary x3
602  for (x3 = 1; x3<n3-1; ++x3) {
603  mind = mval; fail = false;
604  v = in[x1][x2][x3];
605  for (int k3 = -1; k3<=1;++k3)
606  for (int k2 = -1; k2<=0;++k2)
607  for (int k1 = -1; k1<=0;++k1)
608  if (k3!=0||k2!=0||k1!=0) {
609  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
610  if (d<mind) mind = d;
611  }
612  if (!fail) { // local min found
613  minima[x1][x2][x3]=mind;
614  minima_found = true;
615  }
616  }
617  }
618  // corner border cases 8 in total
619  // corner 000
620  x3=0; x2=0; x1=0;
621  v = in[x1][x2][x3];
622  mind = mval; fail = false;
623  for (int k3 = 0; k3<=1;++k3)
624  for (int k2 = 0; k2<=1;++k2)
625  for (int k1 = 0; k1<=1;++k1)
626  if (k3!=0||k2!=0||k1!=0) {
627  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
628  if (d<mind) mind = d;
629  }
630  if (!fail) {
631  minima[x1][x2][x3]=mind;
632  minima_found = true;
633  }
634  // corner 001
635  x3=n3-1; x2=0; x1=0;
636  v = in[x1][x2][x3];
637  mind = mval; fail = false;
638  for (int k3 = -1; k3<=0;++k3)
639  for (int k2 = 0; k2<=1;++k2)
640  for (int k1 = 0; k1<=1;++k1)
641  if (k3!=0||k2!=0||k1!=0) {
642  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
643  if (d<mind) mind = d;
644  }
645  if (!fail) {
646  minima[x1][x2][x3]=mind;
647  minima_found = true;
648  }
649  // corner 010
650  x3=0; x2=n2-1; x1=0;
651  v = in[x1][x2][x3];
652  mind = mval; fail = false;
653  for (int k3 = 0; k3<=1;++k3)
654  for (int k2 = -1; k2<=0;++k2)
655  for (int k1 = 0; k1<=1;++k1)
656  if (k3!=0||k2!=0||k1!=0) {
657  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
658  if (d<mind) mind = d;
659  }
660  if (!fail) {
661  minima[x1][x2][x3]=mind;
662  minima_found = true;
663  }
664  // corner 011
665  x3=n3-1; x2=n2-1; x1=0;
666  v = in[x1][x2][x3];
667  mind = mval; fail = false;
668  for (int k3 = -1; k3<=0;++k3)
669  for (int k2 = -1; k2<=0;++k2)
670  for (int k1 = 0; k1<=1;++k1)
671  if (k3!=0||k2!=0||k1!=0) {
672  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
673  if (d<mind) mind = d;
674  }
675  if (!fail) {
676  minima[x1][x2][x3]=mind;
677  minima_found = true;
678  }
679  // corner 100
680  x3=0; x2=0; x1=n1-1;
681  v = in[x1][x2][x3];
682  mind = mval; fail = false;
683  for (int k3 = 0; k3<=1;++k3)
684  for (int k2 = 0; k2<=1;++k2)
685  for (int k1 = -1; k1<=0;++k1)
686  if (k3!=0||k2!=0||k1!=0) {
687  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
688  if (d<mind) mind = d;
689  }
690  if (!fail) {
691  minima[x1][x2][x3]=mind;
692  minima_found = true;
693  }
694  // corner 101
695  x3=n3-1; x2=0; x1=n1-1;
696  v = in[x1][x2][x3];
697  mind = mval; fail = false;
698  for (int k3 = -1; k3<=0;++k3)
699  for (int k2 = 0; k2<=1;++k2)
700  for (int k1 = -1; k1<=0;++k1)
701  if (k3!=0||k2!=0||k1!=0) {
702  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
703  if (d<mind) mind = d;
704  }
705  if (!fail) {
706  minima[x1][x2][x3]=mind;
707  minima_found = true;
708  } // corner 110
709  mind = mval; fail = false;
710  x3=0; x2=n2-1; x1=n1-1;
711  v = in[x1][x2][x3];
712  for (int k3 = 0; k3<=1;++k3)
713  for (int k2 = -1; k2<=0;++k2)
714  for (int k1 = -1; k1<=0;++k1)
715  if (k3!=0||k2!=0||k1!=0) {
716  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
717  if (d<mind) mind = d;
718  }
719  if (!fail) {
720  minima[x1][x2][x3]=mind;
721  minima_found = true;
722  } // corner 111
723  mind = mval; fail = false;
724  x3=n3-1; x2=n2-1; x1=n1-1;
725  v = in[x1][x2][x3];
726  for (int k3 = -1; k3<=0;++k3)
727  for (int k2 = -1; k2<=0;++k2)
728  for (int k1 = -1; k1<=0;++k1)
729  if (k3!=0||k2!=0||k1!=0) {
730  d = in[x1+k1][x2+k2][x3+k3]-v; if (d<=thresh) {fail=true; break;}
731  if (d<mind) mind = d;
732  }
733  if (!fail) {
734  minima[x1][x2][x3]=mind;
735  minima_found = true;
736  }
737  return minima_found;
738 }
739 
740 #define VBL_LOCAL_MINIMA_INSTANTIATE(T) \
741 template vbl_array_1d<T > vbl_local_minima(vbl_array_1d<T >const&, T); \
742 template vbl_array_2d<T > vbl_local_minima(vbl_array_2d<T >const&, T); \
743 template vbl_array_3d<T > vbl_local_minima(vbl_array_3d<T >const&, T); \
744 template bool local_minima(vbl_array_1d<T >const&, vbl_array_1d<T >&, T); \
745 template bool local_minima(vbl_array_2d<T >const&, vbl_array_2d<T >&, T); \
746 template bool local_minima(vbl_array_3d<T >const&, vbl_array_3d<T >&, T)
747 
748 #endif // VBL_LOCAL_MINIMA_TXX_
vbl_array_2d
simple 2D array.
Definition: vbl_array_2d.h:25
vbl_array_1d::size
size_type size() const
Definition: vbl_array_1d.h:148
vbl_array_1d
A simple container.
Definition: vbl_array_1d.h:28
vbl_array_3d::get_row3_count
size_type get_row3_count() const
Definition: vbl_array_3d.h:129
vbl_array_2d::rows
size_type rows() const
Return number of rows.
Definition: vbl_array_2d.h:122
vbl_array_2d::cols
size_type cols() const
Return number of columns.
Definition: vbl_array_2d.h:125
vbl_local_minima.h
Find local minima in arrays.
local_minima
bool local_minima(vbl_array_1d< T > const &in, vbl_array_1d< T > &minima, T thresh)
DEPRECATED.
Definition: vbl_local_minima.hxx:13
vbl_array_3d
Templated 3-dimensional array.
Definition: vbl_array_3d.h:38
vbl_array_3d::get_row2_count
size_type get_row2_count() const
Definition: vbl_array_3d.h:128
vbl_array_3d::get_row1_count
size_type get_row1_count() const
Definition: vbl_array_3d.h:127
Generated by   doxygen 1.8.15