Visual Servoing Platform  version 3.6.1 under development (2023-12-02)
vpPlane.cpp
1 /****************************************************************************
2  *
3  * ViSP, open source Visual Servoing Platform software.
4  * Copyright (C) 2005 - 2023 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 https://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  * Plane geometrical structure.
33  *
34 *****************************************************************************/
35 
42 #include <visp3/core/vpPlane.h>
43 
44 #include <cmath> // std::fabs
45 #include <limits> // numeric_limits
46 
51 {
52  A = p.A;
53  B = p.B;
54  C = p.C;
55  D = p.D;
56 
57  return *this;
58 }
59 
63 vpPlane::vpPlane() : A(0), B(0), C(0), D(0) {}
64 
75 vpPlane::vpPlane(double a, double b, double c, double d) : A(a), B(b), C(c), D(d) {}
76 
80 vpPlane::vpPlane(const vpPlane &P) : A(0), B(0), C(0), D(0)
81 {
82  setA(P.getA());
83  setB(P.getB());
84  setC(P.getC());
85  setD(P.getD());
86 }
87 
107 vpPlane::vpPlane(const vpPoint &P, const vpColVector &n, vpPlaneFrame frame) : A(0), B(0), C(0), D(0)
108 {
109  // Equation of the plane is given by:
110  A = n[0];
111  B = n[1];
112  C = n[2];
113 
114  if (frame == vpPlane::camera_frame)
115  D = -(A * P.get_X() + B * P.get_Y() + C * P.get_Z());
116  else
117  D = -(A * P.get_oX() + B * P.get_oY() + C * P.get_oZ());
118 }
119 
125 void vpPlane::init(const vpPlane &P)
126 {
127  setA(P.getA());
128  setB(P.getB());
129  setC(P.getC());
130  setD(P.getD());
131 }
132 
144 void vpPlane::init(const vpColVector &P, const vpColVector &n)
145 {
146  // Equation of the plane is given by:
147  A = n[0];
148  B = n[1];
149  C = n[2];
150 
151  D = -(A * P[0] + B * P[1] + C * P[2]);
152 }
153 
165 void vpPlane::init(const vpPoint &P, const vpPoint &Q, const vpPoint &R, vpPlaneFrame frame)
166 {
167  vpColVector a(3);
168  vpColVector b(3);
169  vpColVector n(3);
170  if (frame == vpPlane::camera_frame) {
171  // Calculate vector corresponding to PQ
172  a[0] = P.get_X() - Q.get_X();
173  a[1] = P.get_Y() - Q.get_Y();
174  a[2] = P.get_Z() - Q.get_Z();
175 
176  // Calculate vector corresponding to PR
177  b[0] = P.get_X() - R.get_X();
178  b[1] = P.get_Y() - R.get_Y();
179  b[2] = P.get_Z() - R.get_Z();
180  } else {
181  // Calculate vector corresponding to PQ
182  a[0] = P.get_oX() - Q.get_oX();
183  a[1] = P.get_oY() - Q.get_oY();
184  a[2] = P.get_oZ() - Q.get_oZ();
185 
186  // Calculate vector corresponding to PR
187  b[0] = P.get_oX() - R.get_oX();
188  b[1] = P.get_oY() - R.get_oY();
189  b[2] = P.get_oZ() - R.get_oZ();
190  }
191  // Calculate normal vector to plane PQ x PR
192  n = vpColVector::cross(a, b);
193 
194  // Equation of the plane is given by:
195  A = n[0];
196  B = n[1];
197  C = n[2];
198  if (frame == vpPlane::camera_frame)
199  D = -(A * P.get_X() + B * P.get_Y() + C * P.get_Z());
200  else
201  D = -(A * P.get_oX() + B * P.get_oY() + C * P.get_oZ());
202 
203  double norm = sqrt(A * A + B * B + C * C);
204  A /= norm;
205  B /= norm;
206  C /= norm;
207  D /= norm;
208 }
209 
222 vpPlane::vpPlane(const vpPoint &P, const vpPoint &Q, const vpPoint &R, vpPlaneFrame frame) : A(0), B(0), C(0), D(0)
223 {
224  init(P, Q, R, frame);
225 }
226 
232 double vpPlane::computeZ(double x, double y) const
233 {
234  return -getD() / (getA() * x + getB() * y + getC());
235 }
236 
246 {
247  vpColVector n(3);
248  n[0] = A;
249  n[1] = B;
250  n[2] = C;
251 
252  return n;
253 }
254 
266 {
267  n.resize(3);
268  n[0] = A;
269  n[1] = B;
270  n[2] = C;
271 }
272 
279 void vpPlane::projectionPointOnPlan(const vpPoint &P, vpPoint &Pproj) const
280 {
281  double x0, y0, z0;
282  double rho;
283 
284  x0 = P.get_X() / P.get_W();
285  y0 = P.get_Y() / P.get_W();
286  z0 = P.get_Z() / P.get_W();
287 
288  rho = -(A * x0 + B * y0 + C * z0 + D) / (A * A + B * B + C * C);
289 
290  Pproj.set_X(x0 + A * rho);
291  Pproj.set_Y(y0 + B * rho);
292  Pproj.set_Z(z0 + C * rho);
293  Pproj.set_W(1);
294 }
295 
296 double vpPlane::rayIntersection(const vpPoint &M0, const vpPoint &M1, vpColVector &H) const
297 {
298 
299  double k, scal;
300 
301  // if(M0.get_X()!=0 || M0.get_Y()!=0 || M0.get_Z()!=0)
302  if (std::fabs(M0.get_X()) > std::numeric_limits<double>::epsilon() ||
303  std::fabs(M0.get_Y()) > std::numeric_limits<double>::epsilon() ||
304  std::fabs(M0.get_Z()) > std::numeric_limits<double>::epsilon()) {
305  double R[3];
306  R[0] = M1.get_X() - M0.get_X();
307  R[1] = M1.get_Y() - M0.get_Y();
308  R[2] = M1.get_Z() - M0.get_Z();
309 
310  scal = getA() * R[0] + getB() * R[1] + getC() * R[2];
311  // if (scal != 0)
312  if (std::fabs(scal) > std::numeric_limits<double>::epsilon())
313  k = -(getA() * M0.get_X() + getB() * M0.get_Y() + getC() * M0.get_Z() + getD()) / scal;
314  else
315  k = 0;
316 
317  H[0] = M0.get_X() + k * R[0];
318  H[1] = M0.get_Y() + k * R[1];
319  H[2] = M0.get_Z() + k * R[2];
320  } else {
321  scal = getA() * M1.get_X() + getB() * M1.get_Y() + getC() * M1.get_Z();
322  // if (scal != 0)
323  if (std::fabs(scal) > std::numeric_limits<double>::epsilon())
324  k = -getD() / scal;
325  else
326  k = 0;
327  H[0] = k * M1.get_X();
328  H[1] = k * M1.get_Y();
329  H[2] = k * M1.get_Z();
330  }
331 
332  return k;
333 }
334 
336 {
337 
338  double k, scal;
339 
340  scal = A * M1[0] + B * M1[1] + C * M1[2];
341  // if (scal != 0)
342  if (std::fabs(scal) > std::numeric_limits<double>::epsilon())
343  k = -getD() / scal;
344  else
345  k = 0;
346  H[0] = k * M1[0];
347  H[1] = k * M1[1];
348  H[2] = k * M1[2];
349 
350  return k;
351 }
352 
362 {
363  // Save current plane parameters
364  double Ao = A, Bo = B, Co = C, Do = D;
365  A = cMo[0][0] * Ao + cMo[0][1] * Bo + cMo[0][2] * Co;
366  B = cMo[1][0] * Ao + cMo[1][1] * Bo + cMo[1][2] * Co;
367  C = cMo[2][0] * Ao + cMo[2][1] * Bo + cMo[2][2] * Co;
368  D = Do - (cMo[0][3] * A + cMo[1][3] * B + cMo[2][3] * C);
369 }
370 
377 VISP_EXPORT std::ostream &operator<<(std::ostream &os, vpPlane &p)
378 {
379  return (os << "(" << p.getA() << "," << p.getB() << "," << p.getC() << "," << p.getD() << ") ");
380 };
friend std::ostream & operator<<(std::ostream &s, const vpArray2D< Type > &A)
Definition: vpArray2D.h:504
Implementation of column vector and the associated operations.
Definition: vpColVector.h:163
static vpColVector cross(const vpColVector &a, const vpColVector &b)
Definition: vpColVector.h:1165
void resize(unsigned int i, bool flagNullify=true)
Definition: vpColVector.h:1049
Implementation of an homogeneous matrix and operations on such kind of matrices.
This class defines the container for a plane geometrical structure.
Definition: vpPlane.h:54
double C
Definition: vpPlane.h:62
vpPlaneFrame
Definition: vpPlane.h:65
@ camera_frame
Definition: vpPlane.h:65
double rayIntersection(const vpPoint &M0, const vpPoint &M1, vpColVector &H) const
Definition: vpPlane.cpp:296
void changeFrame(const vpHomogeneousMatrix &cMo)
Definition: vpPlane.cpp:361
void setA(double a)
Definition: vpPlane.h:80
void projectionPointOnPlan(const vpPoint &P, vpPoint &Pproj) const
Definition: vpPlane.cpp:279
void setD(double d)
Definition: vpPlane.h:86
double A
Definition: vpPlane.h:62
double D
Definition: vpPlane.h:62
double computeZ(double x, double y) const
Definition: vpPlane.cpp:232
double getD() const
Definition: vpPlane.h:106
void setC(double c)
Definition: vpPlane.h:84
double B
Definition: vpPlane.h:62
vpColVector getNormal() const
Definition: vpPlane.cpp:245
double getA() const
Definition: vpPlane.h:100
double getC() const
Definition: vpPlane.h:104
double getB() const
Definition: vpPlane.h:102
vpPlane & operator=(const vpPlane &f)
Definition: vpPlane.cpp:50
vpPlane()
Definition: vpPlane.cpp:63
void init(const vpPoint &P, const vpPoint &Q, const vpPoint &R, vpPlaneFrame frame=camera_frame)
Definition: vpPlane.cpp:165
void setB(double b)
Definition: vpPlane.h:82
double getIntersection(const vpColVector &M1, vpColVector &H) const
Definition: vpPlane.cpp:335
Class that defines a 3D point in the object frame and allows forward projection of a 3D point in the ...
Definition: vpPoint.h:77
double get_oX() const
Get the point oX coordinate in the object frame.
Definition: vpPoint.cpp:450
void set_W(double cW)
Set the point cW coordinate in the camera frame.
Definition: vpPoint.cpp:488
double get_Y() const
Get the point cY coordinate in the camera frame.
Definition: vpPoint.cpp:443
double get_oZ() const
Get the point oZ coordinate in the object frame.
Definition: vpPoint.cpp:454
void set_X(double cX)
Set the point cX coordinate in the camera frame.
Definition: vpPoint.cpp:482
double get_W() const
Get the point cW coordinate in the camera frame.
Definition: vpPoint.cpp:447
void set_Y(double cY)
Set the point cY coordinate in the camera frame.
Definition: vpPoint.cpp:484
double get_Z() const
Get the point cZ coordinate in the camera frame.
Definition: vpPoint.cpp:445
void set_Z(double cZ)
Set the point cZ coordinate in the camera frame.
Definition: vpPoint.cpp:486
double get_oY() const
Get the point oY coordinate in the object frame.
Definition: vpPoint.cpp:452
double get_X() const
Get the point cX coordinate in the camera frame.
Definition: vpPoint.cpp:441