Visual Servoing Platform  version 3.3.1 under development (2020-08-12)
vpImageTools.cpp
1 /****************************************************************************
2  *
3  * ViSP, open source Visual Servoing Platform software.
4  * Copyright (C) 2005 - 2019 by Inria. All rights reserved.
5  *
6  * This software is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  * See the file LICENSE.txt at the root directory of this source
11  * distribution for additional information about the GNU GPL.
12  *
13  * For using ViSP with software that can not be combined with the GNU
14  * GPL, please contact Inria about acquiring a ViSP Professional
15  * Edition License.
16  *
17  * See http://visp.inria.fr for more information.
18  *
19  * This software was developed at:
20  * Inria Rennes - Bretagne Atlantique
21  * Campus Universitaire de Beaulieu
22  * 35042 Rennes Cedex
23  * France
24  *
25  * If you have questions regarding the use of this file, please contact
26  * Inria at visp@inria.fr
27  *
28  * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
29  * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
30  *
31  * Description:
32  * Image tools.
33  *
34  * Authors:
35  * Fabien Spindler
36  *
37  *****************************************************************************/
38 
39 #include <visp3/core/vpCPUFeatures.h>
40 #include <visp3/core/vpImageConvert.h>
41 #include <visp3/core/vpImageTools.h>
42 
43 #include <Simd/SimdLib.hpp>
44 
45 #if defined __SSE2__ || defined _M_X64 || (defined _M_IX86_FP && _M_IX86_FP >= 2)
46 #include <emmintrin.h>
47 #define VISP_HAVE_SSE2 1
48 
49 #if defined __SSE3__ || (defined _MSC_VER && _MSC_VER >= 1500)
50 #include <pmmintrin.h>
51 #define VISP_HAVE_SSE3 1
52 #endif
53 #if defined __SSSE3__ || (defined _MSC_VER && _MSC_VER >= 1500)
54 #include <tmmintrin.h>
55 #define VISP_HAVE_SSSE3 1
56 #endif
57 #endif
58 
115 void vpImageTools::changeLUT(vpImage<unsigned char> &I, unsigned char A, unsigned char A_star, unsigned char B,
116  unsigned char B_star)
117 {
118  // Test if input values are valid
119  if (B <= A) {
120  vpERROR_TRACE("Bad gray levels");
122  }
123  unsigned char v;
124 
125  double factor = (double)(B_star - A_star) / (double)(B - A);
126 
127  for (unsigned int i = 0; i < I.getHeight(); i++)
128  for (unsigned int j = 0; j < I.getWidth(); j++) {
129  v = I[i][j];
130 
131  if (v <= A)
132  I[i][j] = A_star;
133  else if (v >= B)
134  I[i][j] = B_star;
135  else
136  I[i][j] = (unsigned char)(A_star + factor * (v - A));
137  }
138 }
139 
153  vpImage<unsigned char> &Idiff)
154 {
155  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
156  throw(vpException(vpException::dimensionError, "The two images have not the same size"));
157  }
158 
159  if ((I1.getHeight() != Idiff.getHeight()) || (I1.getWidth() != Idiff.getWidth())) {
160  Idiff.resize(I1.getHeight(), I1.getWidth());
161  }
162 
163  bool checkSSSE3 = vpCPUFeatures::checkSSSE3();
164 #if !VISP_HAVE_SSSE3
165  checkSSSE3 = false;
166 #endif
167 
168  unsigned int i = 0;
169  if (checkSSSE3) {
170 #if VISP_HAVE_SSSE3
171  if (I1.getSize() >= 16) {
172  const __m128i mask1 = _mm_set_epi8(-1, 14, -1, 12, -1, 10, -1, 8, -1, 6, -1, 4, -1, 2, -1, 0);
173  const __m128i mask2 = _mm_set_epi8(-1, 15, -1, 13, -1, 11, -1, 9, -1, 7, -1, 5, -1, 3, -1, 1);
174 
175  const __m128i mask_out2 = _mm_set_epi8(14, -1, 12, -1, 10, -1, 8, -1, 6, -1, 4, -1, 2, -1, 0, -1);
176 
177  for (; i <= I1.getSize()-16; i+= 16) {
178  const __m128i vdata1 = _mm_loadu_si128(reinterpret_cast<const __m128i *>(I1.bitmap + i));
179  const __m128i vdata2 = _mm_loadu_si128(reinterpret_cast<const __m128i *>(I2.bitmap + i));
180 
181  __m128i vdata1_reorg = _mm_shuffle_epi8(vdata1, mask1);
182  __m128i vdata2_reorg = _mm_shuffle_epi8(vdata2, mask1);
183 
184  const __m128i vshift = _mm_set1_epi16(128);
185  __m128i vdata_diff = _mm_add_epi16(_mm_sub_epi16(vdata1_reorg, vdata2_reorg), vshift);
186 
187  const __m128i v255 = _mm_set1_epi16(255);
188  const __m128i vzero = _mm_setzero_si128();
189  const __m128i vdata_diff_min_max1 = _mm_max_epi16(_mm_min_epi16(vdata_diff, v255), vzero);
190 
191  vdata1_reorg = _mm_shuffle_epi8(vdata1, mask2);
192  vdata2_reorg = _mm_shuffle_epi8(vdata2, mask2);
193 
194  vdata_diff = _mm_add_epi16(_mm_sub_epi16(vdata1_reorg, vdata2_reorg), vshift);
195  const __m128i vdata_diff_min_max2 = _mm_max_epi16(_mm_min_epi16(vdata_diff, v255), vzero);
196 
197  _mm_storeu_si128(reinterpret_cast<__m128i *>(Idiff.bitmap + i), _mm_or_si128(_mm_shuffle_epi8(vdata_diff_min_max1, mask1),
198  _mm_shuffle_epi8(vdata_diff_min_max2, mask_out2)));
199  }
200  }
201 #endif
202  }
203 
204  for (; i < I1.getSize(); i++) {
205  int diff = I1.bitmap[i] - I2.bitmap[i] + 128;
206  Idiff.bitmap[i] = static_cast<unsigned char>(vpMath::maximum(vpMath::minimum(diff, 255), 0));
207  }
208 }
209 
224 {
225  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
226  throw(vpException(vpException::dimensionError, "Cannot compute image difference. The two images "
227  "(%ux%u) and (%ux%u) have not the same size",
228  I1.getWidth(), I1.getHeight(), I2.getWidth(), I2.getHeight()));
229  }
230 
231  if ((I1.getHeight() != Idiff.getHeight()) || (I1.getWidth() != Idiff.getWidth())) {
232  Idiff.resize(I1.getHeight(), I1.getWidth());
233  }
234 
235  bool checkSSSE3 = vpCPUFeatures::checkSSSE3();
236 #if !VISP_HAVE_SSSE3
237  checkSSSE3 = false;
238 #endif
239 
240  unsigned int i = 0;
241  if (checkSSSE3) {
242 #if VISP_HAVE_SSSE3
243  if (I1.getSize() >= 4) {
244  const __m128i mask1 = _mm_set_epi8(-1, 14, -1, 12, -1, 10, -1, 8, -1, 6, -1, 4, -1, 2, -1, 0);
245  const __m128i mask2 = _mm_set_epi8(-1, 15, -1, 13, -1, 11, -1, 9, -1, 7, -1, 5, -1, 3, -1, 1);
246 
247  const __m128i mask_out2 = _mm_set_epi8(14, -1, 12, -1, 10, -1, 8, -1, 6, -1, 4, -1, 2, -1, 0, -1);
248 
249  for (; i <= I1.getSize()-4; i+= 4) {
250  const __m128i vdata1 = _mm_loadu_si128(reinterpret_cast<const __m128i *>(I1.bitmap + i));
251  const __m128i vdata2 = _mm_loadu_si128(reinterpret_cast<const __m128i *>(I2.bitmap + i));
252 
253  __m128i vdata1_reorg = _mm_shuffle_epi8(vdata1, mask1);
254  __m128i vdata2_reorg = _mm_shuffle_epi8(vdata2, mask1);
255 
256  const __m128i vshift = _mm_set1_epi16(128);
257  __m128i vdata_diff = _mm_add_epi16(_mm_sub_epi16(vdata1_reorg, vdata2_reorg), vshift);
258 
259  const __m128i v255 = _mm_set1_epi16(255);
260  const __m128i vzero = _mm_setzero_si128();
261  const __m128i vdata_diff_min_max1 = _mm_max_epi16(_mm_min_epi16(vdata_diff, v255), vzero);
262 
263  vdata1_reorg = _mm_shuffle_epi8(vdata1, mask2);
264  vdata2_reorg = _mm_shuffle_epi8(vdata2, mask2);
265 
266  vdata_diff = _mm_add_epi16(_mm_sub_epi16(vdata1_reorg, vdata2_reorg), vshift);
267  const __m128i vdata_diff_min_max2 = _mm_max_epi16(_mm_min_epi16(vdata_diff, v255), vzero);
268 
269  _mm_storeu_si128(reinterpret_cast<__m128i *>(Idiff.bitmap + i), _mm_or_si128(_mm_shuffle_epi8(vdata_diff_min_max1, mask1),
270  _mm_shuffle_epi8(vdata_diff_min_max2, mask_out2)));
271  }
272  }
273 #endif
274  }
275 
276  for (; i < I1.getSize(); i++) {
277  int diffR = I1.bitmap[i].R - I2.bitmap[i].R + 128;
278  int diffG = I1.bitmap[i].G - I2.bitmap[i].G + 128;
279  int diffB = I1.bitmap[i].B - I2.bitmap[i].B + 128;
280  int diffA = I1.bitmap[i].A - I2.bitmap[i].A + 128;
281  Idiff.bitmap[i].R = static_cast<unsigned char>(vpMath::maximum(vpMath::minimum(diffR, 255), 0));
282  Idiff.bitmap[i].G = static_cast<unsigned char>(vpMath::maximum(vpMath::minimum(diffG, 255), 0));
283  Idiff.bitmap[i].B = static_cast<unsigned char>(vpMath::maximum(vpMath::minimum(diffB, 255), 0));
284  Idiff.bitmap[i].A = static_cast<unsigned char>(vpMath::maximum(vpMath::minimum(diffA, 255), 0));
285  }
286 }
287 
299  vpImage<unsigned char> &Idiff)
300 {
301  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
302  throw(vpException(vpException::dimensionError, "The two images do not have the same size"));
303  }
304 
305  if ((I1.getHeight() != Idiff.getHeight()) || (I1.getWidth() != Idiff.getWidth())) {
306  Idiff.resize(I1.getHeight(), I1.getWidth());
307  }
308 
309  unsigned int n = I1.getHeight() * I1.getWidth();
310  for (unsigned int b = 0; b < n; b++) {
311  int diff = I1.bitmap[b] - I2.bitmap[b];
312  Idiff.bitmap[b] = static_cast<unsigned char>(vpMath::abs(diff));
313  }
314 }
315 
324 {
325  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
326  throw(vpException(vpException::dimensionError, "The two images do not have the same size"));
327  }
328 
329  if ((I1.getHeight() != Idiff.getHeight()) || (I1.getWidth() != Idiff.getWidth())) {
330  Idiff.resize(I1.getHeight(), I1.getWidth());
331  }
332 
333  unsigned int n = I1.getHeight() * I1.getWidth();
334  for (unsigned int b = 0; b < n; b++) {
335  Idiff.bitmap[b] = vpMath::abs(I1.bitmap[b] - I2.bitmap[b]);
336  }
337 }
338 
353 {
354  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
355  throw(vpException(vpException::dimensionError, "The two images do not have the same size"));
356  }
357 
358  if ((I1.getHeight() != Idiff.getHeight()) || (I1.getWidth() != Idiff.getWidth())) {
359  Idiff.resize(I1.getHeight(), I1.getWidth());
360  }
361 
362  unsigned int n = I1.getHeight() * I1.getWidth();
363  for (unsigned int b = 0; b < n; b++) {
364  int diffR = I1.bitmap[b].R - I2.bitmap[b].R;
365  int diffG = I1.bitmap[b].G - I2.bitmap[b].G;
366  int diffB = I1.bitmap[b].B - I2.bitmap[b].B;
367  // int diffA = I1.bitmap[b].A - I2.bitmap[b].A;
368  Idiff.bitmap[b].R = static_cast<unsigned char>(vpMath::abs(diffR));
369  Idiff.bitmap[b].G = static_cast<unsigned char>(vpMath::abs(diffG));
370  Idiff.bitmap[b].B = static_cast<unsigned char>(vpMath::abs(diffB));
371  // Idiff.bitmap[b].A = diffA;
372  Idiff.bitmap[b].A = 0;
373  }
374 }
375 
390  vpImage<unsigned char> &Ires, bool saturate)
391 {
392  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
393  throw(vpException(vpException::dimensionError, "The two images do not have the same size"));
394  }
395 
396  if ((I1.getHeight() != Ires.getHeight()) || (I1.getWidth() != Ires.getWidth())) {
397  Ires.resize(I1.getHeight(), I1.getWidth());
398  }
399 
400  typedef Simd::View<Simd::Allocator> View;
401  View img1(I1.getWidth(), I1.getHeight(), I1.getWidth(), View::Gray8, I1.bitmap);
402  View img2(I2.getWidth(), I2.getHeight(), I2.getWidth(), View::Gray8, I2.bitmap);
403  View imgAdd(Ires.getWidth(), Ires.getHeight(), Ires.getWidth(), View::Gray8, Ires.bitmap);
404 
405  Simd::OperationBinary8u(img1, img2, imgAdd, saturate ? SimdOperationBinary8uSaturatedAddition : SimdOperationBinary8uAddition);
406 }
407 
422  vpImage<unsigned char> &Ires, bool saturate)
423 {
424  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
425  throw(vpException(vpException::dimensionError, "The two images do not have the same size"));
426  }
427 
428  if ((I1.getHeight() != Ires.getHeight()) || (I1.getWidth() != Ires.getWidth())) {
429  Ires.resize(I1.getHeight(), I1.getWidth());
430  }
431 
432  typedef Simd::View<Simd::Allocator> View;
433  View img1(I1.getWidth(), I1.getHeight(), I1.getWidth(), View::Gray8, I1.bitmap);
434  View img2(I2.getWidth(), I2.getHeight(), I2.getWidth(), View::Gray8, I2.bitmap);
435  View imgAdd(Ires.getWidth(), Ires.getHeight(), Ires.getWidth(), View::Gray8, Ires.bitmap);
436 
437  Simd::OperationBinary8u(img1, img2, imgAdd, saturate ? SimdOperationBinary8uSaturatedSubtraction : SimdOperationBinary8uSubtraction);
438 }
439 
451 void vpImageTools::initUndistortMap(const vpCameraParameters &cam, unsigned int width, unsigned int height,
452  vpArray2D<int> &mapU, vpArray2D<int> &mapV,
453  vpArray2D<float> &mapDu, vpArray2D<float> &mapDv)
454 {
455  mapU.resize(height, width, false, false);
456  mapV.resize(height, width, false, false);
457  mapDu.resize(height, width, false, false);
458  mapDv.resize(height, width, false, false);
459 
461  bool is_KannalaBrandt = (projModel==vpCameraParameters::ProjWithKannalaBrandtDistortion); // Check the projection model used
462 
463  float u0 = static_cast<float>(cam.get_u0());
464  float v0 = static_cast<float>(cam.get_v0());
465  float px = static_cast<float>(cam.get_px());
466  float py = static_cast<float>(cam.get_py());
467  float kud;
468  std::vector<double> dist_coefs;
469 
470  if(!is_KannalaBrandt)
471  kud = static_cast<float>(cam.get_kud());
472 
473  else
474  dist_coefs = cam.getKannalaBrandtDistortionCoefficients();
475 
476  if (!is_KannalaBrandt && std::fabs(static_cast<double>(kud)) <= std::numeric_limits<double>::epsilon()) {
477  // There is no need to undistort the image (Perpective projection)
478  for (unsigned int i = 0; i < height; i++) {
479  for (unsigned int j = 0; j < width; j++) {
480  mapU[i][j] = static_cast<int>(j);
481  mapV[i][j] = static_cast<int>(i);
482  mapDu[i][j] = 0;
483  mapDv[i][j] = 0;
484  }
485  }
486 
487  return;
488  }
489 
490  float invpx, invpy;
491  float kud_px2 = 0., kud_py2 = 0., deltau_px, deltav_py;
492  float fr1, fr2;
493  float deltav, deltau;
494  float u_float, v_float;
495  int u_round, v_round;
496  double r, scale;
497  double theta, theta_d;
498  double theta2, theta4, theta6, theta8;
499 
500  invpx = 1.0f / px;
501  invpy = 1.0f / py;
502 
503  if(!is_KannalaBrandt)
504  {
505  kud_px2 = kud * invpx * invpx;
506  kud_py2 = kud * invpy * invpy;
507  }
508 
509  for (unsigned int v = 0; v < height; v++) {
510  deltav = v - v0;
511 
512  if(!is_KannalaBrandt)
513  fr1 = 1.0f + kud_py2 * deltav * deltav;
514  else
515  deltav_py = deltav * invpy;
516 
517  for (unsigned int u = 0; u < width; u++) {
518  // computation of u,v : corresponding pixel coordinates in I.
519  deltau = u - u0;
520  if(!is_KannalaBrandt)
521  {
522  fr2 = fr1 + kud_px2 * deltau * deltau;
523 
524  u_float = deltau * fr2 + u0;
525  v_float = deltav * fr2 + v0;
526  }
527 
528  else
529  {
530  deltau_px = deltau * invpx;
531  r = sqrt(vpMath::sqr(deltau_px) + vpMath::sqr(deltav_py));
532  theta = atan(r);
533 
534  theta2 = vpMath::sqr(theta);
535  theta4 = vpMath::sqr(theta2);
536  theta6 = theta2 * theta4;
537  theta8 = vpMath::sqr(theta4);
538 
539  theta_d = theta * (1 + dist_coefs[0]*theta2 + dist_coefs[1]*theta4 +
540  dist_coefs[2]*theta6 + dist_coefs[3]*theta8);
541 
542  //scale = (r == 0) ? 1.0 : theta_d / r;
543  scale = (std::fabs(r) < std::numeric_limits<double>::epsilon()) ? 1.0 : theta_d / r;
544  u_float = static_cast<float>(deltau*scale + u0);
545  v_float = static_cast<float>(deltav*scale + v0);
546  }
547 
548  u_round = static_cast<int>(u_float);
549  v_round = static_cast<int>(v_float);
550 
551  mapU[v][u] = u_round;
552  mapV[v][u] = v_round;
553 
554  mapDu[v][u] = u_float - u_round;
555  mapDv[v][u] = v_float - v_round;
556  }
557  }
558 }
559 
572 {
573  if (I.getSize() == 0) {
574  std::cerr << "Error, input image is empty." << std::endl;
575  return;
576  }
577 
578  II.resize(I.getHeight() + 1, I.getWidth() + 1, 0.0);
579  IIsq.resize(I.getHeight() + 1, I.getWidth() + 1, 0.0);
580 
581  for (unsigned int i = 1; i < II.getHeight(); i++) {
582  for (unsigned int j = 1; j < II.getWidth(); j++) {
583  II[i][j] = I[i - 1][j - 1] + II[i - 1][j] + II[i][j - 1] - II[i - 1][j - 1];
584  IIsq[i][j] = vpMath::sqr(I[i - 1][j - 1]) + IIsq[i - 1][j] + IIsq[i][j - 1] - IIsq[i - 1][j - 1];
585  }
586  }
587 }
588 
597 #if VISP_HAVE_SSE2
598  bool useOptimized)
599 #else
600  const bool)
601 #endif
602 {
603  if ((I1.getHeight() != I2.getHeight()) || (I1.getWidth() != I2.getWidth())) {
604  throw vpException(vpException::dimensionError, "Error: in vpImageTools::normalizedCorrelation(): "
605  "image dimension mismatch between I1=%ux%u and I2=%ux%u",
606  I1.getHeight(), I1.getWidth(), I2.getHeight(), I2.getWidth());
607  }
608 
609  const double a = I1.getMeanValue();
610  const double b = I2.getMeanValue();
611 
612  double ab = 0.0;
613  double a2 = 0.0;
614  double b2 = 0.0;
615 
616  unsigned int cpt = 0;
617 
618 #if VISP_HAVE_SSE2
619  if (vpCPUFeatures::checkSSE2() && I1.getSize() >= 2 && useOptimized) {
620  const double *ptr_I1 = I1.bitmap;
621  const double *ptr_I2 = I2.bitmap;
622 
623  const __m128d v_mean_a = _mm_set1_pd(a);
624  const __m128d v_mean_b = _mm_set1_pd(b);
625  __m128d v_ab = _mm_setzero_pd();
626  __m128d v_a2 = _mm_setzero_pd();
627  __m128d v_b2 = _mm_setzero_pd();
628 
629  for (; cpt <= I1.getSize() - 2; cpt += 2, ptr_I1 += 2, ptr_I2 += 2) {
630  const __m128d v1 = _mm_loadu_pd(ptr_I1);
631  const __m128d v2 = _mm_loadu_pd(ptr_I2);
632  const __m128d norm_a = _mm_sub_pd(v1, v_mean_a);
633  const __m128d norm_b = _mm_sub_pd(v2, v_mean_b);
634  v_ab = _mm_add_pd(v_ab, _mm_mul_pd(norm_a, norm_b));
635  v_a2 = _mm_add_pd(v_a2, _mm_mul_pd(norm_a, norm_a));
636  v_b2 = _mm_add_pd(v_b2, _mm_mul_pd(norm_b, norm_b));
637  }
638 
639  double v_res_ab[2], v_res_a2[2], v_res_b2[2];
640  _mm_storeu_pd(v_res_ab, v_ab);
641  _mm_storeu_pd(v_res_a2, v_a2);
642  _mm_storeu_pd(v_res_b2, v_b2);
643 
644  ab = v_res_ab[0] + v_res_ab[1];
645  a2 = v_res_a2[0] + v_res_a2[1];
646  b2 = v_res_b2[0] + v_res_b2[1];
647  }
648 #endif
649 
650  for (; cpt < I1.getSize(); cpt++) {
651  ab += (I1.bitmap[cpt] - a) * (I2.bitmap[cpt] - b);
652  a2 += vpMath::sqr(I1.bitmap[cpt] - a);
653  b2 += vpMath::sqr(I2.bitmap[cpt] - b);
654  }
655 
656  return ab / sqrt(a2 * b2);
657 }
658 
666 {
667  unsigned int height = I.getHeight(), width = I.getWidth();
668  V.resize(width); // resize and nullify
669 
670  for (unsigned int i = 0; i < height; ++i)
671  for (unsigned int j = 0; j < width; ++j)
672  V[j] += I[i][j];
673  for (unsigned int j = 0; j < width; ++j)
674  V[j] /= height;
675 }
676 
682 {
683  double s = I.getSum();
684  for (unsigned int i = 0; i < I.getHeight(); ++i)
685  for (unsigned int j = 0; j < I.getWidth(); ++j)
686  I(i, j, I(i, j) / s);
687 }
688 
697  const vpImageInterpolationType &method)
698 {
699  switch (method) {
701  return I(vpMath::round(point.get_i()), vpMath::round(point.get_j()));
702  case INTERPOLATION_LINEAR: {
703  int x1 = (int)floor(point.get_i());
704  int x2 = (int)ceil(point.get_i());
705  int y1 = (int)floor(point.get_j());
706  int y2 = (int)ceil(point.get_j());
707  double v1, v2;
708  if (x1 == x2) {
709  v1 = I(x1, y1);
710  v2 = I(x1, y2);
711  } else {
712  v1 = (x2 - point.get_i()) * I(x1, y1) + (point.get_i() - x1) * I(x2, y1);
713  v2 = (x2 - point.get_i()) * I(x1, y2) + (point.get_i() - x1) * I(x2, y2);
714  }
715  if (y1 == y2)
716  return v1;
717  return (y2 - point.get_j()) * v1 + (point.get_j() - y1) * v2;
718  }
719  case INTERPOLATION_CUBIC: {
721  "vpImageTools::interpolate(): bi-cubic interpolation is not implemented.");
722  }
723  default: {
724  throw vpException(vpException::notImplementedError, "vpImageTools::interpolate(): invalid interpolation type");
725  }
726  }
727 }
728 
736 {
737  unsigned int x_d = vpMath::round(r.getHeight());
738  unsigned int y_d = vpMath::round(r.getWidth());
739  double x1 = r.getTopLeft().get_i();
740  double y1 = r.getTopLeft().get_j();
741  double t = r.getOrientation();
742  Dst.resize(x_d, y_d);
743  for (unsigned int x = 0; x < x_d; ++x) {
744  for (unsigned int y = 0; y < y_d; ++y) {
745  Dst(x, y,
746  (unsigned char)interpolate(Src, vpImagePoint(x1 + x * cos(t) + y * sin(t), y1 - x * sin(t) + y * cos(t)),
748  }
749  }
750 }
751 
759 {
760  unsigned int x_d = vpMath::round(r.getHeight());
761  unsigned int y_d = vpMath::round(r.getWidth());
762  double x1 = r.getTopLeft().get_i();
763  double y1 = r.getTopLeft().get_j();
764  double t = r.getOrientation();
765  Dst.resize(x_d, y_d);
766  for (unsigned int x = 0; x < x_d; ++x) {
767  for (unsigned int y = 0; y < y_d; ++y) {
768  Dst(x, y, interpolate(Src, vpImagePoint(x1 + x * cos(t) + y * sin(t), y1 - x * sin(t) + y * cos(t)),
770  }
771  }
772 }
773 
789  vpImage<double> &I_score, unsigned int step_u, unsigned int step_v,
790  bool useOptimized)
791 {
792  if (I.getSize() == 0) {
793  std::cerr << "Error, input image is empty." << std::endl;
794  return;
795  }
796 
797  if (I_tpl.getSize() == 0) {
798  std::cerr << "Error, template image is empty." << std::endl;
799  return;
800  }
801 
802  if (I_tpl.getHeight() > I.getHeight() || I_tpl.getWidth() > I.getWidth()) {
803  std::cerr << "Error, template image is bigger than input image." << std::endl;
804  return;
805  }
806 
807  vpImage<double> I_double, I_tpl_double;
808  vpImageConvert::convert(I, I_double);
809  vpImageConvert::convert(I_tpl, I_tpl_double);
810 
811  unsigned int height_tpl = I_tpl.getHeight(), width_tpl = I_tpl.getWidth();
812  I_score.resize(I.getHeight() - height_tpl, I.getWidth() - width_tpl, 0.0);
813 
814  if (useOptimized) {
815  vpImage<double> II, IIsq;
816  integralImage(I, II, IIsq);
817 
818  vpImage<double> II_tpl, IIsq_tpl;
819  integralImage(I_tpl, II_tpl, IIsq_tpl);
820 
821  // zero-mean template image
822  const double sum2 = (II_tpl[height_tpl][width_tpl] + II_tpl[0][0] - II_tpl[0][width_tpl] - II_tpl[height_tpl][0]);
823  const double mean2 = sum2 / I_tpl.getSize();
824  for (unsigned int cpt = 0; cpt < I_tpl_double.getSize(); cpt++) {
825  I_tpl_double.bitmap[cpt] -= mean2;
826  }
827 
828 #if defined _OPENMP && _OPENMP >= 200711 // OpenMP 3.1
829 #pragma omp parallel for schedule(dynamic)
830  for (unsigned int i = 0; i < I.getHeight() - height_tpl; i += step_v) {
831  for (unsigned int j = 0; j < I.getWidth() - width_tpl; j += step_u) {
832  I_score[i][j] = normalizedCorrelation(I_double, I_tpl_double, II, IIsq, II_tpl, IIsq_tpl, i, j);
833  }
834  }
835 #else
836  // error C3016: 'i': index variable in OpenMP 'for' statement must have signed integral type
837  int end = (int)((I.getHeight() - height_tpl) / step_v) + 1;
838  std::vector<unsigned int> vec_step_v((size_t)end);
839  for (unsigned int cpt = 0, idx = 0; cpt < I.getHeight() - height_tpl; cpt += step_v, idx++) {
840  vec_step_v[(size_t)idx] = cpt;
841  }
842 #if defined _OPENMP // only to disable warning: ignoring #pragma omp parallel [-Wunknown-pragmas]
843 #pragma omp parallel for schedule(dynamic)
844 #endif
845  for (int cpt = 0; cpt < end; cpt++) {
846  for (unsigned int j = 0; j < I.getWidth() - width_tpl; j += step_u) {
847  I_score[vec_step_v[cpt]][j] =
848  normalizedCorrelation(I_double, I_tpl_double, II, IIsq, II_tpl, IIsq_tpl, vec_step_v[cpt], j);
849  }
850  }
851 #endif
852  } else {
853  vpImage<double> I_cur;
854 
855  for (unsigned int i = 0; i < I.getHeight() - height_tpl; i += step_v) {
856  for (unsigned int j = 0; j < I.getWidth() - width_tpl; j += step_u) {
857  vpRect roi(vpImagePoint(i, j), vpImagePoint(i + height_tpl - 1, j + width_tpl - 1));
858  vpImageTools::crop(I_double, roi, I_cur);
859 
860  I_score[i][j] = vpImageTools::normalizedCorrelation(I_cur, I_tpl_double, useOptimized);
861  }
862  }
863  }
864 }
865 
866 // Reference:
867 // http://blog.demofox.org/2015/08/15/resizing-images-with-bicubic-interpolation/
868 // t is a value that goes from 0 to 1 to interpolate in a C1 continuous way
869 // across uniformly sampled data points. when t is 0, this will return B.
870 // When t is 1, this will return C. In between values will return an
871 // interpolation between B and C. A and B are used to calculate the slopes at
872 // the edges.
873 float vpImageTools::cubicHermite(const float A, const float B, const float C, const float D, const float t)
874 {
875  float a = (-A + 3.0f * B - 3.0f * C + D) / 2.0f;
876  float b = A + 2.0f * C - (5.0f * B + D) / 2.0f;
877  float c = (-A + C) / 2.0f;
878  float d = B;
879 
880  return a * t * t * t + b * t * t + c * t + d;
881 }
882 
883 int vpImageTools::coordCast(double x)
884 {
885  return x < 0 ? -1 : static_cast<int>(x);
886 }
887 
888 double vpImageTools::lerp(double A, double B, double t) {
889  return A * (1.0 - t) + B * t;
890 }
891 
892 float vpImageTools::lerp(float A, float B, float t) {
893  return A * (1.0f - t) + B * t;
894 }
895 
896 int64_t vpImageTools::lerp2(int64_t A, int64_t B, int64_t t, int64_t t_1) {
897  return A * t_1 + B * t;
898 }
899 
901  const vpImage<double> &II, const vpImage<double> &IIsq,
902  const vpImage<double> &II_tpl, const vpImage<double> &IIsq_tpl,
903  unsigned int i0, unsigned int j0)
904 {
905  double ab = 0.0;
906 #if VISP_HAVE_SSE2
907  bool use_sse_version = true;
908  if (vpCPUFeatures::checkSSE2() && I2.getWidth() >= 2) {
909  const double *ptr_I1 = I1.bitmap;
910  const double *ptr_I2 = I2.bitmap;
911 
912  __m128d v_ab = _mm_setzero_pd();
913 
914  for (unsigned int i = 0; i < I2.getHeight(); i++) {
915  unsigned int j = 0;
916  ptr_I1 = &I1.bitmap[(i0 + i) * I1.getWidth() + j0];
917 
918  for (; j <= I2.getWidth() - 2; j += 2, ptr_I1 += 2, ptr_I2 += 2) {
919  const __m128d v1 = _mm_loadu_pd(ptr_I1);
920  const __m128d v2 = _mm_loadu_pd(ptr_I2);
921  v_ab = _mm_add_pd(v_ab, _mm_mul_pd(v1, v2));
922  }
923 
924  for (; j < I2.getWidth(); j++) {
925  ab += (I1[i0 + i][j0 + j]) * I2[i][j];
926  }
927  }
928 
929  double v_res_ab[2];
930  _mm_storeu_pd(v_res_ab, v_ab);
931 
932  ab += v_res_ab[0] + v_res_ab[1];
933  } else {
934  use_sse_version = false;
935  }
936 #else
937  bool use_sse_version = false;
938 #endif
939 
940  if (!use_sse_version) {
941  for (unsigned int i = 0; i < I2.getHeight(); i++) {
942  for (unsigned int j = 0; j < I2.getWidth(); j++) {
943  ab += (I1[i0 + i][j0 + j]) * I2[i][j];
944  }
945  }
946  }
947 
948  unsigned int height_tpl = I2.getHeight(), width_tpl = I2.getWidth();
949  const double sum1 =
950  (II[i0 + height_tpl][j0 + width_tpl] + II[i0][j0] - II[i0][j0 + width_tpl] - II[i0 + height_tpl][j0]);
951  const double sum2 = (II_tpl[height_tpl][width_tpl] + II_tpl[0][0] - II_tpl[0][width_tpl] - II_tpl[height_tpl][0]);
952 
953  double a2 = ((IIsq[i0 + I2.getHeight()][j0 + I2.getWidth()] + IIsq[i0][j0] - IIsq[i0][j0 + I2.getWidth()] -
954  IIsq[i0 + I2.getHeight()][j0]) -
955  (1.0 / I2.getSize()) * vpMath::sqr(sum1));
956 
957  double b2 = ((IIsq_tpl[I2.getHeight()][I2.getWidth()] + IIsq_tpl[0][0] - IIsq_tpl[0][I2.getWidth()] -
958  IIsq_tpl[I2.getHeight()][0]) -
959  (1.0 / I2.getSize()) * vpMath::sqr(sum2));
960  return ab / sqrt(a2 * b2);
961 }
962 
974  const vpArray2D<float> &mapDu, const vpArray2D<float> &mapDv, vpImage<unsigned char> &Iundist)
975 {
976  Iundist.resize(I.getHeight(), I.getWidth());
977 
978 #if defined _OPENMP // only to disable warning: ignoring #pragma omp parallel [-Wunknown-pragmas]
979 #pragma omp parallel for schedule(dynamic)
980 #endif
981  for (int i_ = 0; i_ < static_cast<int>(I.getHeight()); i_++) {
982  const unsigned int i = static_cast<unsigned int>(i_);
983  for (unsigned int j = 0; j < I.getWidth(); j++) {
984 
985  int u_round = mapU[i][j];
986  int v_round = mapV[i][j];
987 
988  float du = mapDu[i][j];
989  float dv = mapDv[i][j];
990 
991  if (0 <= u_round && 0 <= v_round && u_round < static_cast<int>(I.getWidth()) - 1
992  && v_round < static_cast<int>(I.getHeight()) - 1) {
993  // process interpolation
994  float col0 = lerp(I[v_round][u_round], I[v_round][u_round + 1], du);
995  float col1 = lerp(I[v_round + 1][u_round], I[v_round + 1][u_round + 1], du);
996  float value = lerp(col0, col1, dv);
997 
998  Iundist[i][j] = static_cast<unsigned char>(value);
999  } else {
1000  Iundist[i][j] = 0;
1001  }
1002  }
1003  }
1004 }
1005 
1016 void vpImageTools::remap(const vpImage<vpRGBa> &I, const vpArray2D<int> &mapU, const vpArray2D<int> &mapV,
1017  const vpArray2D<float> &mapDu, const vpArray2D<float> &mapDv, vpImage<vpRGBa> &Iundist)
1018 {
1019  Iundist.resize(I.getHeight(), I.getWidth());
1020 
1021  bool checkSSE2 = vpCPUFeatures::checkSSE2();
1022 #if !VISP_HAVE_SSE2
1023  checkSSE2 = false;
1024 #endif
1025 
1026  if (checkSSE2) {
1027 #if defined VISP_HAVE_SSE2
1028 #if defined _OPENMP // only to disable warning: ignoring #pragma omp parallel [-Wunknown-pragmas]
1029 #pragma omp parallel for schedule(dynamic)
1030 #endif
1031  for (int i_ = 0; i_ < static_cast<int>(I.getHeight()); i_++) {
1032  const unsigned int i = static_cast<unsigned int>(i_);
1033  for (unsigned int j = 0; j < I.getWidth(); j++) {
1034 
1035  int u_round = mapU[i][j];
1036  int v_round = mapV[i][j];
1037 
1038  const __m128 vdu = _mm_set1_ps(mapDu[i][j]);
1039  const __m128 vdv = _mm_set1_ps(mapDv[i][j]);
1040 
1041  if (0 <= u_round && 0 <= v_round && u_round < static_cast<int>(I.getWidth()) - 1
1042  && v_round < static_cast<int>(I.getHeight()) - 1) {
1043  #define VLERP(va, vb, vt) _mm_add_ps(va, _mm_mul_ps(_mm_sub_ps(vb, va), vt));
1044 
1045  // process interpolation
1046  const __m128 vdata1 =
1047  _mm_set_ps(static_cast<float>(I[v_round][u_round].A), static_cast<float>(I[v_round][u_round].B),
1048  static_cast<float>(I[v_round][u_round].G), static_cast<float>(I[v_round][u_round].R));
1049 
1050  const __m128 vdata2 =
1051  _mm_set_ps(static_cast<float>(I[v_round][u_round + 1].A), static_cast<float>(I[v_round][u_round + 1].B),
1052  static_cast<float>(I[v_round][u_round + 1].G), static_cast<float>(I[v_round][u_round + 1].R));
1053 
1054  const __m128 vdata3 =
1055  _mm_set_ps(static_cast<float>(I[v_round + 1][u_round].A), static_cast<float>(I[v_round + 1][u_round].B),
1056  static_cast<float>(I[v_round + 1][u_round].G), static_cast<float>(I[v_round + 1][u_round].R));
1057 
1058  const __m128 vdata4 = _mm_set_ps(
1059  static_cast<float>(I[v_round + 1][u_round + 1].A), static_cast<float>(I[v_round + 1][u_round + 1].B),
1060  static_cast<float>(I[v_round + 1][u_round + 1].G), static_cast<float>(I[v_round + 1][u_round + 1].R));
1061 
1062  const __m128 vcol0 = VLERP(vdata1, vdata2, vdu);
1063  const __m128 vcol1 = VLERP(vdata3, vdata4, vdu);
1064  const __m128 vvalue = VLERP(vcol0, vcol1, vdv);
1065 
1066  #undef VLERP
1067 
1068  float values[4];
1069  _mm_storeu_ps(values, vvalue);
1070  Iundist[i][j].R = static_cast<unsigned char>(values[0]);
1071  Iundist[i][j].G = static_cast<unsigned char>(values[1]);
1072  Iundist[i][j].B = static_cast<unsigned char>(values[2]);
1073  Iundist[i][j].A = static_cast<unsigned char>(values[3]);
1074  } else {
1075  Iundist[i][j] = 0;
1076  }
1077  }
1078  }
1079 #endif
1080  } else {
1081 #if defined _OPENMP // only to disable warning: ignoring #pragma omp parallel [-Wunknown-pragmas]
1082 #pragma omp parallel for schedule(dynamic)
1083 #endif
1084  for (int i_ = 0; i_ < static_cast<int>(I.getHeight()); i_++) {
1085  const unsigned int i = static_cast<unsigned int>(i_);
1086  for (unsigned int j = 0; j < I.getWidth(); j++) {
1087 
1088  int u_round = mapU[i][j];
1089  int v_round = mapV[i][j];
1090 
1091  float du = mapDu[i][j];
1092  float dv = mapDv[i][j];
1093 
1094  if (0 <= u_round && 0 <= v_round && u_round < static_cast<int>(I.getWidth()) - 1
1095  && v_round < static_cast<int>(I.getHeight()) - 1) {
1096  // process interpolation
1097  float col0 = lerp(I[v_round][u_round].R, I[v_round][u_round + 1].R, du);
1098  float col1 = lerp(I[v_round + 1][u_round].R, I[v_round + 1][u_round + 1].R, du);
1099  float value = lerp(col0, col1, dv);
1100 
1101  Iundist[i][j].R = static_cast<unsigned char>(value);
1102 
1103  col0 = lerp(I[v_round][u_round].G, I[v_round][u_round + 1].G, du);
1104  col1 = lerp(I[v_round + 1][u_round].G, I[v_round + 1][u_round + 1].G, du);
1105  value = lerp(col0, col1, dv);
1106 
1107  Iundist[i][j].G = static_cast<unsigned char>(value);
1108 
1109  col0 = lerp(I[v_round][u_round].B, I[v_round][u_round + 1].B, du);
1110  col1 = lerp(I[v_round + 1][u_round].B, I[v_round + 1][u_round + 1].B, du);
1111  value = lerp(col0, col1, dv);
1112 
1113  Iundist[i][j].B = static_cast<unsigned char>(value);
1114 
1115  col0 = lerp(I[v_round][u_round].A, I[v_round][u_round + 1].A, du);
1116  col1 = lerp(I[v_round + 1][u_round].A, I[v_round + 1][u_round + 1].A, du);
1117  value = lerp(col0, col1, dv);
1118 
1119  Iundist[i][j].A = static_cast<unsigned char>(value);
1120  } else {
1121  Iundist[i][j] = 0;
1122  }
1123  }
1124  }
1125  }
1126 }
1127 
1128 void vpImageTools::resizeSimdlib(const vpImage<vpRGBa> &Isrc, unsigned int resizeWidth,
1129  unsigned int resizeHeight, vpImage<vpRGBa> &Idst,
1130  int method)
1131 {
1132  Idst.resize(resizeHeight, resizeWidth);
1133 
1134  typedef Simd::View<Simd::Allocator> View;
1135  View src(Isrc.getWidth(), Isrc.getHeight(), Isrc.getWidth() * sizeof(vpRGBa), View::Bgra32, Isrc.bitmap);
1136  View dst(Idst.getWidth(), Idst.getHeight(), Idst.getWidth() * sizeof(vpRGBa), View::Bgra32, Idst.bitmap);
1137 
1138  Simd::Resize(src, dst, method == INTERPOLATION_LINEAR ? SimdResizeMethodBilinear : SimdResizeMethodArea);
1139 }
1140 
1141 void vpImageTools::resizeSimdlib(const vpImage<unsigned char> &Isrc, unsigned int resizeWidth,
1142  unsigned int resizeHeight, vpImage<unsigned char> &Idst,
1143  int method)
1144 {
1145  Idst.resize(resizeHeight, resizeWidth);
1146 
1147  typedef Simd::View<Simd::Allocator> View;
1148  View src(Isrc.getWidth(), Isrc.getHeight(), Isrc.getWidth(), View::Gray8, Isrc.bitmap);
1149  View dst(Idst.getWidth(), Idst.getHeight(), Idst.getWidth(), View::Gray8, Idst.bitmap);
1150 
1151  Simd::Resize(src, dst, method == INTERPOLATION_LINEAR ? SimdResizeMethodBilinear : SimdResizeMethodArea);
1152 }
1153 
1154 bool vpImageTools::checkFixedPoint(unsigned int x, unsigned int y, const vpMatrix &T, bool affine)
1155 {
1156  double a0 = T[0][0]; double a1 = T[0][1]; double a2 = T[0][2];
1157  double a3 = T[1][0]; double a4 = T[1][1]; double a5 = T[1][2];
1158  double a6 = affine ? 0.0 : T[2][0];
1159  double a7 = affine ? 0.0 : T[2][1];
1160  double a8 = affine ? 1.0 : T[2][2];
1161 
1162  double w = a6 * x + a7 * y + a8;
1163  double x2 = (a0 * x + a1 * y + a2) / w;
1164  double y2 = (a3 * x + a4 * y + a5) / w;
1165 
1166  const double limit = 1 << 15;
1167  return (vpMath::abs(x2) < limit) && (vpMath::abs(y2) < limit);
1168 }
Implementation of a matrix and operations on matrices.
Definition: vpMatrix.h:156
static void extract(const vpImage< unsigned char > &Src, vpImage< unsigned char > &Dst, const vpRectOriented &r)
double get_i() const
Definition: vpImagePoint.h:203
static void imageDifferenceAbsolute(const vpImage< unsigned char > &I1, const vpImage< unsigned char > &I2, vpImage< unsigned char > &Idiff)
void resize(unsigned int nrows, unsigned int ncols, bool flagNullify=true, bool recopy_=true)
Definition: vpArray2D.h:305
static void imageAdd(const vpImage< unsigned char > &I1, const vpImage< unsigned char > &I2, vpImage< unsigned char > &Ires, bool saturate=false)
vpCameraParametersProjType get_projModel() const
static void convert(const vpImage< unsigned char > &src, vpImage< vpRGBa > &dest)
void resize(unsigned int h, unsigned int w)
resize the image : Image initialization
Definition: vpImage.h:800
unsigned char B
Blue component.
Definition: vpRGBa.h:150
Implementation of row vector and the associated operations.
Definition: vpRowVector.h:115
Type * bitmap
points toward the bitmap
Definition: vpImage.h:143
double getOrientation() const
Get the rectangle orientation (rad).
#define vpERROR_TRACE
Definition: vpDebug.h:393
error that can be emited by ViSP classes.
Definition: vpException.h:71
Implementation of a generic 2D array used as base class for matrices and vectors. ...
Definition: vpArray2D.h:131
double getSum() const
Definition: vpImage.h:1578
Error that can be emited by the vpImage class and its derivates.
static void normalize(vpImage< double > &I)
unsigned char G
Green component.
Definition: vpRGBa.h:149
static void imageDifference(const vpImage< unsigned char > &I1, const vpImage< unsigned char > &I2, vpImage< unsigned char > &Idiff)
static double interpolate(const vpImage< unsigned char > &I, const vpImagePoint &point, const vpImageInterpolationType &method=INTERPOLATION_NEAREST)
Definition: vpRGBa.h:66
static Type abs(const Type &x)
Definition: vpMath.h:160
static double normalizedCorrelation(const vpImage< double > &I1, const vpImage< double > &I2, bool useOptimized=true)
static Type maximum(const Type &a, const Type &b)
Definition: vpMath.h:145
double getWidth() const
Get the rectangle width.
static void remap(const vpImage< unsigned char > &I, const vpArray2D< int > &mapU, const vpArray2D< int > &mapV, const vpArray2D< float > &mapDu, const vpArray2D< float > &mapDv, vpImage< unsigned char > &Iundist)
static double sqr(double x)
Definition: vpMath.h:116
VISP_EXPORT bool checkSSE2()
double get_j() const
Definition: vpImagePoint.h:214
unsigned char A
Additionnal component.
Definition: vpRGBa.h:151
Generic class defining intrinsic camera parameters.
std::vector< double > getKannalaBrandtDistortionCoefficients() const
static Type minimum(const Type &a, const Type &b)
Definition: vpMath.h:153
VISP_EXPORT bool checkSSSE3()
static void templateMatching(const vpImage< unsigned char > &I, const vpImage< unsigned char > &I_tpl, vpImage< double > &I_score, unsigned int step_u, unsigned int step_v, bool useOptimized=true)
Type getMeanValue() const
Return the mean value of the bitmap.
Definition: vpImage.h:907
static int round(double x)
Definition: vpMath.h:245
static void initUndistortMap(const vpCameraParameters &cam, unsigned int width, unsigned int height, vpArray2D< int > &mapU, vpArray2D< int > &mapV, vpArray2D< float > &mapDu, vpArray2D< float > &mapDv)
static void columnMean(const vpImage< double > &I, vpRowVector &result)
vpImagePoint getTopLeft() const
Get the top-left corner.
unsigned int getHeight() const
Definition: vpImage.h:188
static void integralImage(const vpImage< unsigned char > &I, vpImage< double > &II, vpImage< double > &IIsq)
double getHeight() const
Get the rectangle height.
static void crop(const vpImage< Type > &I, double roi_top, double roi_left, unsigned int roi_height, unsigned int roi_width, vpImage< Type > &crop, unsigned int v_scale=1, unsigned int h_scale=1)
Definition: vpImageTools.h:305
unsigned int getSize() const
Definition: vpImage.h:227
unsigned char R
Red component.
Definition: vpRGBa.h:148
double get_kud() const
Defines a rectangle in the plane.
Definition: vpRect.h:79
static void imageSubtract(const vpImage< unsigned char > &I1, const vpImage< unsigned char > &I2, vpImage< unsigned char > &Ires, bool saturate=false)
Class that defines a 2D point in an image. This class is useful for image processing and stores only ...
Definition: vpImagePoint.h:87
static void changeLUT(vpImage< unsigned char > &I, unsigned char A, unsigned char newA, unsigned char B, unsigned char newB)
unsigned int getWidth() const
Definition: vpImage.h:246
void resize(unsigned int i, bool flagNullify=true)
Definition: vpRowVector.h:271
Defines an oriented rectangle in the plane.