Basic Image AlgorithmS Library  2.8.0
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Gauss.cpp
1 /*
2 This file is part of the BIAS library (Basic ImageAlgorithmS).
3 
4 Copyright (C) 2003-2009 (see file CONTACT for details)
5  Multimediale Systeme der Informationsverarbeitung
6  Institut fuer Informatik
7  Christian-Albrechts-Universitaet Kiel
8 
9 
10 BIAS is free software; you can redistribute it and/or modify
11 it under the terms of the GNU Lesser General Public License as published by
12 the Free Software Foundation; either version 2.1 of the License, or
13 (at your option) any later version.
14 
15 BIAS is distributed in the hope that it will be useful,
16 but WITHOUT ANY WARRANTY; without even the implied warranty of
17 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 GNU Lesser General Public License for more details.
19 
20 You should have received a copy of the GNU Lesser General Public License
21 along with BIAS; if not, write to the Free Software
22 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 */
24 
25 #include "Gauss.hh"
26 #include <Base/Image/ImageConvert.hh>
27 
28 using namespace BIAS;
29 using namespace std;
30 
31 //////////////////////////////////////////////////////////////////////////
32 // implementation
33 //////////////////////////////////////////////////////////////////////////
34 
35 
36 template <class InputStorageType, class OutputStorageType>
39  : Convolution<InputStorageType, OutputStorageType>()
40 {
41  _GaussSigma = 0.7;
42  _GaussRatio = 0.01;
43  _LastRatio = _LastSigma = -1.0;
44  _GaussIO = NULL;
45  _GaussOO = NULL;
46 }
47 
48 template <class InputStorageType, class OutputStorageType>
51  : _GaussSigma(other._GaussSigma),
52  _GaussRatio(other._GaussRatio),
53  _LastSigma(other._LastSigma),
54  _LastRatio(other._LastRatio)
55 {
56  if (other._GaussIO!=NULL) _GaussIO = other._GaussIO->Clone();
57  else _GaussIO = NULL;
58  if (other._GaussOO!=NULL) _GaussOO = other._GaussOO->Clone();
59  else _GaussOO = NULL;
60 }
61 
62 template <class StorageType, class KernelType>
64 {
65  if (_GaussIO!=NULL) delete _GaussIO;
66  if (_GaussOO!=NULL) delete _GaussOO;
67 }
68 
69 
70 void DecreaseROI(const ROI& RoiIn, ROI& RoiOut,
71  int hws) {
72  unsigned minx, miny, maxx, maxy;
73  RoiIn.GetCorners(minx, miny, maxx, maxy);
74 
75  minx += hws;
76  miny += hws;
77  maxx -= hws;
78  maxy -= hws;
79  if (minx>maxx) minx = maxx = (minx+maxx)/2;
80  if (miny>maxy) miny = maxy = (miny+maxy)/2;
81  RoiOut.SetCorners(minx, miny, maxx, maxy);
82 }
83 
84 
85 ////////////////////////////////////////////////////////////////////////////
86 
87 template <class InputStorageType, class OutputStorageType>
90 {
91  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter\n");
92  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
93  _CalculateKernels(_GaussSigma, _GaussRatio);
94 
95  if (_fm.IsSize1x1()){ // do not filter if filter mask is 1x1
96  return ImageConvert::ConvertST(src, dst, dst.GetStorageType());
97  }
98  if ((src.GetChannelCount()==1) && (this->GetBorderHandling() == TBH_valid))
99  return FilterGrey(src, dst);
101 }
102 
103 template <class InputStorageType, class OutputStorageType>
106 {
107  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter\n");
108  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
109  _CalculateKernels(_GaussSigma, _GaussRatio);
111 }
112 
113 template <class InputStorageType, class OutputStorageType>
116 {
117  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter\n");
118  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
119  _CalculateKernels(_GaussSigma, _GaussRatio);
121 }
122 
123 
124 template <class InputStorageType, class OutputStorageType>
128  bool ResetROI)
129 {
130  if (_GaussIO==NULL)
132  if (_GaussOO==NULL)
134 
135  // this is the final target sigma:
136  double smoothsigma = _GaussSigma;
137 
138  // compute maximum sigma for the different kernel sizes
139  // approximate rule:
140  // if you want to neglect filter mask entries below 10% of the maximum,
141  // your gauss'es half window size should be of size 5*sigma
142  const double sqrtlogratio = 1.0 / sqrt(-2.0 * log(_GaussRatio));
143  const double maxsigma3x3 = sqrtlogratio;
144  const double maxsigma5x5 = 2.0 * sqrtlogratio;
145  const double maxsigma7x7 = 3.0 * sqrtlogratio;
146  const double maxsigma9x9 = 4.0 * sqrtlogratio;
147  const double maxsigma11x11 = 5.0 * sqrtlogratio;
148  const double maxsigma13x13 = 6.0 * sqrtlogratio;
149  const double maxsigmaperpass = maxsigma13x13;
150 
151  // compute how many passes we would need with the largest filter
152  int MultiPass = (int)rint(smoothsigma * smoothsigma /
153  (maxsigmaperpass*maxsigmaperpass));
154  if (MultiPass>1) {
155  // we need multiple passes
156  smoothsigma /= sqrtf(float(MultiPass));
157  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass="<<MultiPass);
158  } else {
159  BIASDOUT(D_GAUSS_FAST_GREY,"Single pass is sufficient.");
160  }
161  _GaussIO->SetSigma(smoothsigma);
162  _GaussIO->SetRatio(_GaussRatio);
163  _GaussOO->SetSigma(smoothsigma);
164  _GaussOO->SetRatio(_GaussRatio);
165  if (_GaussSigma<maxsigma3x3) {
166  BIASDOUT(D_GAUSS_FAST_GREY,"Run 3x3");
167  _GaussIO->Filter3x3Grey(src, dst);
168  } else if (_GaussSigma<maxsigma5x5) {
169  BIASDOUT(D_GAUSS_FAST_GREY,"Run 5x5");
170  _GaussIO->Filter5x5Grey(src, dst);
171  } else if (_GaussSigma<maxsigma7x7) {
172  BIASDOUT(D_GAUSS_FAST_GREY,"Run 7x7");
173  _GaussIO->Filter7x7Grey(src,dst);
174  } else if (_GaussSigma<maxsigma9x9){
175  BIASDOUT(D_GAUSS_FAST_GREY,"Run 9x9");
176  _GaussIO->Filter9x9Grey(src, dst);
177  } else if (_GaussSigma<maxsigma11x11){
178  BIASDOUT(D_GAUSS_FAST_GREY,"Run 11x11");
179  _GaussIO->Filter11x11Grey(src, dst);
180  } else {
181  BIASDOUT(D_GAUSS_FAST_GREY,"Run 13x13");
182  _GaussIO->Filter13x13Grey(src, dst);
183  }
184  while (MultiPass-->1) {
185  if (ResetROI) dst.UnsetROI();
186  if (_GaussSigma<maxsigma3x3) {
187  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass: Run 3x3");
188  _GaussOO->Filter3x3Grey(dst, dst);
189  } else if (_GaussSigma<maxsigma5x5) {
190  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass: Run 5x5");
191  _GaussOO->Filter5x5Grey(dst, dst);
192  } else if (_GaussSigma<maxsigma7x7) {
193  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass: Run 7x7");
194  _GaussOO->Filter7x7Grey(dst, dst);
195  } else if (_GaussSigma<maxsigma9x9) {
196  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass: Run 9x9");
197  _GaussOO->Filter9x9Grey(dst, dst);
198  } else if (_GaussSigma<maxsigma11x11){
199  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass: Run 11x11");
200  _GaussOO->Filter11x11Grey(dst, dst);
201  } else {
202  BIASDOUT(D_GAUSS_FAST_GREY,"MultiPass: Run 13x13");
203  _GaussOO->Filter13x13Grey(dst, dst);
204  }
205  }
206  if (ResetROI) dst.UnsetROI();
207  return 0;
208 }
209 
210 template <class InputStorageType, class OutputStorageType>
214  float threshold) {
215  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter3x3\n");
216  SetHalfWinSize(1, false);
217  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
218  _CalculateKernels(_GaussSigma, _GaussRatio);
219 
220  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
221  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
222  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
223  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
224 
225 #ifdef BIAS_DEBUG
226  if (src.GetChannelCount()!=1) {
227  BIASERR("image with multiple channels. stopping.");
228  BIASABORT;
229  }
230 #endif
231  if(!dst.SamePixelAndChannelCount(src)) {
232  if (!dst.IsEmpty()) dst.Release();
233  dst.Init(src.GetWidth(), src.GetHeight(), 1);
234  }
235  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
236  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
237  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
239  }
240 
241  // horizontal convolution ignoring wrap-arounds
242  const InputStorageType* pS = src.GetImageData();
243  CONV_FLOAT* pD = _tmpFloat.GetImageData();
244  // copy first byte
245  *pD++ = CONV_FLOAT(*pS);
246  const int width = src.GetWidth(), height = src.GetHeight();
247  const InputStorageType* pEndH = src.GetImageData() + width*height-1;
248  while (++pS<pEndH) {
249  if(pS[0] <= threshold ||
250  pS[-1] <= threshold ||
251  pS[1] <= threshold)
252  *pD++ = CONV_FLOAT(pS[0]);
253  else
254  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1]));
255  }
256  // copy last byte
257  *pD = CONV_FLOAT(*pS);
258 
259  // vertical convolution
260  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-1);
261  CONV_FLOAT* pI = _tmpFloat.GetImageData();
262  OutputStorageType* pO = dst.GetImageData();
263  const CONV_FLOAT* pEndFirst2 = _tmpFloat.GetImageData() + width;
264  const CONV_FLOAT* pEndLast2 = _tmpFloat.GetImageData() + width*height;
265  // copy first two lines
266  while (pI<pEndFirst2) *pO++ = OutputStorageType(*pI++);
267  // convolve main
268  pI--;
269  while (++pI<pEndV) {
270  if(pI[0] <= threshold ||
271  pI[width] <= threshold ||
272  pI[-width] <= threshold) {
273  *pO++ = OutputStorageType(pI[0]);
274  } else
275  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width]));
276  }
277  // copy last line
278  while (pI<pEndLast2) *pO++ = OutputStorageType(*pI++);
279 
280  // update roi
281  DecreaseROI(*src.GetROI(), *dst.GetROI(), 1);
282 
283  return 0;
284 }
285 
286 
287 template <class InputStorageType, class OutputStorageType>
291  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter3x3\n");
292  SetHalfWinSize(1, false);
293  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
294  _CalculateKernels(_GaussSigma, _GaussRatio);
295 
296  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
297  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
298  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
299  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
300 
301 #ifdef BIAS_DEBUG
302  if (src.GetChannelCount()!=1) {
303  BIASERR("image with multiple channels. stopping.");
304  BIASABORT;
305  }
306 #endif
307  if(!dst.SamePixelAndChannelCount(src)) {
308  if (!dst.IsEmpty()) dst.Release();
309  dst.Init(src.GetWidth(), src.GetHeight(), 1);
310  }
311  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
312  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
313  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
315  }
316 
317  // horizontal convolution ignoring wrap-arounds
318  const InputStorageType* pS = src.GetImageData();
319  CONV_FLOAT* pD = _tmpFloat.GetImageData();
320  // copy first byte
321  *pD++ = CONV_FLOAT(*pS);
322  const int width = src.GetWidth(), height = src.GetHeight();
323  const InputStorageType* pEndH = src.GetImageData() + width*height-1;
324  while (++pS<pEndH) {
325  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1]));
326  }
327  // copy last byte
328  *pD = CONV_FLOAT(*pS);
329 
330  // vertical convolution
331  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-1);
332  CONV_FLOAT* pI = _tmpFloat.GetImageData();
333  OutputStorageType* pO = dst.GetImageData();
334  const CONV_FLOAT* pEndFirst2 = _tmpFloat.GetImageData() + width;
335  const CONV_FLOAT* pEndLast2 = _tmpFloat.GetImageData() + width*height;
336  // copy first two lines
337  while (pI<pEndFirst2) *pO++ = OutputStorageType(*pI++);
338  // convolve main
339  pI--;
340  while (++pI<pEndV) {
341  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width]));
342  }
343  // copy last line
344  while (pI<pEndLast2) *pO++ = OutputStorageType(*pI++);
345 
346  // update roi
347  DecreaseROI(*src.GetROI(), *dst.GetROI(), 1);
348 
349 
350  return 0;
351 }
352 
353 template <class InputStorageType, class OutputStorageType>
357  float threshold) {
358  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter5x5\n");
359  SetHalfWinSize(2, false);
360  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
361  _CalculateKernels(_GaussSigma, _GaussRatio);
362 
363  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
364  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
365  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
366  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
367  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
368  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
369 
370 
371 #ifdef BIAS_DEBUG
372  if (src.GetChannelCount()!=1) {
373  BIASERR("image with multiple channels. stopping.");
374  BIASABORT;
375  }
376 #endif
377  if(!dst.SamePixelAndChannelCount(src)) {
378  if (!dst.IsEmpty()) dst.Release();
379  dst.Init(src.GetWidth(), src.GetHeight(), 1);
380  }
381  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
382  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
383  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
385  }
386 
387  // horizontal convolution ignoring wrap-arounds
388  const InputStorageType* pS = src.GetImageData();
389  CONV_FLOAT* pD = _tmpFloat.GetImageData();
390  // copy first two bytes
391  *pD++ = CONV_FLOAT(*pS++);
392  *pD++ = CONV_FLOAT(*pS);
393  const int width = src.GetWidth(), height = src.GetHeight();
394  const int width2 = 2*width;
395  const InputStorageType* pEndH = src.GetImageData() + width*height-2;
396  while (++pS<pEndH) {
397  if(pS[0] <= threshold ||
398  pS[-1] <= threshold ||
399  pS[1] <= threshold ||
400  pS[-2] <= threshold ||
401  pS[2] <= threshold )
402  *pD++ = CONV_FLOAT(pS[0]);
403  else
404  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1])+h2*(pS[-2]+pS[2]));
405  }
406  // copy last two bytes
407  *pD++ = CONV_FLOAT(*pS++);
408  *pD = CONV_FLOAT(*pS);
409 
410  // vertical convolution
411  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-2);
412  CONV_FLOAT* pI = _tmpFloat.GetImageData();
413  OutputStorageType* pO = dst.GetImageData();
414  const CONV_FLOAT* pEndFirst2 = _tmpFloat.GetImageData() + 2*width;
415  const CONV_FLOAT* pEndLast2 = _tmpFloat.GetImageData() + width*height;
416  // copy first two lines
417  while (pI<pEndFirst2) *pO++ = OutputStorageType(*pI++);
418  // convolve main
419  pI--;
420  while (++pI<pEndV) {
421  if(pI[0] <= threshold ||
422  pI[width] <= threshold ||
423  pI[-width] <= threshold ||
424  pI[width2] <= threshold ||
425  pI[-width2] <= threshold )
426  *pO++ = OutputStorageType(pI[0]);
427  else
428  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width])
429  + v2*(pI[-width2]+pI[width2]));
430 
431  }
432  // copy last two lines
433  while (pI<pEndLast2) *pO++ = OutputStorageType(*pI++);
434 
435  // update roi
436  DecreaseROI(*src.GetROI(), *dst.GetROI(), 2);
437 
438  return 0;
439 }
440 
441 template <class InputStorageType, class OutputStorageType>
445  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter5x5\n");
446  SetHalfWinSize(2, false);
447  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
448  _CalculateKernels(_GaussSigma, _GaussRatio);
449 
450  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
451  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
452  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
453  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
454  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
455  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
456 
457 
458 #ifdef BIAS_DEBUG
459  if (src.GetChannelCount()!=1) {
460  BIASERR("image with multiple channels. stopping.");
461  BIASABORT;
462  }
463 #endif
464  if(!dst.SamePixelAndChannelCount(src)) {
465  if (!dst.IsEmpty()) dst.Release();
466  dst.Init(src.GetWidth(), src.GetHeight(), 1);
467  }
468  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
469  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
470  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
472  }
473 
474  // horizontal convolution ignoring wrap-arounds
475  const InputStorageType* pS = src.GetImageData();
476  CONV_FLOAT* pD = _tmpFloat.GetImageData();
477  // copy first two bytes
478  *pD++ = CONV_FLOAT(*pS++);
479  *pD++ = CONV_FLOAT(*pS);
480  const int width = src.GetWidth(), height = src.GetHeight();
481  const int width2 = 2*width;
482  const InputStorageType* pEndH = src.GetImageData() + width*height-2;
483  while (++pS<pEndH) {
484  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1])+h2*(pS[-2]+pS[2]));
485  }
486  // copy last two bytes
487  *pD++ = CONV_FLOAT(*pS++);
488  *pD = CONV_FLOAT(*pS);
489 
490  // vertical convolution
491  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-2);
492  CONV_FLOAT* pI = _tmpFloat.GetImageData();
493  OutputStorageType* pO = dst.GetImageData();
494  const CONV_FLOAT* pEndFirst2 = _tmpFloat.GetImageData() + 2*width;
495  const CONV_FLOAT* pEndLast2 = _tmpFloat.GetImageData() + width*height;
496  // copy first two lines
497  while (pI<pEndFirst2) *pO++ = OutputStorageType(*pI++);
498  // convolve main
499  pI--;
500  while (++pI<pEndV) {
501  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width])
502  + v2*(pI[-width2]+pI[width2]));
503 
504  }
505  // copy last two lines
506  while (pI<pEndLast2) *pO++ = OutputStorageType(*pI++);
507 
508  // update roi
509  DecreaseROI(*src.GetROI(), *dst.GetROI(), 2);
510 
511  return 0;
512 }
513 
514 template <class InputStorageType, class OutputStorageType>
518 {
519  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter7x7\n");
520  SetHalfWinSize(3, false);
521  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
522  _CalculateKernels(_GaussSigma, _GaussRatio);
523 
524  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
525  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
526  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
527  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
528  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
529  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
530  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
531  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
532 
533 #ifdef BIAS_DEBUG
534  if (src.GetChannelCount()!=1) {
535  BIASERR("image with multiple channels. stopping.");
536  BIASABORT;
537  }
538 #endif
539  if(!dst.SamePixelAndChannelCount(src)) {
540  if (!dst.IsEmpty()) dst.Release();
541  dst.Init(src.GetWidth(), src.GetHeight(), 1);
542  }
543  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
544  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
545  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
547  }
548 
549  // horizontal convolution ignoring wrap-arounds
550  const InputStorageType* pS = src.GetImageData();
551  CONV_FLOAT* pD = _tmpFloat.GetImageData();
552  // copy first three bytes
553  *pD++ = CONV_FLOAT(*pS++);
554  *pD++ = CONV_FLOAT(*pS++);
555  *pD++ = CONV_FLOAT(*pS);
556  const int width = src.GetWidth(), height = src.GetHeight();
557  const int width2 = 2*width, width3=3*width;
558  const InputStorageType* pEndH = src.GetImageData() + width*height-3;
559  while (++pS<pEndH) {
560  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1])
561  +h2*(pS[-2]+pS[2])+h3*(pS[-3]+pS[3]));
562  }
563  // copy last three bytes
564  *pD++ = CONV_FLOAT(*pS++);
565  *pD++ = CONV_FLOAT(*pS++);
566  *pD = CONV_FLOAT(*pS);
567 
568  // vertical convolution
569  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-3);
570  CONV_FLOAT* pI = _tmpFloat.GetImageData();
571  OutputStorageType* pO = dst.GetImageData();
572  const CONV_FLOAT* pEndFirst3 = _tmpFloat.GetImageData() + 3*width;
573  const CONV_FLOAT* pEndLast3 = _tmpFloat.GetImageData() + width*height;
574  // copy first three lines
575  while (pI<pEndFirst3) *pO++ = OutputStorageType(*pI++);
576  // convolve main
577  pI--;
578  while (++pI<pEndV) {
579  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width])
580  + v2*(pI[-width2]+pI[width2])
581  + v3*(pI[-width3]+pI[width3]));
582  }
583  // copy last three lines
584  while (pI<pEndLast3) *pO++ = OutputStorageType(*pI++);
585 
586  // update roi
587  DecreaseROI(*src.GetROI(), *dst.GetROI(), 3);
588 
589  return 0;
590 }
591 
592 template <class InputStorageType, class OutputStorageType>
596 {
597  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter7x7\n");
598  SetHalfWinSize(3, false);
599  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
600  _CalculateKernels(_GaussSigma, _GaussRatio);
601 
602  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
603  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
604  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
605  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
606  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
607  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
608  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
609  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
610 
611 #ifdef BIAS_DEBUG
612  if ((src.GetChannelCount()!=3) || !src.IsInterleaved()) {
613  BIASERR("image is not 3-channel interleaved! stopping."
614  <<src.GetChannelCount()<<" "<<src.IsInterleaved());
615  BIASABORT;
616  }
617 #endif
618  if(!dst.SamePixelAndChannelCount(src)) {
619  if (!dst.IsEmpty()) dst.Release();
620  dst.Init(src.GetWidth(), src.GetHeight(), 3);
621  }
622  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
623  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
624  _tmpFloat.Init(src.GetWidth(),src.GetHeight(), 3,
626  }
627 
628  // horizontal convolution ignoring wrap-arounds
629  const InputStorageType* pS = src.GetImageData();
630  CONV_FLOAT* pD = _tmpFloat.GetImageData();
631  // copy first three bytes
632  *pD++ = CONV_FLOAT(*pS++);
633  *pD++ = CONV_FLOAT(*pS++);
634  *pD++ = CONV_FLOAT(*pS++);
635 
636  *pD++ = CONV_FLOAT(*pS++);
637  *pD++ = CONV_FLOAT(*pS++);
638  *pD++ = CONV_FLOAT(*pS++);
639 
640  *pD++ = CONV_FLOAT(*pS++);
641  *pD++ = CONV_FLOAT(*pS++);
642  *pD++ = CONV_FLOAT(*pS);
643 
644  const int width = src.GetWidth(), height = src.GetHeight();
645 
646  const int widthstep=3*width, widthstep2 = 6*width, widthstep3=9*width;
647  const InputStorageType* pEndH = src.GetImageData() + width*height*3-9;
648  while (++pS<pEndH) {
649  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-3]+pS[3])
650  +h2*(pS[-6]+pS[6])+h3*(pS[-9]+pS[9]));
651  }
652  // copy last three pixels
653  *pD++ = CONV_FLOAT(*pS++);
654  *pD++ = CONV_FLOAT(*pS++);
655  *pD++ = CONV_FLOAT(*pS++);
656  *pD++ = CONV_FLOAT(*pS++);
657  *pD++ = CONV_FLOAT(*pS++);
658  *pD++ = CONV_FLOAT(*pS++);
659  *pD++ = CONV_FLOAT(*pS++);
660  *pD++ = CONV_FLOAT(*pS++);
661  *pD = CONV_FLOAT(*pS);
662 
663  // vertical convolution
664  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-9);
665  CONV_FLOAT* pI = _tmpFloat.GetImageData();
666  OutputStorageType* pO = dst.GetImageData();
667  const CONV_FLOAT* pEndFirst3 = _tmpFloat.GetImageData() + 9*width;
668  const CONV_FLOAT* pEndLast3 = _tmpFloat.GetImageData() + 3*width*height;
669  // copy first three lines
670  while (pI<pEndFirst3) *pO++ = OutputStorageType(*pI++);
671  // convolve main
672  pI--;
673  while (++pI<pEndV) {
674  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-widthstep]+pI[widthstep])
675  + v2*(pI[-widthstep2]+pI[widthstep2])
676  + v3*(pI[-widthstep3]+pI[widthstep3]));
677  }
678  // copy last three lines
679  while (pI<pEndLast3) *pO++ = OutputStorageType(*pI++);
680 
681  // update roi
682  DecreaseROI(*src.GetROI(), *dst.GetROI(), 3);
683 
684  return 0;
685 }
686 
687 template <class InputStorageType, class OutputStorageType>
691 {
692  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter7x7\n");
693  SetHalfWinSize(3, false);
694  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
695  _CalculateKernels(_GaussSigma, _GaussRatio);
696 
697  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
698  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
699  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
700  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
701  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
702  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
703  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
704  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
705 
706 #ifdef BIAS_DEBUG
707  if ((src.GetChannelCount()!=4) || !src.IsInterleaved()) {
708  BIASERR("image is not 3-channel interleaved! stopping."
709  <<src.GetChannelCount()<<" "<<src.IsInterleaved());
710  BIASABORT;
711  }
712 #endif
713  if(!dst.SamePixelAndChannelCount(src)) {
714  if (!dst.IsEmpty()) dst.Release();
715  dst.Init(src.GetWidth(), src.GetHeight(), 4);
716  }
717  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
718  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
719  _tmpFloat.Init(src.GetWidth(),src.GetHeight(), 4,
721  }
722 
723  // horizontal convolution ignoring wrap-arounds
724  const InputStorageType* pS = src.GetImageData();
725  CONV_FLOAT* pD = _tmpFloat.GetImageData();
726  // copy first three pixels
727  *pD++ = CONV_FLOAT(*pS++);
728  *pD++ = CONV_FLOAT(*pS++);
729  *pD++ = CONV_FLOAT(*pS++);
730 
731  *pD++ = CONV_FLOAT(*pS++);
732  *pD++ = CONV_FLOAT(*pS++);
733  *pD++ = CONV_FLOAT(*pS++);
734 
735  *pD++ = CONV_FLOAT(*pS++);
736  *pD++ = CONV_FLOAT(*pS++);
737  *pD++ = CONV_FLOAT(*pS++);
738 
739  *pD++ = CONV_FLOAT(*pS++);
740  *pD++ = CONV_FLOAT(*pS++);
741  *pD++ = CONV_FLOAT(*pS);
742 
743  const int width = src.GetWidth(), height = src.GetHeight();
744 
745  const int widthstep=4*width, widthstep2 = 8*width, widthstep3=12*width;
746  const InputStorageType* pEndH = src.GetImageData() + width*height*4-12;
747  while (++pS<pEndH) {
748  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-4]+pS[4])
749  +h2*(pS[-8]+pS[8])+h3*(pS[-12]+pS[12]));
750  }
751  // copy last three pixels
752  *pD++ = CONV_FLOAT(*pS++);
753  *pD++ = CONV_FLOAT(*pS++);
754  *pD++ = CONV_FLOAT(*pS++);
755  *pD++ = CONV_FLOAT(*pS++);
756  *pD++ = CONV_FLOAT(*pS++);
757  *pD++ = CONV_FLOAT(*pS++);
758  *pD++ = CONV_FLOAT(*pS++);
759  *pD++ = CONV_FLOAT(*pS++);
760  *pD++ = CONV_FLOAT(*pS++);
761  *pD++ = CONV_FLOAT(*pS++);
762  *pD++ = CONV_FLOAT(*pS++);
763  *pD = CONV_FLOAT(*pS);
764 
765  // vertical convolution
766  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-12);
767  CONV_FLOAT* pI = _tmpFloat.GetImageData();
768  OutputStorageType* pO = dst.GetImageData();
769  const CONV_FLOAT* pEndFirst3 = _tmpFloat.GetImageData() + 12*width;
770  const CONV_FLOAT* pEndLast3 = _tmpFloat.GetImageData() + 4*width*height;
771  // copy first three lines
772  while (pI<pEndFirst3) *pO++ = OutputStorageType(*pI++);
773  // convolve main
774  pI--;
775  while (++pI<pEndV) {
776  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-widthstep]+pI[widthstep])
777  + v2*(pI[-widthstep2]+pI[widthstep2])
778  + v3*(pI[-widthstep3]+pI[widthstep3]));
779  }
780  // copy last three lines
781  while (pI<pEndLast3) *pO++ = OutputStorageType(*pI++);
782 
783  // update roi
784  DecreaseROI(*src.GetROI(), *dst.GetROI(), 3);
785 
786  return 0;
787 }
788 
789 
790 template <class InputStorageType, class OutputStorageType>
794 {
795  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter9x9\n");
796  SetHalfWinSize(4, false);
797  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
798  _CalculateKernels(_GaussSigma, _GaussRatio);
799  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[4]);
800  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
801  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
802  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
803  const FM_FLOAT h4 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
804  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[4]);
805  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
806  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
807  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
808  const FM_FLOAT v4 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
809 
810 #ifdef BIAS_DEBUG
811  if (src.GetChannelCount()!=1) {
812  BIASERR("image with multiple channels. stopping.");
813  BIASABORT;
814  }
815 #endif
816  if(!dst.SamePixelAndChannelCount(src)) {
817  if (!dst.IsEmpty()) dst.Release();
818  dst.Init(src.GetWidth(), src.GetHeight(), 1);
819  }
820  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
821  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
822  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
824  }
825 
826 
827  // horizontal convolution ignoring wrap-arounds
828  const InputStorageType* pS = src.GetImageData();
829  CONV_FLOAT* pD = _tmpFloat.GetImageData();
830  // copy first four bytes
831  *pD++ = CONV_FLOAT(*pS++);
832  *pD++ = CONV_FLOAT(*pS++);
833  *pD++ = CONV_FLOAT(*pS++);
834  *pD++ = CONV_FLOAT(*pS);
835 
836  const int width = src.GetWidth(), height = src.GetHeight();
837  const int width2 = 2*width, width3=3*width, width4=4*width;
838  const InputStorageType* pEndH = src.GetImageData() + width*height-4;
839  while (++pS<pEndH) {
840  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1])+h2*(pS[-2]+pS[2])
841  +h3*(pS[-3]+pS[3])+h4*(pS[-4]+pS[4]));
842  }
843  // copy last four bytes
844  *pD++ = CONV_FLOAT(*pS++);
845  *pD++ = CONV_FLOAT(*pS++);
846  *pD++ = CONV_FLOAT(*pS++);
847  *pD = CONV_FLOAT(*pS);
848 
849  // vertical convolution
850  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-4);
851  CONV_FLOAT* pI = _tmpFloat.GetImageData();
852  OutputStorageType* pO = dst.GetImageData();
853  const CONV_FLOAT* pEndFirst4 = _tmpFloat.GetImageData() + 4*width;
854  const CONV_FLOAT* pEndLast4 = _tmpFloat.GetImageData() + width*height;
855  // copy first four lines
856  while (pI<pEndFirst4) *pO++ = OutputStorageType(*pI++);
857  // convolve main
858  pI--;
859  while (++pI<pEndV) {
860  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width])
861  + v2*(pI[-width2]+pI[width2])
862  + v3*(pI[-width3]+pI[width3])
863  + v4*(pI[-width4]+pI[width4]));
864  }
865  // copy last four lines
866  while (pI<pEndLast4) *pO++ = OutputStorageType(*pI++);
867 
868  // update roi
869  DecreaseROI(*src.GetROI(), *dst.GetROI(), 4);
870 
871  return 0;
872 }
873 
874 template <class InputStorageType, class OutputStorageType>
878 {
879  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter11x11\n");
880  SetHalfWinSize(5, false);
881  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
882  _CalculateKernels(_GaussSigma, _GaussRatio);
883  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[5]);
884  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[4]);
885  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
886  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
887  const FM_FLOAT h4 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
888  const FM_FLOAT h5 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
889 
890  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[5]);
891  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[4]);
892  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
893  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
894  const FM_FLOAT v4 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
895  const FM_FLOAT v5 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
896 
897 #ifdef BIAS_DEBUG
898  if (src.GetChannelCount()!=1) {
899  BIASERR("image with multiple channels. stopping.");
900  BIASABORT;
901  }
902 #endif
903  if(!dst.SamePixelAndChannelCount(src)) {
904  if (!dst.IsEmpty()) dst.Release();
905  dst.Init(src.GetWidth(), src.GetHeight(), 1);
906  }
907  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
908  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
909  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
911  }
912 
913  // horizontal convolution ignoring wrap-arounds
914  const InputStorageType* pS = src.GetImageData();
915  CONV_FLOAT* pD = _tmpFloat.GetImageData();
916  // copy first five bytes
917  *pD++ = CONV_FLOAT(*pS++);
918  *pD++ = CONV_FLOAT(*pS++);
919  *pD++ = CONV_FLOAT(*pS++);
920  *pD++ = CONV_FLOAT(*pS++);
921  *pD++ = CONV_FLOAT(*pS);
922 
923  const int width = src.GetWidth(), height = src.GetHeight();
924  const int width2 = 2*width, width3=3*width, width4=4*width, width5=5*width;
925  const InputStorageType* pEndH = src.GetImageData() + width*height-5;
926  while (++pS<pEndH) {
927  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1])+h2*(pS[-2]+pS[2])
928  +h3*(pS[-3]+pS[3])+h4*(pS[-4]+pS[4])+h5*(pS[-5]+pS[5]));
929  }
930  // copy last five bytes
931  *pD++ = CONV_FLOAT(*pS++);
932  *pD++ = CONV_FLOAT(*pS++);
933  *pD++ = CONV_FLOAT(*pS++);
934  *pD++ = CONV_FLOAT(*pS++);
935  *pD = CONV_FLOAT(*pS);
936 
937  // vertical convolution
938  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-5);
939  CONV_FLOAT* pI = _tmpFloat.GetImageData();
940  OutputStorageType* pO = dst.GetImageData();
941  const CONV_FLOAT* pEndFirst5 = _tmpFloat.GetImageData() + 5*width;
942  const CONV_FLOAT* pEndLast5 = _tmpFloat.GetImageData() + width*height;
943  // copy first four lines
944  while (pI<pEndFirst5) *pO++ = OutputStorageType(*pI++);
945  // convolve main
946  pI--;
947  while (++pI<pEndV) {
948  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width])
949  + v2*(pI[-width2]+pI[width2])
950  + v3*(pI[-width3]+pI[width3])
951  + v4*(pI[-width4]+pI[width4])
952  + v5*(pI[-width5]+pI[width5]));
953  }
954  // copy last four lines
955  while (pI<pEndLast5) *pO++ = OutputStorageType(*pI++);
956 
957  // update roi
958  DecreaseROI(*src.GetROI(), *dst.GetROI(), 5);
959 
960  return 0;
961 }
962 
963 template <class InputStorageType, class OutputStorageType>
967 {
968  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter13x13\n");
969  SetHalfWinSize(6, false);
970  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
971  _CalculateKernels(_GaussSigma, _GaussRatio);
972  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[6]);
973  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[5]);
974  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[4]);
975  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
976  const FM_FLOAT h4 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
977  const FM_FLOAT h5 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
978  const FM_FLOAT h6 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
979 
980  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[6]);
981  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[5]);
982  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[4]);
983  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
984  const FM_FLOAT v4 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
985  const FM_FLOAT v5 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
986  const FM_FLOAT v6 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
987 
988 #ifdef BIAS_DEBUG
989  if (src.GetChannelCount()!=1) {
990  BIASERR("image with multiple channels. stopping.");
991  BIASABORT;
992  }
993 #endif
994  if(!dst.SamePixelAndChannelCount(src)) {
995  if (!dst.IsEmpty()) dst.Release();
996  dst.Init(src.GetWidth(), src.GetHeight(), 1);
997  }
998  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
999  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
1000  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
1002  }
1003 
1004  // horizontal convolution ignoring wrap-arounds
1005  const InputStorageType* pS = src.GetImageData();
1006  CONV_FLOAT* pD = _tmpFloat.GetImageData();
1007  // copy first five bytes
1008  *pD++ = CONV_FLOAT(*pS++);
1009  *pD++ = CONV_FLOAT(*pS++);
1010  *pD++ = CONV_FLOAT(*pS++);
1011  *pD++ = CONV_FLOAT(*pS++);
1012  *pD++ = CONV_FLOAT(*pS++);
1013  *pD++ = CONV_FLOAT(*pS);
1014 
1015  const int width = src.GetWidth(), height = src.GetHeight();
1016  const int width2 = 2*width, width3=3*width, width4=4*width,
1017  width5=5*width, width6=6*width;
1018  const InputStorageType* pEndH = src.GetImageData() + width*height-6;
1019  while (++pS<pEndH) {
1020  *pD++ = CONV_FLOAT(h0*pS[0]+h1*(pS[-1]+pS[1])+h2*(pS[-2]+pS[2])
1021  +h3*(pS[-3]+pS[3])+h4*(pS[-4]+pS[4])
1022  +h5*(pS[-5]+pS[5])+h6*(pS[-6]+pS[6]));
1023  }
1024  // copy last five bytes
1025  *pD++ = CONV_FLOAT(*pS++);
1026  *pD++ = CONV_FLOAT(*pS++);
1027  *pD++ = CONV_FLOAT(*pS++);
1028  *pD++ = CONV_FLOAT(*pS++);
1029  *pD++ = CONV_FLOAT(*pS++);
1030  *pD = CONV_FLOAT(*pS);
1031 
1032  // vertical convolution
1033  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-6);
1034  CONV_FLOAT* pI = _tmpFloat.GetImageData();
1035  OutputStorageType* pO = dst.GetImageData();
1036  const CONV_FLOAT* pEndFirst6 = _tmpFloat.GetImageData() + 6*width;
1037  const CONV_FLOAT* pEndLast6 = _tmpFloat.GetImageData() + width*height;
1038  // copy first four lines
1039  while (pI<pEndFirst6) *pO++ = OutputStorageType(*pI++);
1040  // convolve main
1041  pI--;
1042  while (++pI<pEndV) {
1043  *pO++ = OutputStorageType(v0*pI[0] + v1*(pI[-width]+pI[width])
1044  + v2*(pI[-width2]+pI[width2])
1045  + v3*(pI[-width3]+pI[width3])
1046  + v4*(pI[-width4]+pI[width4])
1047  + v5*(pI[-width5]+pI[width5])
1048  + v6*(pI[-width6]+pI[width6]));
1049  }
1050  // copy last four lines
1051  while (pI<pEndLast6) *pO++ = OutputStorageType(*pI++);
1052 
1053  // update roi
1054  DecreaseROI(*src.GetROI(), *dst.GetROI(), 6);
1055 
1056  return 0;
1057 }
1058 
1059 ////////////////////////////////////////////////////////////////////////////
1060 
1061 template <class InputStorageType, class OutputStorageType>
1063 _CalculateKernels(double Sigma, double Ratio)
1064 {
1065  const int hws = (int)(floor(sqrt(-2.0*Sigma*Sigma*log(Ratio))));
1066  const int size = (hws<<1)+1;
1067 
1068  Vector<FM_FLOAT> H(size), V(size);
1069  BIASCDOUT(D_CONV_KERNEL, "Sigma = "<<Sigma<<"\tRatio = "<<Ratio<<endl);
1070  BIASCDOUT(D_CONV_KERNEL, "size = "<<size<<"\thws = "<<hws<<endl);
1071 
1072  double sum=0.0;
1073  int i, ip, im;
1074  if (Sigma==0.0){
1075  H[0]=1.0; V[0]=1.0;
1076  } else {
1077  double fac=1.0/(2.0*Sigma*Sigma);
1078  register double g;
1079  // calculate the filter coefficients as gauss function
1080  for (i=0; i<=hws; i++){
1081  ip = hws+i;
1082  im = hws-i;
1083  g = expf((float)(-(double)(i*i)*fac));
1084  V[ip] = V[im] = H[im] = H[ip] = (FM_FLOAT)g;
1085  //BIASCDOUT(D_CONV_KERNEL, i<<"v:"<<this->_VVec[im]<<" "<<endl);
1086  sum += g*2.0;
1087  }
1088  sum -= V[hws]; // added twice
1089  BIASCDOUT(D_CONV_KERNEL, "sum = "<<sum<<endl);
1090 
1091  // normalize filter, so that the sum is 1.0
1092  sum = 1.0/sum;
1093  BIASCDOUT(D_CONV_KERNEL, "sum = "<<sum<<endl);
1094  for (i=0; i<=hws; i++){
1095  ip=hws+i;
1096  im=hws-i;
1097  H[ip] = H[im] = V[ip] = V[im] = V[im] * (FM_FLOAT)sum;
1098  BIASCDOUT(D_CONV_KERNEL, H[im]<<" "<<endl);
1099  }
1100  }
1101  _fm.Init(H,V);
1102  // compute the int approximation
1103  int thebits = _fm.ComputeIntPrecisionBits(sizeof(InputStorageType)*8,
1104  sizeof(CONV_INT)*8-1);
1105  _fm.CreateIntFilter(thebits,thebits,thebits);
1106 
1107  _LastSigma = Sigma;
1108  _LastRatio = Ratio;
1109 
1110 #ifdef BIAS_DEBUG
1111  if (this->DebugLevelIsSet(D_CONV_KERNEL)){
1112  sum = 0.0;
1113  for (i=-hws; i<=hws; i++){
1114  sum += H[i+hws];
1115  BIASCDOUT(D_CONV_KERNEL, H[i+hws]<<" ");
1116  }
1117  BIASCDOUT(D_CONV_KERNEL, " sum = "<<sum<<endl);
1118  }
1119 #endif
1120 
1121  //#define COMPARE
1122 #ifdef COMPARE
1123  fprintf(stdout, "gauss: \n");
1124  for (i=-hws; i<=hws; i++){
1125  fprintf(stdout, "%f ",HVec[i+hws]);
1126  }
1127  fprintf(stdout, "\n");
1128 #endif
1129 }
1130 
1131 #define MIN_NORM 0.7
1132 
1133 template <class InputStorageType, class OutputStorageType>
1137  const InputStorageType& thresh)
1138 {
1139  BIASCDOUT(D_FILTERBASE_CALLSTACK, "Gauss::Filter7x7\n");
1140  SetHalfWinSize(3, false);
1141  if ((_GaussSigma != _LastSigma)||(_GaussRatio != _LastRatio))
1142  _CalculateKernels(_GaussSigma, _GaussRatio);
1143 
1144  const FM_FLOAT h0 = (const FM_FLOAT)((*_fm.GetSepfh())[3]);
1145  const FM_FLOAT h1 = (const FM_FLOAT)((*_fm.GetSepfh())[2]);
1146  const FM_FLOAT h2 = (const FM_FLOAT)((*_fm.GetSepfh())[1]);
1147  const FM_FLOAT h3 = (const FM_FLOAT)((*_fm.GetSepfh())[0]);
1148  const FM_FLOAT v0 = (const FM_FLOAT)((*_fm.GetSepfv())[3]);
1149  const FM_FLOAT v1 = (const FM_FLOAT)((*_fm.GetSepfv())[2]);
1150  const FM_FLOAT v2 = (const FM_FLOAT)((*_fm.GetSepfv())[1]);
1151  const FM_FLOAT v3 = (const FM_FLOAT)((*_fm.GetSepfv())[0]);
1152 
1153 
1154  // cout<<"Kernel is "<<h0<<" "<<h1<<" "<<h2<<" "<<h3<<" "<<v0<<" "<<v1<<" "<<v2<<" "<<v3<<" "<<endl;
1155 
1156 #ifdef BIAS_DEBUG
1157  if (src.GetChannelCount()!=1) {
1158  BIASERR("image with multiple channels. stopping.");
1159  BIASABORT;
1160  }
1161 #endif
1162  if(!dst.SamePixelAndChannelCount(src)) {
1163  if (!dst.IsEmpty()) dst.Release();
1164  dst.Init(src.GetWidth(), src.GetHeight(), 1);
1165  }
1166  if(!_tmpFloat.SamePixelAndChannelCount(src)) {
1167  if (!_tmpFloat.IsEmpty()) _tmpFloat.Release();
1168  _tmpFloat.Init(src.GetWidth(),src.GetHeight(),1,
1170  }
1171 
1172  // horizontal convolution ignoring wrap-arounds
1173  const InputStorageType* pS = src.GetImageData();
1174  CONV_FLOAT* pD = _tmpFloat.GetImageData();
1175  // copy first three bytes
1176  *pD++ = CONV_FLOAT(thresh);
1177  pS++;
1178  *pD++ = CONV_FLOAT(thresh);
1179  pS++;
1180  *pD++ = CONV_FLOAT(thresh);
1181 
1182  const int width = src.GetWidth(), height = src.GetHeight();
1183  const int width2 = 2*width, width3=3*width;
1184  const InputStorageType* pEndH = src.GetImageData() + width*height-3;
1185 
1186  while (++pS<pEndH) {
1187  register CONV_FLOAT normalization = 0;
1188  register CONV_FLOAT accumulator = 0;
1189  if (pS[-3]>thresh) {
1190  normalization += h3;
1191  accumulator += h3*CONV_FLOAT(pS[-3]);
1192  }
1193  if (pS[-2]>thresh) {
1194  normalization += h2;
1195  accumulator += h2*CONV_FLOAT(pS[-2]);
1196  }
1197  if (pS[-1]>thresh) {
1198  normalization += h1;
1199  accumulator += h1*CONV_FLOAT(pS[-1]);
1200  }
1201  if (pS[0]>thresh) {
1202  normalization += h0;
1203  accumulator += h0*CONV_FLOAT(pS[0]);
1204  }
1205  if (pS[1]>thresh) {
1206  normalization += h1;
1207  accumulator += h1*CONV_FLOAT(pS[1]);
1208  }
1209 
1210  if (pS[2]>thresh) {
1211  normalization += h2;
1212  accumulator += h2*CONV_FLOAT(pS[2]);
1213  }
1214  if (pS[3]>thresh) {
1215  normalization += h3;
1216  accumulator += h3*CONV_FLOAT(pS[3]);
1217  }
1218  if (normalization>MIN_NORM)
1219  *pD++ = CONV_FLOAT(accumulator / normalization);
1220  else
1221  *pD++ = CONV_FLOAT(thresh);
1222  }
1223  // copy last three bytes
1224  *pD++ = CONV_FLOAT(thresh);pS++;
1225  *pD++ = CONV_FLOAT(thresh);pS++;
1226  *pD = CONV_FLOAT(thresh);pS++;
1227 
1228  const CONV_FLOAT cfthresh = CONV_FLOAT(thresh);
1229 
1230  // vertical convolution
1231  const CONV_FLOAT* pEndV = _tmpFloat.GetImageData() + width*(height-3);
1232  CONV_FLOAT* pI = _tmpFloat.GetImageData();
1233  OutputStorageType* pO = dst.GetImageData();
1234  const CONV_FLOAT* pEndFirst3 = _tmpFloat.GetImageData() + 3*width;
1235  const CONV_FLOAT* pEndLast3 = _tmpFloat.GetImageData() + width*height;
1236  // copy first three lines
1237  while (pI<pEndFirst3) {
1238  *pO++ = OutputStorageType(thresh);
1239  pI++;
1240  }
1241  // convolve main
1242  pI--;
1243  while (++pI<pEndV) {
1244  register CONV_FLOAT normalization = 0;
1245  register CONV_FLOAT accumulator = 0;
1246 
1247  if (pI[-width3]>cfthresh) {
1248  normalization += v3;
1249  accumulator += v3*CONV_FLOAT(pI[-width3]);
1250  }
1251  if (pI[-width2]>cfthresh) {
1252  normalization += v2;
1253  accumulator += v2*CONV_FLOAT(pI[-width2]);
1254  }
1255  if (pI[-width]>cfthresh) {
1256  normalization += v1;
1257  accumulator += v1*CONV_FLOAT(pI[-width]);
1258  }
1259  if (pI[0]>cfthresh) {
1260  normalization += v0;
1261  accumulator += v0*CONV_FLOAT(pI[0]);
1262  }
1263  if (pI[width]>cfthresh) {
1264  normalization += v1;
1265  accumulator += v1*CONV_FLOAT(pI[width]);
1266  }
1267 
1268  if (pI[width2]>cfthresh) {
1269  normalization += v2;
1270  accumulator += v2*CONV_FLOAT(pI[width2]);
1271  }
1272  if (pI[width3]>cfthresh) {
1273  normalization += v3;
1274  accumulator += v3*CONV_FLOAT(pI[width3]);
1275  }
1276 
1277  if (normalization>MIN_NORM)
1278  *pO++ = OutputStorageType(accumulator / normalization);
1279  else
1280  *pO++ = OutputStorageType(thresh);
1281 
1282  }
1283  // copy last three lines
1284  while (pI<pEndLast3) {
1285  *pO++ = OutputStorageType(thresh);
1286  pI++;
1287  }
1288 
1289  // update roi
1290  DecreaseROI(*src.GetROI(), *dst.GetROI(), 3);
1291 
1292  return 0;
1293 }
1294 
1295 
1296 
1297 
1298 //////////////////////////////////////////////////////////////////////////
1299 // instantiation
1300 //////////////////////////////////////////////////////////////////////////
1301 namespace BIAS{
1302 #define FILTER_INSTANTIATION_CLASS Gauss
1303 #include "Filterinst.hh"
1304 }
void Release()
reimplemented from ImageBase
Definition: Image.cpp:1579
double _GaussSigma
the parameter sigma of gaussian kernel
Definition: Gauss.hh:189
virtual int Filter13x13Grey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:965
virtual int Filter7x7Grey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:516
virtual int Filter11x11Grey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:876
class for handling different region of interest (ROI) representations...
Definition: ROI.hh:118
Gauss< InputStorageType, OutputStorageType > * _GaussIO
Definition: Gauss.hh:212
generic convolution class.
Definition: Convolution.hh:66
int SetCorners(unsigned UpperLeftX, unsigned UpperLeftY, unsigned LowerRightX, unsigned LowerRightY)
Sets a rectangular region of interest.
Definition: ROI.cpp:287
void GetCorners(unsigned &UpperLeftX, unsigned &UpperLeftY, unsigned &LowerRightX, unsigned &LowerRightY) const
Return the region of interest, by saving the coordinates within the variables defined by the paramete...
Definition: ROI.hh:443
virtual int FilterFloat(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
do the convolution using floating point calculations
bool IsEmpty() const
check if ImageData_ points to allocated image buffer or not
Definition: ImageBase.hh:245
virtual int FilterGrey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst, bool ResetROI=false)
wrapper around FilterNxNGrey, which decides on sigma and cutoffratio which filter to apply and how of...
Definition: Gauss.cpp:126
bool IsInterleaved() const
Definition: ImageBase.hh:491
int Filter7x7GreyIgnoreBelowThreshold(const Image< InputStorageType > &src, Image< OutputStorageType > &dst, const InputStorageType &thresh)
7x7 gauss filtering, values below threshold are ignored useful for depth map filtering ...
Definition: Gauss.cpp:1135
virtual int Filter(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
decides on the image types which FilterFunction should be called
virtual int Filter3x3GreyThreshold(const Image< InputStorageType > &src, Image< OutputStorageType > &dst, float threshold=0.0)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored Keep src ...
Definition: Gauss.cpp:212
double _GaussRatio
minimum ratio of 1D kernel center and border, ignore smaller entries
Definition: Gauss.hh:197
float image storage type
Definition: ImageBase.hh:118
virtual int Filter3x3Grey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:289
unsigned int GetWidth() const
Definition: ImageBase.hh:312
virtual int Filter7x7RGB(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:594
smoothing with gaussian kernel
Definition: Gauss.hh:51
ROI * GetROI()
Returns a pointer to the roi object.
Definition: ImageBase.hh:615
virtual int Filter5x5GreyThreshold(const Image< InputStorageType > &src, Image< OutputStorageType > &dst, float threshold=0.0)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored...
Definition: Gauss.cpp:355
void _CalculateKernels(double Sigma, double Ratio)
calculates the kernel
Definition: Gauss.cpp:1063
virtual Gauss< InputStorageType, OutputStorageType > * Clone() const
Definition: Gauss.hh:157
static int ConvertST(const BIAS::ImageBase &source, BIAS::ImageBase &dest, ImageBase::EStorageType targetST)
Function to convert the storage type of images e.g.
unsigned int GetChannelCount() const
returns the number of Color channels, e.g.
Definition: ImageBase.hh:382
virtual int FilterInt(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
sets gauss kernel if params changed and calls convolution
Definition: Gauss.cpp:115
unsigned int GetHeight() const
Definition: ImageBase.hh:319
bool SamePixelAndChannelCount(const ImageBase &Image) const
checks if data area has same &quot;size&quot; as Image of other type
Definition: ImageBase.hh:73
Gauss< OutputStorageType, OutputStorageType > * _GaussOO
to allow for iterated gauss convolution saving the intermediate image we need a different instance...
Definition: Gauss.hh:211
virtual int Filter(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
sets gauss kernel if params changed and calls convolution or fast grey implementation if possible ...
Definition: Gauss.cpp:89
void Init(unsigned int Width, unsigned int Height, unsigned int channels=1, enum EStorageType storageType=ST_unsignedchar, const bool interleaved=true)
calls Init from ImageBase storageType is ignored, just dummy argument
Definition: Image.cpp:421
const StorageType * GetImageData() const
overloaded GetImageData() from ImageBase
Definition: Image.hh:137
virtual int Filter7x7RGBA(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:689
double _LastRatio
Definition: Gauss.hh:200
double _LastSigma
Definition: Gauss.hh:200
enum EStorageType GetStorageType() const
Definition: ImageBase.hh:414
virtual int Filter9x9Grey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:792
virtual int Filter5x5Grey(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
fast, approximate, direct implementation, ignoring wrap-arounds, roi is updated but ignored ...
Definition: Gauss.cpp:443
virtual int FilterFloat(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
sets gauss kernel if params changed and calls convolution
Definition: Gauss.cpp:105
virtual int FilterInt(const Image< InputStorageType > &src, Image< OutputStorageType > &dst)
do the convolution using integer arithmetics and shifts
void UnsetROI()
deprecated, use GetROI()-&gt;UnsetROI()
Definition: ImageBase.cpp:1057