Basic Image AlgorithmS Library  2.8.0
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
clfTVL1Flow.cpp
1 /*
2  This file is part of the BIAS library (Basic ImageAlgorithmS).
3 
4  Copyright (C) 2003, 2004 (see file CONTACTS 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 <OpenCLFramework/Algorithm/clfTVL1Flow.hh>
26 
27 //#define WRITE_DBG_IMAGES
28 #undef WRITE_DBG_IMAGES
29 
30 using namespace std;
31 using namespace BIAS;
32 using namespace BIAS;
33 
34 clfTVL1Flow::clfTVL1Flow(clfContext *ctx) {
35 
36  { // default parameter set
37  sigma_ = 1.0f/sqrt(8.0f); // was 1.0f/sqrt(8.0f);
38  tau_ = 1.0f/sqrt(8.0f); // 1.0f/sqrt(8.0f);
39  lambda_ = 20; // was 50
40  gamma_ = 0.02f; //was 0.02f
41  warps_ = 1; //was 1
42  maxIterations_ = 20; // was 10
43  }
44 
45  { // initialize all the cl stuff
46  context_ = ctx;
47  maxComputeUnits_ = context_->GetDeviceInfo().maxComputeUnits;
48 
49  colorConvCL_ = new clfColorConversion(context_);
50  downsampleCL_ = new clfSimpleFilter(context_, "DownsampleBy2");
51  upsampleCL_ = new clfSimpleFilter(context_, "UpsampleBy2");
52 
53  clearFilter_ = new clfFilter<float, float>(context_);
54 
55  programCL_ = context_->CreateProgram();
56 
57  // add code to our program
58  programCL_->AddSource("/OpenCLFramework/Algorithm/cl/tvl1flow.cl");
59  programCL_->Build();
60  programCL_->AddKernel("Warp");
61  programCL_->AddKernel("UvwGrad");
62  programCL_->AddKernel("UpdateP");
63  programCL_->AddKernel("UpdateDiv");
64  programCL_->AddKernel("UpdateUVW");
65  programCL_->AddKernel("UpdateUVWold");
66  programCL_->AddKernel("WeightedGradX");
67  programCL_->AddKernel("WeightedGradY");
68  programCL_->AddKernel("CreateColorCoded");
69  }
70  { // initialize base images and variables
71  numIm_ = 0;
72  image1_ = context_->CreateImage2D();
73  image2_ = context_->CreateImage2D();
74  debug_ = context_->CreateImage2D();
75  levels_ = 0;
76  }
77 }
78 
79 clfTVL1Flow::~clfTVL1Flow() {
80  // free all the memory
81  for (unsigned int i=0;i<pyramid1_.size();i++) {
82  delete pyramid1_[i];
83  delete pyramid2_[i];
84  delete uvw_[i];
85  delete uvwX_[i];
86  delete uvwY_[i];
87  delete uvw0_[i];
88  delete uvwOld_[i];
89  delete image2warpedt_[i];
90  delete imageDiv_[i];
91  delete P1_[i];
92  delete P2_[i];
93  delete imageCopy_[i];
94  delete imageCopy2_[i];
95  delete Ixyg_[i];
96  }
97  delete programCL_;
98  delete colorConvCL_;
99  delete downsampleCL_;
100  delete upsampleCL_;
101  delete clearFilter_;
102 
103  delete image1_;
104  delete image2_;
105  delete debug_;
106 }
107 
108 void clfTVL1Flow::SetNumberOfIterations(unsigned int iternumber)
109 {
110  maxIterations_ = iternumber;
111 }
112 
113 void clfTVL1Flow::SetSmoothness(const float& smoothness)
114 {
115  lambda_ = smoothness;
116 }
117 
118 void clfTVL1Flow::GetUVW(BIAS::Image<float> &image) {
119  uvw_[0]->CopyToBiasImage(image);
120 }
121 
122 void clfTVL1Flow::GetColorCoded(BIAS::Image<float> &image) {
123  // run the color coding kernel on biggest uvw image
124  programCL_->KernelSetArgument("CreateColorCoded", 0, *uvw_[0]);
125  programCL_->KernelSetArgument("CreateColorCoded", 1, *imageCopy_[0]);
126  unsigned int w,h;
127  imageCopy_[0]->GetImageDim(w,h);
128  context_->RunOn2DRange(*programCL_, "CreateColorCoded", w,h,16,16);
129  imageCopy_[0]->CopyToBiasImage(image);
130 }
131 
132 
133 void clfTVL1Flow::AddImage(BIAS::Image<unsigned char> image) {
134  if (numIm_ == 0) {
135  // on first image added, reserve all the memory needed
136  image1_->AllocateFromBiasImage(image, true, false, true);
137  image2_->AllocateFromBiasImage(image, true, false, true);
138  int curWidth = image.GetWidth();
139  int curHeight = image.GetHeight();
140  while (curHeight % 16 == 0) {
141  pyramid1_.push_back(newImage_(curWidth, curHeight));
142  pyramid2_.push_back(newImage_(curWidth, curHeight));
143  uvw_.push_back(newImage_(curWidth, curHeight,true));
144  uvwX_.push_back(newImage_(curWidth, curHeight));
145  uvwY_.push_back(newImage_(curWidth, curHeight));
146  uvw0_.push_back(newImage_(curWidth, curHeight));
147  uvwOld_.push_back(newImage_(curWidth, curHeight));
148  image2warpedt_.push_back(newImage_(curWidth, curHeight));
149  imageDiv_.push_back(newImage_(curWidth, curHeight));
150  P1_.push_back(newImage_(curWidth, curHeight,true));
151  P2_.push_back(newImage_(curWidth, curHeight,true));
152  Ixyg_.push_back(newImage_(curWidth, curHeight));
153  imageCopy_.push_back(newImage_(curWidth, curHeight));
154  imageCopy2_.push_back(newImage_(curWidth, curHeight));
155  curWidth /= 2;
156  curHeight /= 2;
157  }
158  numIm_ = 1;
159  } else if (numIm_ == 1) {
160  // on second image, store the new image data
161  image2_->WriteToImage( image.GetImageData() );
162  // convert both images to grey, float, normalized by 1/255
163  if (image.GetColorModel() == ImageBase::CM_RGB) {
164  colorConvCL_->UCharRGBToFloatGreyNormalized( image1_, pyramid1_[0]);
165  colorConvCL_->UCharRGBToFloatGreyNormalized( image2_, pyramid2_[0]);
166  } else {
167  colorConvCL_->UCharGreyToFloatGrey( image1_, pyramid1_[0], true);
168  colorConvCL_->UCharGreyToFloatGrey( image2_, pyramid2_[0], true);
169  }
170  // create pyramids
171  CreatePyramid_(pyramid1_);
172  CreatePyramid_(pyramid2_);
173  levels_ = (unsigned int)pyramid1_.size();
174  numIm_ = 2;
175  } else {
176  // on every following image
177  // switch pointers
178  clfImage2D *tmp = image1_;
179  image1_ = image2_;
180  image2_ = tmp;
181  for (unsigned int i=0;i<levels_;i++) {
182  tmp = pyramid1_[i];
183  pyramid1_[i] = pyramid2_[i];
184  pyramid2_[i] = tmp;
185  }
186  // upload new image and convert
187  image2_->WriteToImage( image.GetImageData() );
188  if (image.GetColorModel() == ImageBase::CM_RGB) {
189  colorConvCL_->UCharRGBToFloatGreyNormalized( image2_, pyramid2_[0]);
190  } else {
191  colorConvCL_->UCharGreyToFloatGrey( image1_, pyramid1_[0], true);
192  }
193  // create new pyramid
194  CreatePyramid_(pyramid2_);
195  // clear necessary lowest level images
196  clearFilter_->Clear(0.0f, uvw_[levels_-1]);
197  clearFilter_->Clear(0.0f, P1_[levels_-1]);
198  clearFilter_->Clear(0.0f, P2_[levels_-1]);
199  numIm_++;
200  }
201 }
202 
203 clfImage2D* clfTVL1Flow::newImage_(unsigned int w, unsigned int h, bool clear) {
204  clfImage2D* pim = context_->CreateImage2D();
205  pim->Allocate(ImageBase::ST_float, ImageBase::CM_RGBA, w, h, 0, false, false);
206  if (clear) {
207  clearFilter_->Clear(0.0f, pim);
208  }
209  return pim;
210 }
211 
212 
213 void clfTVL1Flow::CreatePyramid_(std::vector<clfImage2D *> &pyramid) {
214  for (unsigned int i=1;i<pyramid.size();i++) {
215  downsampleCL_->FilterScale( pyramid[i-1], pyramid[i], 0.5f);
216  }
217 }
218 
219 void clfTVL1Flow::Upsample_(clfImage2D *src, clfImage2D *dst, float scale) {
220  upsampleCL_->GetProgram()->KernelSetArgument("UpsampleBy2", 2, scale);
221  upsampleCL_->Filter(src,dst);
222 }
223 
224 void clfTVL1Flow::Compute() {
225  if (numIm_ < 2) return;
226  // from bottom to top in pyramid
227  for (int level = (int)pyramid1_.size() -1;level>=0;level--) {
228  if (level != (int) pyramid1_.size() -1) {
229  // upsample results from previous level to current level
230  Upsample_(uvw_[level+1], uvw_[level], 2.0f);
231  Upsample_(P1_[level+1], P1_[level]);
232  Upsample_(P2_[level+1], P2_[level]);
233  }
234  // compute flow for this level
235  ComputeFlow_(level);
236 
237 #ifdef WRITE_DBG_IMAGES
238  Image<float> debug;
239  pyramid1_[level]->CopyToBiasImage(debug);
240  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-im0.mip", debug);
241  pyramid2_[level]->CopyToBiasImage(debug);
242  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-im1.mip", debug);
243  uvw_[level]->CopyToBiasImage(debug);
244  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-uvw.mip", debug);
245  uvwX_[level]->CopyToBiasImage(debug);
246  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-uvwX.mip", debug);
247  uvwY_[level]->CopyToBiasImage(debug);
248  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-uvwY.mip", debug);
249  image2warpedt_[level]->CopyToBiasImage(debug);
250  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-image2warpedt_.mip", debug);
251  Ixyg_[level]->CopyToBiasImage(debug);
252  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-Ixyg_.mip", debug);
253  P1_[level]->CopyToBiasImage(debug);
254  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-p1_.mip", debug);
255  P2_[level]->CopyToBiasImage(debug);
256  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-p2_.mip", debug);
257  imageDiv_[level]->CopyToBiasImage(debug);
258  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-div_.mip", debug);
259  Ixyg_[level]->CopyToBiasImage(debug);
260  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-grad.mip", debug);
261 #endif // WRITE_DBG_IMAGES
262  }
263 }
264 
265 void clfTVL1Flow::Warp_(unsigned int level) {
266  clfImage2D *Ixyg = Ixyg_[level];
267  clfImage2D *image1 = pyramid1_[level];
268  clfImage2D *image2 = pyramid2_[level];
269  clfImage2D *uvw = uvw_[level];
270  clfImage2D *warped = image2warpedt_[level];
271  clfImage2D *tmp = imageCopy_[level];
272 
273  unsigned int w,h;
274  image1->GetImageDim(w,h);
275 
276  // create gradient image in x dir
277  programCL_->KernelSetArgument("WeightedGradX", 0, *image1);
278  programCL_->KernelSetArgument("WeightedGradX", 1, *tmp);
279  context_->RunOn2DRange(*programCL_, "WeightedGradX", w,h, 16,16);
280 
281  // create gradient image in y dir
282  // the kernel will also compute squared grad
283  programCL_->KernelSetArgument("WeightedGradY", 0, *image1);
284  programCL_->KernelSetArgument("WeightedGradY", 1, *tmp);
285  programCL_->KernelSetArgument("WeightedGradY", 2, *Ixyg);
286  programCL_->KernelSetArgument("WeightedGradY", 3, gamma_);
287  context_->RunOn2DRange(*programCL_, "WeightedGradY", w,h, 16,16);
288 
289  // warp images, kernel will create difference image (It in channel 3 of warped image
290  programCL_->KernelSetArgument("Warp", 0, *image1); // read
291  programCL_->KernelSetArgument("Warp", 1, *image2); // read
292  programCL_->KernelSetArgument("Warp", 2, *uvw); // read
293  programCL_->KernelSetArgument("Warp", 3, *warped); // write
294  context_->RunOn2DRange(*programCL_, "Warp", w,h,16,16);
295 
296 }
297 
298 void clfTVL1Flow::ComputeFlow_(unsigned int level) {
299  // copy pointers for more readability
300  clfImage2D *uvw = uvw_[level];
301  clfImage2D *uvwold = uvwOld_[level];
302  clfImage2D *uvw0 = uvw0_[level];
303  clfImage2D *Ixyg = Ixyg_[level];
304  clfImage2D *uvwX = uvwX_[level];
305  clfImage2D *uvwY = uvwY_[level];
306  clfImage2D *p1 = P1_[level];
307  clfImage2D *p2 = P2_[level];
308  clfImage2D *div = imageDiv_[level];
309  clfImage2D *tmp = imageCopy_[level];
310  clfImage2D *tmp2 = imageCopy2_[level];
311  clfImage2D *warped = image2warpedt_[level];
312 
313  // set all the parameters
314  programCL_->KernelSetArgument("UvwGrad", 0, *uvwold); // read
315  programCL_->KernelSetArgument("UvwGrad", 1, *uvwX); // write
316  programCL_->KernelSetArgument("UvwGrad", 2, *uvwY); // write
317 
318  programCL_->KernelSetArgument("UpdateP", 0, *tmp); //read
319  programCL_->KernelSetArgument("UpdateP", 1, *tmp2); //read
320  programCL_->KernelSetArgument("UpdateP", 2, *uvwX); //read
321  programCL_->KernelSetArgument("UpdateP", 3, *uvwY); //read
322  programCL_->KernelSetArgument("UpdateP", 4, *p1); //write
323  programCL_->KernelSetArgument("UpdateP", 5, *p2); //write
324  programCL_->KernelSetArgument("UpdateP", 6, sigma_);
325 
326  programCL_->KernelSetArgument("UpdateDiv", 0, *p1); //read
327  programCL_->KernelSetArgument("UpdateDiv", 1, *p2); //read
328  programCL_->KernelSetArgument("UpdateDiv", 2, *div);//write
329 
330  programCL_->KernelSetArgument("UpdateUVW", 0, *uvwold); // read
331  programCL_->KernelSetArgument("UpdateUVW", 1, *div); //read
332  programCL_->KernelSetArgument("UpdateUVW", 2, *uvw0); //read
333  programCL_->KernelSetArgument("UpdateUVW", 3, *Ixyg); //read
334  programCL_->KernelSetArgument("UpdateUVW", 4, *warped); //read
335 
336  programCL_->KernelSetArgument("UpdateUVW", 5, *uvw); // write
337  programCL_->KernelSetArgument("UpdateUVW", 6, tau_);
338  programCL_->KernelSetArgument("UpdateUVW", 7, lambda_);
339  programCL_->KernelSetArgument("UpdateUVW", 8, gamma_);
340 
341  programCL_->KernelSetArgument("UpdateUVWold", 0, *tmp);
342  programCL_->KernelSetArgument("UpdateUVWold", 1, *uvw);
343  programCL_->KernelSetArgument("UpdateUVWold", 2, *uvwold);
344 
345  // get pyramid level size
346  unsigned int w,h;
347  uvw->GetImageDim(w,h);
348 
349  // backup uvw
350  uvw->CopyToImage(*uvwold);
351 
352  // default value is to do one warp only
353  for (unsigned int j=0;j<warps_;j++) {
354  // keep state before iterating
355  uvw->CopyToImage(*uvw0);
356  // warp image2 according to previous result
357  Warp_(level);
358 
359  for (unsigned int k=0;k<maxIterations_;k++) {
360  // gradients of direction
361  context_->RunOn2DRange(*programCL_, "UvwGrad", w,h,16,16);
362 
363  // update statistics
364  p1->CopyToImage( *tmp );
365  p2->CopyToImage( *tmp2 );
366  context_->RunOn2DRange(*programCL_, "UpdateP", w,h,16,16);
367  // update div
368  context_->RunOn2DRange(*programCL_, "UpdateDiv", w,h,16,16);
369 
370  uvw->CopyToImage(*uvwold);
371  // generate new uvw
372  context_->RunOn2DRange(*programCL_, "UpdateUVW", w,h,16,16);
373  // update for next iteration
374  uvwold->CopyToImage(*tmp);
375  context_->RunOn2DRange(*programCL_, "UpdateUVWold", w,h, 16,16);
376 
377 #ifdef WRITE_DBG_IMAGES
378  if (level==4) {
379  Image<float> debug;
380  uvwold->CopyToBiasImage(debug);
381  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-uvwold.mip", debug);
382  uvw->CopyToBiasImage(debug);
383  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-uvw.mip", debug);
384  p1->CopyToBiasImage(debug);
385  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-p1.mip", debug);
386  p2->CopyToBiasImage(debug);
387  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-p2.mip", debug);
388  div->CopyToBiasImage(debug);
389  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-div.mip", debug);
390  Ixyg->CopyToBiasImage(debug);
391  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-Ixyg.mip", debug);
392  uvwX->CopyToBiasImage(debug);
393  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-uvwX.mip", debug);
394  uvwY->CopyToBiasImage(debug);
395  ImageIO::Save("cbg-"+FileHandling::LeadingZeroString(level,2)+"-"+FileHandling::LeadingZeroString(k,2)+"-uvwY.mip", debug);
396  }
397 #endif // WRITE_DBG_IMAGES
398  }
399  }
400 }
401 
BIAS::clfImage2D::GetImageDim
void GetImageDim(unsigned int &width, unsigned int &height) const
Definition: clfImage2D.hh:110
BIAS::clfImage2D::Allocate
void Allocate(BIAS::ImageBase::EStorageType st, BIAS::ImageBase::EColorModel cm, unsigned int width, unsigned int height, unsigned int stride=0, bool readonly=false, bool writeonly=false, const void *hostptr=NULL, bool copy=false)
Allocation of a memory buffer as 2D image, either call directly or use wrapper for BIAS::ImageBase...
Definition: clfImage2D.cpp:62
BIAS::clfImage2D
OpenCL Image2D wrapper.
Definition: clfImage2D.hh:46
BIAS::ImageBase::GetWidth
unsigned int GetWidth() const
Definition: ImageBase.hh:312
BIAS::clfContext
OpenCL Context wrapper.
Definition: clfContext.hh:49
BIAS::clfContext::GetDeviceInfo
clfDeviceInfo GetDeviceInfo(unsigned int device=0)
Definition: clfContext.hh:141
BIAS::ImageBase::GetHeight
unsigned int GetHeight() const
Definition: ImageBase.hh:319
BIAS::Image::GetImageData
const StorageType * GetImageData() const
overloaded GetImageData() from ImageBase
Definition: Image.hh:137
BIAS::ImageBase::GetColorModel
enum EColorModel GetColorModel() const
Definition: ImageBase.hh:407
BIAS::clfImage2D::CopyToBiasImage
void CopyToBiasImage(BIAS::ImageBase &image, unsigned int originX=0, unsigned int originY=0, unsigned int regionX=0, unsigned int regionY=0)
Definition: clfImage2D.cpp:283
BIAS::clfImage2D::WriteToImage
void WriteToImage(const void *data, unsigned int originX=0, unsigned int originY=0, unsigned int regionX=0, unsigned int regionY=0, bool blocking=true)
write from host memory to image
Definition: clfImage2D.cpp:195
BIAS::clfImage2D::CopyToImage
void CopyToImage(clfImage2D &outputimage, unsigned int srcoriginX=0, unsigned int srcoriginY=0, unsigned int dstoriginX=0, unsigned int dstoriginY=0, unsigned int regionX=0, unsigned int regionY=0, bool blocking=true)
Definition: clfImage2D.cpp:233
BIAS::clfFilter::Build
virtual int Build(unsigned int size)
Definition: clfFilter.cpp:87
Generated on Tue Feb 6 2018 02:01:42 for Basic Image AlgorithmS Library by   doxygen 1.8.5