Basic Image AlgorithmS Library  2.8.0
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
TrackerBaseInterface.hh
1 /* This file is part of the BIAS library (Basic ImageAlgorithmS).
2 
3  Copyright (C) 2003-2009 (see file CONTACT for details)
4  Multimediale Systeme der Informationsverarbeitung
5  Institut fuer Informatik
6  Christian-Albrechts-Universitaet Kiel
7 
8  BIAS is free software; you can redistribute it and/or modify
9  it under the terms of the GNU Lesser General Public License as published by
10  the Free Software Foundation; either version 2.1 of the License, or
11  (at your option) any later version.
12 
13  BIAS is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU Lesser General Public License for more details.
17 
18  You should have received a copy of the GNU Lesser General Public License
19  along with BIAS; if not, write to the Free Software
20  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA*/
21 #ifndef __TrackerBaseInterface__hh__
22 #define __TrackerBaseInterface__hh__
23 
24 #include <Base/Common/BIASpragmaStart.hh>
25 
26 #include <Base/Debug/Debug.hh>
27 #include <Base/Geometry/HomgPoint2D.hh>
28 #include <Base/Math/Matrix2x2.hh>
29 #include <Base/Image/Image.hh>
30 #include <MathAlgo/Median1D.hh>
31 #include <Filter/FilterMask.hh>
32 
33 #include <vector>
34 
35 #define KLT_TYPE double
36 
37 #define D_TRACKERB_INIT 0x00000001
38 #define D_TRACKERB_KLT 0x00000002
39 #define D_TRACKERB_BILINEAR 0x00000004
40 #define D_TRACKERB_INPUT 0x00000008
41 #define D_TRACKERB_REJECT 0x00000010
42 #define D_TRACKERB_AFFINE 0x00000020
43 
44 #define TRACKERBASE_INVALID_PIXEL -1e10
45 
46 namespace BIAS {
47 
49 
51  TRACKER_SUCCESS = 0, ///< success (error < maxerror)
52  TRACKER_NOSTRUCTURE = -1, ///< no spatial structure is present
53  TRACKER_OUTOFIMAGE = -2, ///< point slid out of image
54  TRACKER_TOOMANYITER = -3, ///< maxiter is reached and error is above maxerror
55  TRACKER_OUTOFROI = -4, ///< a point lies outside of valid region in images
56  TRACKER_MEANABSDIF = -5, ///< point was rejected by MAD criterion
57  TRACKER_IGNOREPOINT = -6, ///< input point was at infinity
58  TRACKER_X84 = -7, ///< point is rejected by X84/X84M
59  TRACKER_ODDAFFINEWARP = -8, ///< affine warp seems degenerated
60  TRACKER_ERRORINCREASED = 1 ///< error increased from last iteration
61  };
62 
63  std::ostream& operator<<(std::ostream& os, const enum ETrackerResult& res);
64 
65  std::ostream& operator<<(std::ostream& os, const enum ERejectionType& rt);
66 
67  /** @class TrackerBaseInterface
68  @ingroup g_tracker
69  @brief Base class for the different tracking algorithms, defining the interfaces
70  for the tracking functions. Also stores some basic options and datas,
71  i.e. the pointer to the image datas. See ExampleTrackerBaseInterface for brief
72  exemplary usage.
73  @author woelk 09/2004, kollmann 07/2006 */
74  template <class StorageType>
75  class BIASMatcher2D_EXPORT TrackerBaseInterface : public Debug
76  {
77  public :
79  virtual ~TrackerBaseInterface();
80 
81  /** @brief Prepare for tracking with prefiltered images.
82 
83  The ROI must be set correctly in gradient images (representing the
84  gradient filter mask size).
85  For good results im1 and im2 should be lowpass filtered images.
86 
87  @param im1 (low-pass filtered) image 1
88  @param im2 (low-pass filtered) image 2
89  @param gradx1 horizontal gradient image of im1
90  @param grady1 vertical gradient image of im1
91  @param gradx2 horizontal gradient image of im2
92  @param grady2 vertical gradient image of im2 */
93  void Init(Image<StorageType>& im1, Image<StorageType>& im2,
94  Image<StorageType>& gradx1, Image<StorageType>& grady1,
95  Image<StorageType>& gradx2, Image<StorageType>& grady2);
96 
97  /** @brief Prepare for tracking with on-demand filtering.
98 
99  Lowpass and gradient filtering is done only where needed, using the
100  given filter masks. All filter masks must have odd size.
101 
102  @param im1 (low-pass filtered) image 1
103  @param im2 (low-pass filtered) image 2
104  @param lowpass_mask filter mask for lowpass filtering
105  @param gradx_mask filter mask for horizontal gradient computation
106  @param grady_mask filter mask for vertical gradient computation */
107  void Init(Image<StorageType>& im1, Image<StorageType>& im2,
108  FilterMask& lowpass_mask, FilterMask& gradx_mask,
109  FilterMask& grady_mask);
110 
111  /** @brief Calculates correspondence from image1 to image2. The two images
112  and corresponding gradients must specified using the function Init()
113  before calling this function.
114 
115  The algorithm searches for best correspondence to point p1 from image1
116  in image2. The search starts at p2 in image 2.
117 
118  All homogenous points must be homogenized (i.e. w=1.0).
119 
120  The algorithms either terminates if the displacement resulting from
121  one iteration drops below *_Maxerror or the maximum number of
122  iterations *_MaxIterations is reached.
123 
124  @param[in] p1 Point in image1
125  @param[in] p2 Algorithm starts searching at this point in image2
126  @param[out] p2tracked Best matching point in image2
127  @param[out] error Change in displacement calculated in last iteration
128  @param[out] iter Number of iterations used
129  @param[out] residuumMAD Mean absolute grey value difference between p1 in
130  image 1 and p2tracked in image 2 over the tracking
131  support window
132  @param[out] residuumMSD Squared grey value difference between p1 in image 1
133  and p2tracked in image 2 summed over the tracking
134  support window
135  @param[out] cov Approximation to the covarinace matrix image1 and
136  result in image2 over the tracking window (= sum of
137  absolute difference divided by number of pixels in
138  tracking value)
139  @param[in] AffinePred Predicted affine matrix (only used in affine tracking)
140  @param[out] AffineResult Result of affine tracking (only used in affine tracking)
141  @return See TRACKER_* */
142  int Track(HomgPoint2D& p1, HomgPoint2D& p2, HomgPoint2D& p2tracked,
143  KLT_TYPE& error, int& iter, KLT_TYPE& residuumMAD,
144  KLT_TYPE& residuumMSD, Matrix<KLT_TYPE>& cov,
145  const Matrix2x2<KLT_TYPE>& AffinePred =
147  Matrix2x2<KLT_TYPE>* AffineResult = NULL);
148 
149  /** @brief Same as Track(), but a correspondence is searched for every
150  point in the vector p1. */
151  void Track(std::vector<HomgPoint2D>& p1,
152  std::vector<HomgPoint2D>& p2,
153  std::vector<HomgPoint2D>& p2tracked,
154  std::vector<KLT_TYPE>& error,
155  std::vector<int>& numiter,
156  std::vector<KLT_TYPE>& residuiMAD,
157  std::vector<KLT_TYPE>& residuiMSD,
158  std::vector<int>& res,
159  std::vector<Matrix<KLT_TYPE> >& cov,
160  const std::vector<Matrix2x2<KLT_TYPE> >& AffinePred =
161  std::vector<Matrix2x2<KLT_TYPE> >(),
162  std::vector<Matrix2x2<KLT_TYPE> >* AffineResult = NULL);
163 
164  inline void SetHalfWinSize(const int hws)
165  { _HalfWinSize=hws; _WinSize=(hws<<1)+1; }
166  inline void SetMaxIterations(const int maxiter)
167  { _MaxIterations=maxiter; }
168  inline void SetMaxError(const KLT_TYPE maxerr)
169  { _MaxError=maxerr; }
170 
171  /** !!!
172  Also used for X84M as maximal residuum. That is features with a
173  residuum <= 'maxres' will be accept automatically!
174  !!!
175  **/
176  inline void SetMaxResiduumMAD(const KLT_TYPE maxres)
177  { _MaxResiduumMAD=maxres; }
178  inline void SetRejectionType(const int rejection_type)
179  { _RejectionType = rejection_type; }
180 
181 
182  /** enable brightness variance and offset invariant computation */
183  virtual void SetAffineBrightnessInvariance(bool bi) {
184  _AffineBrightnessInvariance = bi;
185  }
186 
187  // these are needed in the highlevel Tracker class
188  inline int GetHalfWinSize() { return _HalfWinSize; }
189  inline int GetMaxIterations() { return _MaxIterations; }
190  inline KLT_TYPE GetMaxError() { return _MaxError; }
191  inline int GetRejectionType() { return _RejectionType; }
192 
193  virtual void SetWeightMatrix(const BIAS::Matrix<KLT_TYPE>& weight) {}
194 
195  protected:
196 
197  ///////////////////////////////////////////////////////////////////////////
198  // params
199  ///////////////////////////////////////////////////////////////////////////
200 
201  int _HalfWinSize; ///< use support window of size 2*hws+1
202  int _MaxIterations; ///< iteration stops after _MaxIterations
203  /// iteration stops if position refinement is less than *_MaxError
204  KLT_TYPE _MaxError;
205  KLT_TYPE _MaxResiduumMAD; ///< outlier rejection via MAD and
206  // acception via X84M respectively
207  int _RejectionType; ///< The rejection type: RT_MAD,RT_X84 or RT_X84M
208 
209 
210  /// compute brightness "invariant"
212 
213  /// mean and std deviation of image patch
214  KLT_TYPE mean1_, mean2_, dev1_, dev2_;
215 
216  ///////////////////////////////////////////////////////////////////////////
217 
218  // pointers to images
219  Image<StorageType> *_im1, *_gradim1x, *_gradim1y;
220  Image<StorageType> *_im2, *_gradim2x, *_gradim2y;
221 
222  // filtered images are computed when set to true
224 
225  // lowpass filter computation type
226  enum LowpassFilter
227  {
230  LOWPASS_FILTERMASK_SEPARABLE
231  };
232 
233  // gradient filter computation type
234  enum GradientFilter
235  {
238  GRAD_FILTERMASK_SEPARABLE
239  };
240 
241  LowpassFilter _LowpassFilter;
242  GradientFilter _GradXFilter, _GradYFilter;
243  FilterMask _LowpassMask, _GradXMask, _GradYMask;
244  KLT_TYPE _LowpassMaskSum, _GradXMaskSum, _GradYMaskSum;
246 
247  // the image borders modulo filter and support windows or ROI
248  // for image 1 and image 2
249  int _maxx1, _maxy1, _minx1, _miny1;
250  int _maxx2, _maxy2, _minx2, _miny2;
251 
252  // remember last support window size (for recalculation of gaussian weight)
253  int _LastHalfWinSize, _WinSize;
254 
255  // used in iteration
257  Vector2<KLT_TYPE> _e, _d, _disp;
258  Matrix<KLT_TYPE> _gx, _gy, _bl1, _bl2, _gx1, _gy1, _gx2, _gy2, _gsx, _gsy;
259  KLT_TYPE _gxx, _gxy, _gyy, _idet;
260  Matrix<KLT_TYPE> _weight; // the weight matrix used in TrackWeighted_
261 
262  ///////////////////////////////////////////////////////////////////
263  // tracking
264  ///////////////////////////////////////////////////////////////////
265 
266  /** @brief Interface of all tracking algorithms implemented in
267  derived classes. */
268  virtual int Track_(Vector2<KLT_TYPE>& p1, Vector2<KLT_TYPE>& p2,
269  Vector2<KLT_TYPE>& result, KLT_TYPE& error,
270  int &iter,
271  const Matrix2x2<KLT_TYPE>& AffinePred,
272  Matrix2x2<KLT_TYPE>& AffineResult) = 0;
273 
274  ///////////////////////////////////////////////////////////////////
275  // calculation of cov, sad and ssd
276  ///////////////////////////////////////////////////////////////////
277 
278  /** Uses _bl1, _bl2, _gxx, _gxy and _gyy to calculate the residui
279  and covariance matrix and must hence be called immediatly after the
280  Track_() function call!!
281  @author woelk */
282  virtual void EvaluateResult_(KLT_TYPE& mad,KLT_TYPE& msd,
283  Matrix<KLT_TYPE>& cov);
284 
285  ///////////////////////////////////////////////////////////////////
286  // outlier rejection
287  ///////////////////////////////////////////////////////////////////
288 
289  /** @author woelk */
290  int RejectMAD_(std::vector<KLT_TYPE>& sad, std::vector<int>& res);
291 
292  /** @author woelk */
293  int RejectX84_(std::vector<KLT_TYPE>& ssd, std::vector<KLT_TYPE>& sad,
294  std::vector<int>& res);
295 
296  ///////////////////////////////////////////////////////////////////
297  // helper functions
298  ///////////////////////////////////////////////////////////////////
299 
300  /** resizes the internal matrices _gx1, _gy1, _bl1, _gx2, _gy2 and _bl2
301  @author woelk */
302  inline void Resize_(const int halfwinsize);
303 
304  ////////////////////////////////////////////////////////////////////
305  // bilinear interpolation functions
306  ////////////////////////////////////////////////////////////////////
307 
308  KLT_TYPE ComputeMaskSum(const FilterMask& mask);
309 
310  /** The lowpass filter functions apply a filter to a window of
311  the source image and store the result in a matrix. */
312  //@{
313 
314  /// Lowpass filtering using a kernel filter mask.
315  inline void FilterLowpass_ByMask(int minX, int minY,
316  int rows, int cols,
317  const Image<StorageType>& image,
318  Matrix<KLT_TYPE>& out,
319  const FilterMask& mask,
320  KLT_TYPE maskSum);
321 
322  /// Lowpass filtering using a separable filter mask.
323  inline void FilterLowpass_BySeparableMask(int minX, int minY,
324  int rows, int cols,
325  const Image<StorageType>& image,
326  Matrix<KLT_TYPE>& out,
327  const FilterMask& mask,
328  KLT_TYPE maskSum);
329 
330  /// Optimized lowpass filtering using the binomial 3x3 filter.
331  inline void FilterLowpass_Binomial3x3(int minX, int minY,
332  int rows, int cols,
333  const Image<StorageType>& image,
334  Matrix<KLT_TYPE>& out);
335  //@}
336 
337  /** The gradient filter functions apply a filter to an image window
338  stored in a matrix and store the result in another matrix. */
339  //@{
340 
341  /// Optimized gradient X sobel3x3 filter.
342  inline void Filter_GradXSobel3x3(const Matrix<KLT_TYPE>& in,
343  Matrix<KLT_TYPE>& out);
344 
345  /// Optimized gradient Y sobel3x3 filter.
346  inline void Filter_GradYSobel3x3(const Matrix<KLT_TYPE>& in,
347  Matrix<KLT_TYPE>& out);
348 
349  /// Generic gradient filter.
350  inline void Filter_ByMask(const Matrix<KLT_TYPE>& in,
351  Matrix<KLT_TYPE>& out,
352  const FilterMask& mask,
353  KLT_TYPE maskSum);
354 
355  /// Generic gradient filter using separable filter mask.
356  inline void Filter_BySeparableMask(const Matrix<KLT_TYPE>& in,
357  Matrix<KLT_TYPE>& out,
358  const FilterMask& mask,
359  KLT_TYPE maskSum);
360 
361  //@}
362 
363  /** Fills the matrices _gx1, _gy1 and _bl1 by interpolating gradients
364  and lowpass filtered image at positions between (x-hws, y-hws) and
365  (x+hws, y+hws). */
366  inline void BilinearRegion1_(KLT_TYPE x, KLT_TYPE y, int hws);
367 
368  /** Fills the matrices _gx2, _gy2 and _bl2 by interpolating gradients
369  and lowpass filtered image at positions between (x-hws, y-hws) and
370  (x+hws, y+hws). */
371  inline void BilinearRegion2_(KLT_TYPE x, KLT_TYPE y, int hws);
372 
373  /** Bilinear filtering using images as source (for prefiltered images) */
374  inline KLT_TYPE
375  BilinearRegionFromImages_(KLT_TYPE x, KLT_TYPE y, int hws,
376  Image<StorageType>& image,
377  Image<StorageType>& gradx,
378  Image<StorageType>& grady,
379  Matrix<KLT_TYPE>& bl,
380  Matrix<KLT_TYPE>& gx,
381  Matrix<KLT_TYPE>& gy);
382 
383  /** Bilinear filtering using matrices as source */
384  inline KLT_TYPE
385  BilinearRegionFromMatrices_(KLT_TYPE x, KLT_TYPE y, int hws,
386  const Matrix<KLT_TYPE>& image,
387  const Matrix<KLT_TYPE>& gradx,
388  const Matrix<KLT_TYPE>& grady,
389  Matrix<KLT_TYPE>& bl,
390  Matrix<KLT_TYPE>& gx,
391  Matrix<KLT_TYPE>& gy);
392 
393  /** bring region to mean zero and stddev 1 for brightness invariance
394  @return scale of region which has been applied */
395  inline KLT_TYPE NormalizeRegion_(Matrix<KLT_TYPE>& bl,
396  Matrix<KLT_TYPE>& gx,
397  Matrix<KLT_TYPE>& gy,
398  const KLT_TYPE& min_value = -1e9);
399 
400  }; // class
401 
402  ////////////////////////////////////////////////////////////////////
403  // inline functions
404  ////////////////////////////////////////////////////////////////////
405 
406  template <class StorageType>
408  Resize_(const int halfwinsize)
409  {
410  _LastHalfWinSize=halfwinsize;
411  _WinSize=(halfwinsize<<1)+1;
412  _gx1.newsize(_WinSize, _WinSize);
413  _gy1.newsize(_WinSize, _WinSize);
414  _gx2.newsize(_WinSize, _WinSize);
415  _gy2.newsize(_WinSize, _WinSize);
416  _gsx.newsize(_WinSize, _WinSize);
417  _gsy.newsize(_WinSize, _WinSize);
418  _bl1.newsize(_WinSize, _WinSize);
419  _bl2.newsize(_WinSize, _WinSize);
420  }
421 
422  template <class StorageType>
424  FilterLowpass_ByMask(int minX, int minY, int rows, int cols,
425  const Image<StorageType>& image,
426  Matrix<KLT_TYPE>& out,
427  const FilterMask& mask,
428  KLT_TYPE maskSum)
429  {
430  BIASASSERT(!mask.IsSeparable());
431 
432  const Matrix<FM_FLOAT>& kernel = *mask.GetKernelf();
433 
434  int borderRows = (kernel.GetRows()-1) / 2;
435  int borderCols = (kernel.GetCols()-1) / 2;
436 
437  out.newsize(rows, cols);
438 
439  const StorageType** ida = image.GetImageDataArray();
440 
441  for (int i = 0; i < rows; i++)
442  {
443  for (int j = 0; j < cols; j++)
444  {
445  out[i][j] = 0;
446  for (int k = 0; k < (int)kernel.GetRows(); k++)
447  {
448  for (int l = 0; l < (int)kernel.GetCols(); l++)
449  {
450  out[i][j] +=
451  kernel[k][l] * ida[i+k-borderRows+minY][j+l-borderCols+minX];
452  }
453  }
454  out[i][j] /= maskSum;
455  }
456  }
457  }
458 
459  template <class StorageType>
461  FilterLowpass_BySeparableMask(int minX, int minY, int rows, int cols,
462  const Image<StorageType>& image,
463  Matrix<KLT_TYPE>& out,
464  const FilterMask& mask,
465  KLT_TYPE maskSum)
466  {
467  BIASASSERT(mask.IsSeparable());
468 
469  const Vector<FM_FLOAT>& vert = *mask.GetSepfv();
470  const Vector<FM_FLOAT>& horz = *mask.GetSepfh();
471 
472  int borderRows = (vert.Size()-1) / 2;
473  int borderCols = (horz.Size()-1) / 2;
474 
475  out.newsize(rows, cols);
476 
477  const StorageType** ida = image.GetImageDataArray();
478 
479  for (int i = 0; i < rows; i++)
480  {
481  for (int j = 0; j < cols; j++)
482  {
483  out[i][j] = 0;
484  for (int k = 0; k < (int)vert.Size(); k++)
485  {
486  for (int l = 0; l < (int)horz.Size(); l++)
487  {
488  out[i][j] +=
489  vert[k] * horz[l]
490  * ida[i+k-borderRows+minY][j+l-borderCols+minX];
491  }
492  }
493  out[i][j] /= maskSum;
494  }
495  }
496  }
497 
498  template <class StorageType>
500  FilterLowpass_Binomial3x3(int minX, int minY, int rows, int cols,
501  const Image<StorageType>& image,
502  Matrix<KLT_TYPE>& out)
503  {
504  out.newsize(rows, cols);
505 
506  const StorageType** ida = image.GetImageDataArray();
507 
508  for (int i = 0; i < rows; i++)
509  {
510  for (int j = 0; j < cols; j++)
511  {
512  int x = j + minX;
513  int y = i + minY;
514 
515  /* 1 2 1
516  2 4 2
517  1 2 1 */
518  out[i][j] =
519  ida[y-1][x-1] + 2*ida[y-1][x] + ida[y-1][x+1]
520  +2*ida[y ][x-1] + 4*ida[y ][x] + 2*ida[y ][x+1]
521  + ida[y+1][x-1] + 2*ida[y+1][x] + ida[y+1][x+1];
522  out[i][j] /= 16.0;
523  }
524  }
525  }
526 
527  template <class StorageType>
530  {
531  int rows = in.GetRows();
532  int cols = in.GetCols();
533 
534  out.newsize(rows, cols);
535 
536  for (int i = 1; i < rows-1; i++)
537  {
538  for (int j = 1; j < cols-1; j++)
539  {
540  /* -1 0 1
541  -2 0 2
542  -1 0 1 */
543  out[i][j] =
544  - in[i-1][j-1] + in[i-1][j+1] // row 0
545  -2*in[i ][j-1] + 2*in[i ][j+1] // row 1
546  - in[i+1][j-1] + in[i+1][j+1]; // row 2
547  out[i][j] /= 8.0;
548  }
549  }
550  }
551 
552  template <class StorageType>
555  {
556  int rows = in.GetRows();
557  int cols = in.GetCols();
558 
559  out.newsize(rows, cols);
560 
561  for (int i = 1; i < rows-1; i++)
562  {
563  for (int j = 1; j < cols-1; j++)
564  {
565  /* -1 -2 -1
566  0 0 0
567  1 2 1 */
568  out[i][j] =
569  -in[i-1][j-1] - 2*in[i-1][j] - in[i-1][j+1] // row 0
570  +in[i+1][j-1] + 2*in[i+1][j] + in[i+1][j+1]; // row 2
571  out[i][j] /= (StorageType)8.0;
572  }
573  }
574  }
575 
576  template <class StorageType>
579  Matrix<KLT_TYPE>& out,
580  const FilterMask& mask,
581  KLT_TYPE maskSum)
582  {
583  BIASASSERT(!mask.IsSeparable());
584 
585  const Matrix<FM_FLOAT>& kernel = *mask.GetKernelf();
586 
587  int borderRows = (kernel.GetRows()-1) / 2;
588  int borderCols = (kernel.GetCols()-1) / 2;
589 
590  int rows = in.GetRows();
591  int cols = in.GetCols();
592 
593  out.newsize(rows, cols);
594 
595  for (int i = borderRows; i < rows-borderRows; i++)
596  {
597  for (int j = borderCols; j < cols-borderCols; j++)
598  {
599  out[i][j] = 0;
600  for (int k = 0; k < (int)kernel.GetRows(); k++)
601  {
602  for (int l = 0; l < (int)kernel.GetCols(); l++)
603  {
604  out[i][j] += kernel[k][l] * in[i+k-borderRows][j+l-borderCols];
605  }
606  }
607  out[i][j] /= maskSum;
608  }
609  }
610  }
611 
612  template <class StorageType>
615  Matrix<KLT_TYPE>& out,
616  const FilterMask& mask,
617  KLT_TYPE maskSum)
618  {
619  BIASASSERT(mask.IsSeparable());
620 
621  const Vector<FM_FLOAT>& vert = *mask.GetSepfv();
622  const Vector<FM_FLOAT>& horz = *mask.GetSepfh();
623 
624  int borderRows = (vert.Size()-1) / 2;
625  int borderCols = (horz.Size()-1) / 2;
626 
627  int rows = in.GetRows();
628  int cols = in.GetCols();
629 
630  out.newsize(rows, cols);
631 
632  for (int i = borderRows; i < rows-borderRows; i++)
633  {
634  for (int j = borderCols; j < cols-borderCols; j++)
635  {
636  out[i][j] = 0;
637  for (int k = 0; k < (int)vert.Size(); k++)
638  {
639  for (int l = 0; l < (int)horz.Size(); l++)
640  {
641  out[i][j] +=
642  vert[k] * horz[l]
643  * in[i+k-borderRows][j+l-borderCols];
644  }
645  }
646  out[i][j] /= maskSum;
647  }
648  }
649  }
650 
651  template <class StorageType>
653  BilinearRegion1_(KLT_TYPE x, KLT_TYPE y, int hws)
654  {
655  if (_ComputeFilteredImages)
656  {
657  Matrix<KLT_TYPE> lowpass;
658  int lowpassHWS = hws + _LowpassWindowExtraSize;
659  int lowpassWS = 2*lowpassHWS+2;
660  int lowpassMinX = (int)floor(x) - lowpassHWS;
661  int lowpassMinY = (int)floor(y) - lowpassHWS;
662 
663  switch (_LowpassFilter)
664  {
665  case LOWPASS_BINOMIAL3X3:
666  FilterLowpass_Binomial3x3(lowpassMinX, lowpassMinY,
667  lowpassWS, lowpassWS,
668  *_im1, lowpass);
669  break;
670  case LOWPASS_FILTERMASK:
671  FilterLowpass_ByMask(lowpassMinX, lowpassMinY,
672  lowpassWS, lowpassWS,
673  *_im1, lowpass,
674  _LowpassMask, _LowpassMaskSum);
675  break;
676  case LOWPASS_FILTERMASK_SEPARABLE:
677  FilterLowpass_BySeparableMask(lowpassMinX, lowpassMinY,
678  lowpassWS, lowpassWS,
679  *_im1, lowpass,
680  _LowpassMask, _LowpassMaskSum);
681  break;
682  }
683 
684  Matrix<KLT_TYPE> gradx;
685  switch (_GradXFilter)
686  {
687  case GRAD_SOBEL3X3:
688  Filter_GradXSobel3x3(lowpass, gradx);
689  break;
690  case GRAD_FILTERMASK:
691  Filter_ByMask(lowpass, gradx, _GradXMask, _GradXMaskSum);
692  break;
693  case GRAD_FILTERMASK_SEPARABLE:
694  Filter_BySeparableMask(lowpass, gradx, _GradXMask, _GradXMaskSum);
695  break;
696  }
697 
698  Matrix<KLT_TYPE> grady;
699  switch (_GradYFilter)
700  {
701  case GRAD_SOBEL3X3:
702  Filter_GradYSobel3x3(lowpass, grady);
703  break;
704  case GRAD_FILTERMASK:
705  Filter_ByMask(lowpass, grady, _GradYMask, _GradYMaskSum);
706  break;
707  case GRAD_FILTERMASK_SEPARABLE:
708  Filter_BySeparableMask(lowpass, grady, _GradYMask, _GradYMaskSum);
709  break;
710  }
711 
712  dev1_ = this->BilinearRegionFromMatrices_(x - lowpassMinX,
713  y - lowpassMinY,
714  hws,
715  lowpass, gradx, grady,
716  _bl1, _gx1, _gy1);
717  }
718  else
719  {
720  dev1_ = this->BilinearRegionFromImages_(x, y, hws,
721  *_im1, *_gradim1x, *_gradim1y,
722  _bl1, _gx1, _gy1);
723  }
724  }
725 
726  template <class StorageType>
728  BilinearRegion2_(KLT_TYPE x, KLT_TYPE y, int hws)
729  {
730 
731  if (_ComputeFilteredImages)
732  {
733  Matrix<KLT_TYPE> lowpass;
734  int lowpassHWS = hws + _LowpassWindowExtraSize;
735  int lowpassWS = 2*lowpassHWS+2;
736  int lowpassMinX = (int)floor(x) - lowpassHWS;
737  int lowpassMinY = (int)floor(y) - lowpassHWS;
738 
739  switch (_LowpassFilter)
740  {
741  case LOWPASS_BINOMIAL3X3:
742  FilterLowpass_Binomial3x3(lowpassMinX, lowpassMinY,
743  lowpassWS, lowpassWS,
744  *_im2, lowpass);
745  break;
746  case LOWPASS_FILTERMASK:
747  FilterLowpass_ByMask(lowpassMinX, lowpassMinY,
748  lowpassWS, lowpassWS,
749  *_im2, lowpass,
750  _LowpassMask, _LowpassMaskSum);
751  break;
752  case LOWPASS_FILTERMASK_SEPARABLE:
753  FilterLowpass_BySeparableMask(lowpassMinX, lowpassMinY,
754  lowpassWS, lowpassWS,
755  *_im2, lowpass,
756  _LowpassMask, _LowpassMaskSum);
757  break;
758  }
759 
760  Matrix<KLT_TYPE> gradx;
761  switch (_GradXFilter)
762  {
763  case GRAD_SOBEL3X3:
764  Filter_GradXSobel3x3(lowpass, gradx);
765  break;
766  case GRAD_FILTERMASK:
767  Filter_ByMask(lowpass, gradx, _GradXMask, _GradXMaskSum);
768  break;
769  case GRAD_FILTERMASK_SEPARABLE:
770  Filter_BySeparableMask(lowpass, gradx, _GradXMask, _GradXMaskSum);
771  break;
772  }
773 
774  Matrix<KLT_TYPE> grady;
775  switch (_GradYFilter)
776  {
777  case GRAD_SOBEL3X3:
778  Filter_GradYSobel3x3(lowpass, grady);
779  break;
780  case GRAD_FILTERMASK:
781  Filter_ByMask(lowpass, grady, _GradYMask, _GradYMaskSum);
782  break;
783  case GRAD_FILTERMASK_SEPARABLE:
784  Filter_BySeparableMask(lowpass, grady, _GradYMask, _GradYMaskSum);
785  break;
786  }
787 
788  dev2_ = this->BilinearRegionFromMatrices_(x - lowpassMinX,
789  y - lowpassMinY,
790  hws,
791  lowpass, gradx, grady,
792  _bl2, _gx2, _gy2);
793  } else {
794  dev2_ = this->BilinearRegionFromImages_(x, y, hws,
795  *_im2, *_gradim2x, *_gradim2y,
796  _bl2, _gx2, _gy2);
797  }
798  }
799 
800  template <class StorageType>
802  BilinearRegionFromImages_(KLT_TYPE x, KLT_TYPE y, int hws,
803  Image<StorageType>& image,
804  Image<StorageType>& gradx,
805  Image<StorageType>& grady,
806  Matrix<KLT_TYPE>& bl,
807  Matrix<KLT_TYPE>& gx,
808  Matrix<KLT_TYPE>& gy) {
809  StorageType** image_ida = image.GetImageDataArray();
810  StorageType** gradx_ida = gradx.GetImageDataArray();
811  StorageType** grady_ida = grady.GetImageDataArray();
812 
813  int x_floor = (int)floor(x);
814  int y_floor = (int)floor(y);
815  int x_floor_end = x_floor+hws+1;
816  int y_floor_end = y_floor+hws+1;
817  KLT_TYPE brdy = y - (KLT_TYPE)(y_floor);
818  KLT_TYPE brdx = x - (KLT_TYPE)(x_floor);
819  KLT_TYPE brdxy = brdx*brdy;
820 
821  int rf, rc, cf, cc, r, c;
822  for (rf=y_floor-hws, rc=rf+1, r=0; rf<y_floor_end; rf++, rc++, r++) {
823  for (cf=x_floor-hws, cc=cf+1, c=0; cf<x_floor_end; cf++, cc++, c++) {
824  KLT_TYPE ul, ur, ll, lr;
825 
826  ul = (KLT_TYPE)image_ida[rf][cf];
827  ur = (KLT_TYPE)image_ida[rf][cc];
828  ll = (KLT_TYPE)image_ida[rc][cf];
829  lr = (KLT_TYPE)image_ida[rc][cc];
830  bl[r][c] = (ul + brdy*(ll - ul) + brdx*(ur - ul) +
831  brdxy*(ul + lr - ll - ur));
832  ul = (KLT_TYPE)(gradx_ida[rf][cf]);
833  ur = (KLT_TYPE)(gradx_ida[rf][cc]);
834  ll = (KLT_TYPE)(gradx_ida[rc][cf]);
835  lr = (KLT_TYPE)(gradx_ida[rc][cc]);
836  gx[r][c] = (ul + brdy*(ll - ul) + brdx*(ur - ul) +
837  brdxy*(ul + lr - ll - ur));
838  ul = (KLT_TYPE)(grady_ida[rf][cf]);
839  ur = (KLT_TYPE)(grady_ida[rf][cc]);
840  ll = (KLT_TYPE)(grady_ida[rc][cf]);
841  lr = (KLT_TYPE)(grady_ida[rc][cc]);
842  gy[r][c] = (ul + brdy*(ll - ul) + brdx*(ur - ul) +
843  brdxy*(ul + lr - ll - ur));
844  }
845  }
846  if (_AffineBrightnessInvariance) return NormalizeRegion_(bl, gx, gy);
847  return 1.0;
848  }
849 
850  template <class StorageType>
852  BilinearRegionFromMatrices_(KLT_TYPE x, KLT_TYPE y, int hws,
853  const Matrix<KLT_TYPE>& image,
854  const Matrix<KLT_TYPE>& gradx,
855  const Matrix<KLT_TYPE>& grady,
856  Matrix<KLT_TYPE>& bl,
857  Matrix<KLT_TYPE>& gx,
858  Matrix<KLT_TYPE>& gy) {
859  int x_floor = (int)floor(x);
860  int y_floor = (int)floor(y);
861  int x_floor_end = x_floor+hws+1;
862  int y_floor_end = y_floor+hws+1;
863  KLT_TYPE brdy = y - (KLT_TYPE)(y_floor);
864  KLT_TYPE brdx = x - (KLT_TYPE)(x_floor);
865  KLT_TYPE brdxy = brdx*brdy;
866 
867  int rf, rc, cf, cc, r, c;
868  for (rf=y_floor-hws, rc=rf+1, r=0; rf<y_floor_end; rf++, rc++, r++) {
869  for (cf=x_floor-hws, cc=cf+1, c=0; cf<x_floor_end; cf++, cc++, c++) {
870  KLT_TYPE ul, ur, ll, lr;
871 
872  ul = (KLT_TYPE)image[rf][cf];
873  ur = (KLT_TYPE)image[rf][cc];
874  ll = (KLT_TYPE)image[rc][cf];
875  lr = (KLT_TYPE)image[rc][cc];
876  bl[r][c] = (ul + brdy*(ll - ul) + brdx*(ur - ul) +
877  brdxy*(ul + lr - ll - ur));
878  ul = (KLT_TYPE)(gradx[rf][cf]);
879  ur = (KLT_TYPE)(gradx[rf][cc]);
880  ll = (KLT_TYPE)(gradx[rc][cf]);
881  lr = (KLT_TYPE)(gradx[rc][cc]);
882  gx[r][c] = (ul + brdy*(ll - ul) + brdx*(ur - ul) +
883  brdxy*(ul + lr - ll - ur));
884  ul = (KLT_TYPE)(grady[rf][cf]);
885  ur = (KLT_TYPE)(grady[rf][cc]);
886  ll = (KLT_TYPE)(grady[rc][cf]);
887  lr = (KLT_TYPE)(grady[rc][cc]);
888  gy[r][c] = (ul + brdy*(ll - ul) + brdx*(ur - ul) +
889  brdxy*(ul + lr - ll - ur));
890  }
891  }
892  if (_AffineBrightnessInvariance) return NormalizeRegion_(bl, gx, gy);
893  return 1.0;
894  }
895 
896  template <class StorageType>
899  Matrix<KLT_TYPE>& gx,
900  Matrix<KLT_TYPE>& gy,
901  const KLT_TYPE& min_value) {
902 
903  // mean value of patch
904  KLT_TYPE mean = 0;
905 
906  int matrixsize = bl.num_cols()*bl.num_rows();
907  // number of pixels above a certain threshold (for mask tracker)
908  int validpixels = 0;
909 
910  // compute mean on not-black pixels
911  KLT_TYPE *pbl = bl.GetData()-1, *pblend = bl.GetData() + matrixsize;
912  while (++pbl<pblend) {
913  if (*pbl > min_value) {
914  mean += *pbl;
915  validpixels++;
916  } else *pbl = min_value;
917  }
918  // normalize mean for stddev computation
919  if (validpixels>0) mean *= 1.0 / KLT_TYPE(validpixels);
920 
921 
922  // estimate stddev
923  KLT_TYPE stddev = 0.0;
924  pbl = bl.GetData()-1;
925  while (++pbl<pblend) {
926  if (*pbl > min_value) {
927  stddev += (*pbl - mean)*(*pbl - mean);
928  }
929  }
930  if (fabs(stddev)>1e-10) {
931  stddev = sqrt(KLT_TYPE(validpixels)/stddev);
932  } else {
933  stddev = 1.0;
934  }
935 
936  // normalize image and gradient values such that new image patch has
937  // zero mean and stddev 1
938  KLT_TYPE *pgx = gx.GetData(),*pgy = gy.GetData();
939  pbl = bl.GetData()-1;
940  while (++pbl<pblend) {
941  if (*pbl > min_value) {
942  *pbl = (*pbl-mean) * stddev;
943  *pgx *= stddev;
944  *pgy *= stddev;
945  } else {
946  *pgx = *pgy = 0;
947  *pbl = TRACKERBASE_INVALID_PIXEL;
948  }
949  pgx++; pgy++;
950  }
951  return stddev;
952  }
953 
954 
955 } // namespace
956 
957 #include <Base/Common/BIASpragmaEnd.hh>
958 
959 #endif // __TrackerBaseInterface__hh__
class HomgPoint2D describes a point with 2 degrees of freedom in projective coordinates.
Definition: HomgPoint2D.hh:67
Subscript num_cols() const
Definition: cmat.h:320
bool IsSeparable() const
checks if the kernel is separable
Definition: FilterMask.hh:136
KLT_TYPE BilinearRegionFromMatrices_(KLT_TYPE x, KLT_TYPE y, int hws, const Matrix< KLT_TYPE > &image, const Matrix< KLT_TYPE > &gradx, const Matrix< KLT_TYPE > &grady, Matrix< KLT_TYPE > &bl, Matrix< KLT_TYPE > &gx, Matrix< KLT_TYPE > &gy)
Bilinear filtering using matrices as source.
point slid out of image
Base class for the different tracking algorithms, defining the interfaces for the tracking functions...
int _RejectionType
The rejection type: RT_MAD,RT_X84 or RT_X84M.
const Matrix< FM_FLOAT > * GetKernelf() const
returns pointer to float filter mask
Definition: FilterMask.hh:116
void BilinearRegion2_(KLT_TYPE x, KLT_TYPE y, int hws)
Fills the matrices _gx2, _gy2 and _bl2 by interpolating gradients and lowpass filtered image at posit...
Matrix< T > & newsize(Subscript M, Subscript N)
Definition: cmat.h:269
void SetMaxResiduumMAD(const KLT_TYPE maxres)
!!! Also used for X84M as maximal residuum.
bool _AffineBrightnessInvariance
compute brightness &quot;invariant&quot;
void Filter_ByMask(const Matrix< KLT_TYPE > &in, Matrix< KLT_TYPE > &out, const FilterMask &mask, KLT_TYPE maskSum)
Generic gradient filter.
void Filter_BySeparableMask(const Matrix< KLT_TYPE > &in, Matrix< KLT_TYPE > &out, const FilterMask &mask, KLT_TYPE maskSum)
Generic gradient filter using separable filter mask.
unsigned int Size() const
length of the vector
Definition: Vector.hh:134
void SetMaxError(const KLT_TYPE maxerr)
void SetHalfWinSize(const int hws)
void FilterLowpass_ByMask(int minX, int minY, int rows, int cols, const Image< StorageType > &image, Matrix< KLT_TYPE > &out, const FilterMask &mask, KLT_TYPE maskSum)
The lowpass filter functions apply a filter to a window of the source image and store the result in a...
void SetMaxIterations(const int maxiter)
const Vector< FM_FLOAT > * GetSepfh() const
returns pointer to horiz float vector for sep.
Definition: FilterMask.hh:128
success (error &lt; maxerror)
unsigned int GetRows() const
Definition: Matrix.hh:202
no spatial structure is present
KLT_TYPE _MaxError
iteration stops if position refinement is less than *_MaxError
affine warp seems degenerated
input point was at infinity
void FilterLowpass_Binomial3x3(int minX, int minY, int rows, int cols, const Image< StorageType > &image, Matrix< KLT_TYPE > &out)
Optimized lowpass filtering using the binomial 3x3 filter.
const Vector< FM_FLOAT > * GetSepfv() const
returns pointer to vertical float vector for sep.
Definition: FilterMask.hh:132
KLT_TYPE BilinearRegionFromImages_(KLT_TYPE x, KLT_TYPE y, int hws, Image< StorageType > &image, Image< StorageType > &gradx, Image< StorageType > &grady, Matrix< KLT_TYPE > &bl, Matrix< KLT_TYPE > &gx, Matrix< KLT_TYPE > &gy)
Bilinear filtering using images as source (for prefiltered images)
The image template class for specific storage types.
Definition: Image.hh:78
T * GetData()
get the pointer to the data array of the matrix (for faster direct memeory access) ...
Definition: Matrix.hh:185
point is rejected by X84/X84M
std::ostream & operator<<(std::ostream &os, const Array2D< T > &arg)
Definition: Array2D.hh:260
void BilinearRegion1_(KLT_TYPE x, KLT_TYPE y, int hws)
Fills the matrices _gx1, _gy1 and _bl1 by interpolating gradients and lowpass filtered image at posit...
int _HalfWinSize
use support window of size 2*hws+1
unsigned int GetCols() const
Definition: Matrix.hh:204
virtual void SetAffineBrightnessInvariance(bool bi)
enable brightness variance and offset invariant computation
void SetRejectionType(const int rejection_type)
int _MaxIterations
iteration stops after _MaxIterations
void Filter_GradYSobel3x3(const Matrix< KLT_TYPE > &in, Matrix< KLT_TYPE > &out)
Optimized gradient Y sobel3x3 filter.
KLT_TYPE _MaxResiduumMAD
outlier rejection via MAD and
Subscript num_rows() const
Definition: cmat.h:319
a point lies outside of valid region in images
point was rejected by MAD criterion
void FilterLowpass_BySeparableMask(int minX, int minY, int rows, int cols, const Image< StorageType > &image, Matrix< KLT_TYPE > &out, const FilterMask &mask, KLT_TYPE maskSum)
Lowpass filtering using a separable filter mask.
A filter mask (or a kernel) used for convolution.
Definition: FilterMask.hh:61
void Filter_GradXSobel3x3(const Matrix< KLT_TYPE > &in, Matrix< KLT_TYPE > &out)
The gradient filter functions apply a filter to an image window stored in a matrix and store the resu...
error increased from last iteration
void Resize_(const int halfwinsize)
resizes the internal matrices _gx1, _gy1, _bl1, _gx2, _gy2 and _bl2
KLT_TYPE NormalizeRegion_(Matrix< KLT_TYPE > &bl, Matrix< KLT_TYPE > &gx, Matrix< KLT_TYPE > &gy, const KLT_TYPE &min_value=-1e9)
bring region to mean zero and stddev 1 for brightness invariance
maxiter is reached and error is above maxerror
virtual void SetWeightMatrix(const BIAS::Matrix< KLT_TYPE > &weight)
const StorageType ** GetImageDataArray() const
overloaded GetImageDataArray() from ImageBase
Definition: Image.hh:152