Visual Servoing Platform  version 3.5.0 under development (2022-02-15)
vpMomentObject Class Reference

#include <visp3/core/vpMomentObject.h>

Public Types

enum  vpObjectType { DENSE_FULL_OBJECT = 0, DENSE_POLYGON = 1, DISCRETE = 2 }
 
enum  vpCameraImgBckGrndType { BLACK = 0, WHITE = 1 }
 

Public Member Functions

 vpMomentObject (unsigned int order)
 
 vpMomentObject (const vpMomentObject &srcobj)
 
virtual ~vpMomentObject ()
 
void fromImage (const vpImage< unsigned char > &image, unsigned char threshold, const vpCameraParameters &cam)
 
void fromImage (const vpImage< unsigned char > &image, const vpCameraParameters &cam, vpCameraImgBckGrndType bg_type, bool normalize_with_pix_size=true)
 
void fromVector (std::vector< vpPoint > &points)
 
const std::vector< double > & get () const
 
double get (unsigned int i, unsigned int j) const
 
vpObjectType getType () const
 
unsigned int getOrder () const
 
void init (unsigned int orderinp)
 
void init (const vpMomentObject &objin)
 
void setType (vpObjectType input_type)
 

Static Public Member Functions

static void printWithIndices (const vpMomentObject &momobj, std::ostream &os)
 
static vpMatrix convertTovpMatrix (const vpMomentObject &momobj)
 

Public Attributes

bool flg_normalize_intensity
 

Protected Member Functions

void set (unsigned int i, unsigned int j, const double &value_ij)
 
void cacheValues (std::vector< double > &cache, double x, double y)
 

Protected Attributes

unsigned int order
 
vpObjectType type
 
std::vector< double > values
 

Friends

VISP_EXPORT std::ostream & operator<< (std::ostream &os, const vpMomentObject &v)
 

Detailed Description

Class for generic objects.

It contains all basic moments often described by $m_{ij}$ of order $i+j$ going from $m_{00}$ to the order used as parameter in vpMomentObject() constructor. All other moments implemented in ViSP (gravity center, alpha orientation, centered moments...) use this moment object as a combination of its different values.

When constructing a vpMomentObject() you need first to specify the maximum used moment order as parameter.

Then there are three ways to initialize a vpMomentObject. Firstly using fromImage() you can considerer a dense object O defined by an image. Secondly, as described in fromVector() you can also define a dense object O by a closed contour. In these two cases, 2D basic moments are defined by:

\[m_{ij} = \int \int_{O} x^i y^j dx dy\]

Lastly, as presented in fromVector() you can consider a discrete set of n points. In that last case, the basic moments are defined by

\[m_{ij} = \sum_{k=1}^{n} x_k^i y_k^j \]

With setType() method you can specify the object type.

The implementation is based on the following references [42], [6], [44], [2].

Warning
Be careful with the object order. When you specify a maximum order in the vpMomentObject::vpMomentObject constructor (see its detailed description), it will compute all moment orders up to the order you specified. If you want to access the values $ m_{ij} $ with the vpMomentObject::get method, you can do object.get()[j*(order+1)+i].

A few tips about which orders to use in different situations:

  • moment based visual servoing: use vpMomentObject(6). This will compute moment values up to order 6 which will enable vpFeatureMoments up to order 5 which is the maximum order required for common moments.
  • computing gravity center: use vpMomentObject(1). You only need $ m_{00},m_{01},m_{10} $. You should compute moments up to order 1.
  • computing gravity center interaction matrix with vpFeatureMomentGravityCenter: use vpMomentObject(2). This will compute moment values till order 2 since they are needed for the interaction matrix of vpFeatureMoments of order 1.

The following example shows how to create a moment object from 4 discrete points locate on a plane one meter in front of the camera. It shows also how to get the basic moments that are computed and how to compute other classical moments such as the gravity center or the centered moments.

#include <visp3/core/vpMomentCommon.h>
#include <visp3/core/vpMomentObject.h>
#include <visp3/core/vpPoint.h>
int main()
{
// Define an object as 4 clockwise points on a plane (Z=0)
std::vector<vpPoint> vec_p; // vector that contains the 4 points
vec_p.push_back( vpPoint(-0.2, 0.1, 0.0) ); // values in meters
vec_p.push_back( vpPoint(+0.3, 0.1, 0.0) ); // values in meters
vec_p.push_back( vpPoint(+0.2,-0.1, 0.0) ); // values in meters
vec_p.push_back( vpPoint(-0.2,-0.15, 0.0) ); // values in meters
// These points are observed by a camera
vpHomogeneousMatrix cMo(0, 0, 1, 0, 0, 0); // We set the camera to be 1m far the object
// ... update cMo from an image processing
// Apply the perspective projection to update the points coordinates in the camera plane
for (unsigned int i=0; i<vec_p.size(); ++i)
vec_p[i].project(cMo);
std::cout << "Considered points: " << std::endl;
for(unsigned int i=0; i<vec_p.size(); ++i)
std::cout << "point " << i << ": " << vec_p[i].get_x() << ", " << vec_p[i].get_y() << std::endl;
// Define an image moment object from the previous points
vpMomentObject obj(5); // use moments up to order 5
obj.setType(vpMomentObject::DISCRETE); // initialize the object as constituted by discrete points
obj.fromVector(vec_p); // init the object from the points
// --- Access the computed moments by querying the moment object
// 1. Getting a vector of doubles
std::vector<double> moment = obj.get();
std::cout << std::endl << "Basic moment available (from vector of doubles)" << std::endl;
for(unsigned int k=0; k<=obj.getOrder(); k++) {
for(unsigned int l=0; l<(obj.getOrder()+1)-k; l++) {
std::cout << "m" << l << k << "=" << moment[k*(momobj.getOrder()+1)+ l] << "\t";
}
std::cout<<std::endl;
}
// 2. Print the contents of moment object directly
std::cout << std::endl << "Basic moment available: ";
std::cout << obj << std::endl;
// 3. Directly indexing the moment object
std::cout << std::endl << "Direct acces to some basic moments: " << std::endl;
std::cout << "m00: " << obj.get(0, 0) << std::endl;
std::cout << "m10: " << obj.get(1, 0) << std::endl;
std::cout << "m01: " << obj.get(0, 1) << std::endl;
std::cout << "m22: " << obj.get(2, 2) << std::endl;
std::cout << "m20: " << obj.get(2, 0) << std::endl;
std::cout << "m02: " << obj.get(0, 2) << std::endl;
// Get common moments computed using basic moments
double m00 = vpMomentCommon::getSurface(obj); // surface = m00
double alpha = vpMomentCommon::getAlpha(obj); // orientation
std::vector<double> mu_3 = vpMomentCommon::getMu3(obj); // centered moment up to 3rd order
std::cout << std::endl << "Common moments computed using basic moments:" << std::endl;
std::cout << "Surface: " << m00 << std::endl;
std::cout << "Alpha: " << alpha << std::endl;
std::cout << "Centered moments (mu03, mu12, mu21, mu30): ";
for(unsigned int i=0; i<mu_3.size(); ++i)
std::cout << mu_3[i] << " ";
std::cout << std::endl;
return 0;
}

This example produces the following results:

Considered points:
point 0: -0.2, 0.1
point 1: 0.3, 0.1
point 2: 0.2, -0.1
point 3: -0.2, -0.15
Basic moment available (from vector of doubles):
m00=4 m10=0.1 m20=0.21 m30=0.019 m40=0.0129 m50=0.00211
m01=-0.05 m11=0.02 m21=0.003 m31=0.0023 m41=0.00057
m02=0.0525 m12=-0.0015 m22=0.0026 m32=9e-05
m03=-0.002375 m13=0.000575 m23=-4.5e-05
m04=0.00080625 m14=-7.125e-05
m05=-6.59375e-05
Basic moment available:
4 0.1 0.21 0.019 0.0129 0.00211
-0.05 0.02 0.003 0.0023 0.00057 x
0.0525 -0.0015 0.0026 9e-05 x x
-0.002375 0.000575 -4.5e-05 x x x
0.00080625 -7.125e-05 x x x x
-6.59375e-05 x x x x x
Direct acces to some basic moments:
m00: 4
m10: 0.1
m01: -0.05
m22: 0.0026
m20: 0.21
m02: 0.0525
Common moments computed using basic moments:
Surface: 0.259375
Alpha: 0.133296
Centered moments (mu03, mu12, mu21, mu30): 0.003375 0.0045625 -0.00228125 -0.000421875

Note that in the continuous case, the moment object $m_{00}$ corresponds to the surface $a$ of the object. In the discrete case, it is the number of discrete points $n$.

Examples:
manServoMomentsSimple.cpp, mbot-apriltag-ibvs.cpp, servoBebop2.cpp, servoMomentImage.cpp, servoMomentPoints.cpp, servoMomentPolygon.cpp, testMomentAlpha.cpp, and tutorial-count-coins.cpp.

Definition at line 222 of file vpMomentObject.h.

Member Enumeration Documentation

◆ vpCameraImgBckGrndType

Type of camera image background.

Enumerator
BLACK 

Black background.

WHITE 

Not functional right now.

Definition at line 240 of file vpMomentObject.h.

◆ vpObjectType

Type of object that will be considered.

Enumerator
DENSE_FULL_OBJECT 

A set of points (typically from an image) which are interpreted as being dense.

DENSE_POLYGON 

A set of points (stored in clockwise order) describing a polygon. It will be treated as dense.

DISCRETE 

A cloud of points. Treated as discrete.

Definition at line 228 of file vpMomentObject.h.

Constructor & Destructor Documentation

◆ vpMomentObject() [1/2]

vpMomentObject::vpMomentObject ( unsigned int  max_order)
explicit

Default constructor. Initializes the object with the maximum used order. You cannot use higher order moments than the order of the moment object. The parameter specified is the highest desired included order. All orders up to this values will be computed. In other words, a vpMomentObject will compute all $ m_{ij} $ moments with $ i+j \in [0..order] $.

Parameters
max_order: Maximum reached order (i+j) to be used. All considered i+j will be of order smaller or equal than this parameter. For example if this parameter is 5, all moment values of order 0 to 5 included will be computed.

Mani : outsourced the constructor work to void init (unsigned int orderinp);

Definition at line 492 of file vpMomentObject.cpp.

References init().

◆ vpMomentObject() [2/2]

vpMomentObject::vpMomentObject ( const vpMomentObject srcobj)

Copy constructor

Definition at line 501 of file vpMomentObject.cpp.

References init().

◆ ~vpMomentObject()

vpMomentObject::~vpMomentObject ( )
virtual

Virtual destructor to allow polymorphic usage. For instance,

vpMomentObject* obj = new vpWeightedMomentObject(weightfunc,ORDER);

where vpWeightedMomentObject is child class of vpMomentObject

Nothing to destruct. This will allow for a polymorphic usage For instance,

vpMomentObject* obj = new vpWeightedMomentObject(weightfunc,ORDER);

where vpWeightedMomentObject is child class of vpMomentObject

Definition at line 661 of file vpMomentObject.cpp.

Member Function Documentation

◆ cacheValues()

void vpMomentObject::cacheValues ( std::vector< double > &  cache,
double  x,
double  y 
)
protected

Caching to avoid redundant multiplications.

Parameters
cache: Lookup table that contains the order by order values. For example, if the order is 3, cache will contain:
1 x x^2
y x*y x^2*y
y^2 x*y^2 x^2*y^2
x,y: Coordinates of a point.

Definition at line 111 of file vpMomentObject.cpp.

References order.

Referenced by fromImage(), and fromVector().

◆ convertTovpMatrix()

vpMatrix vpMomentObject::convertTovpMatrix ( const vpMomentObject momobj)
static

Converts the raw moments contained in vpMomentObject to a vpMatrix

Parameters
momobj: A vpMomentObject

This function returns a vpMatrix of size (order+1, order+1).

obj.fromImageWeighted(I, cam, vpMomentObject::BLACK); // cam should have the camera parameters

Instead of accessing the moment m21 as obj.get(2,1), you can now do Mpq[2][1]. This is useful when you want to use the functions available in vpMatrix. One use case i see now is to copy the contents of the matrix to a file or std::cout. For instance, like

// Print to console
Mpq.maplePrint(std::cout);
// Or write to a file
std::ofstream fileMpq("Mpq.csv");
Mpq.maplePrint(fileMpq);

The output can be copied and pasted to MAPLE as a matrix.

Warning
The moments that are not calculated have zeros. For instance, for a vpMomentObject of order 8, the moment m[7,2] is not calculated. It will have 0 by default. User discretion is advised.

Definition at line 640 of file vpMomentObject.cpp.

References get(), getOrder(), and order.

◆ fromImage() [1/2]

void vpMomentObject::fromImage ( const vpImage< unsigned char > &  image,
unsigned char  threshold,
const vpCameraParameters cam 
)

Computes basic moments from an image based on this reference [2].

There is no assumption made about whether the input is dense or discrete but it's more common to use vpMomentObject::DENSE_FULL_OBJECT with this method.

Parameters
image: Image to consider.
threshold: Pixels with a luminance lower than this threshold will be considered.
cam: Camera parameters used to convert pixels coordinates in meters in the image plane.

The code below shows how to use this function.

#include <visp3/core/vpImage.h>
#include <visp3/core/vpMomentObject.h>
int main()
{
vpCameraParameters cam; // Camera parameters used for pixel to
meter conversion vpImage<unsigned char> I(288, 384); // Image used to define
the object
// ... Initialize the image
unsigned char threshold = 128; // Gray level used to define which part of
the image belong to the dense object
vpMomentObject obj(3); // Create an image moment object with 3 as maximum
order obj.fromImage(I, threshold, cam); // Initialize the object from the
image
return 0;
}
Examples:
servoMomentImage.cpp, and testMomentAlpha.cpp.

Definition at line 287 of file vpMomentObject.cpp.

References cacheValues(), vpPixelMeterConversion::convertPoint(), vpCameraParameters::get_px(), vpCameraParameters::get_py(), vpImage< Type >::getCols(), vpImage< Type >::getRows(), order, and values.

◆ fromImage() [2/2]

void vpMomentObject::fromImage ( const vpImage< unsigned char > &  image,
const vpCameraParameters cam,
vpCameraImgBckGrndType  bg_type,
bool  normalize_with_pix_size = true 
)

Computes basic moments from an image based on this reference [2].

Intended to be used by vpMomentObject with DENSE_FULL_OBJECT object type, see setType().

Parameters
image: Grayscale image
cam: Camera parameters (to change to )
bg_type: White/Black background surrounding the image
normalize_with_pix_size: When this flag if set, the moments, after calculation are normalized w.r.t pixel size available from camera parameters.

Definition at line 370 of file vpMomentObject.cpp.

References cacheValues(), vpPixelMeterConversion::convertPoint(), flg_normalize_intensity, vpCameraParameters::get_px(), vpCameraParameters::get_py(), vpImage< Type >::getCols(), vpImage< Type >::getRows(), order, values, and WHITE.

◆ fromVector()

void vpMomentObject::fromVector ( std::vector< vpPoint > &  points)

Computes basic moments from a vector of points. There are two cases:

  1. Dense case specified by setType(vpMomentObject::DENSE_POLYGON):
    • The points are the vertices describing a closed contour polygon.
    • They must be stored in a clockwise order.
    • The first and the last points should be the same to close the contour.
  2. Discrete case specified by setType(vpMomentObject::DISCRETE)
    • The points will be interpreted as a discrete point cloud.
Parameters
points: Vector of points.

The code below shows how to use this function to consider a dense object defined by a closed contour.

#include <visp3/core/vpMomentObject.h>
#include <visp3/core/vpPoint.h>
int main()
{
// Define the contour of an object by a 5 clockwise vertices on a plane
std::vector<vpPoint> vec_p; // vector that contains the vertices of the contour polygon
p.set_x(-0.2); p.set_y(0.1); // coordinates in meters in the image plane (vertex 1)
vec_p.push_back(p);
p.set_x(+0.3); p.set_y(0.1); // coordinates in meters in the image plane (vertex 2)
vec_p.push_back(p);
p.set_x(+0.2); p.set_y(-0.1); // coordinates in meters in the image plane (vertex 3)
vec_p.push_back(p);
p.set_x(-0.2); p.set_y(-0.15); // coordinates in meters in the image plane (vertex 4)
vec_p.push_back(p);
p.set_x(-0.2); p.set_y(0.1); // close the contour (vertex 5 = vertex 1)
vec_p.push_back(p);
vpMomentObject obj(4); // Create an image moment object with 4 as maximum order
obj.setType(vpMomentObject::DENSE_POLYGON); // The object is defined by a countour polygon
obj.fromVector(vec_p); // Init the dense object with the polygon
return 0;
}

This other example shows how to consider an object as a discrete set of four points.

#include <visp3/core/vpMomentObject.h>
#include <visp3/core/vpPoint.h>
int main()
{
// Define 4 discrete points on a plane
std::vector<vpPoint> vec_p; // vector that contains the 4 points
p.set_x(-0.2); p.set_y(0.1); // coordinates in meters in the image plane (point 1)
vec_p.push_back(p);
p.set_x(+0.3); p.set_y(0.1); // coordinates in meters in the image plane (point 2)
vec_p.push_back(p);
p.set_x(+0.2); p.set_y(-0.1); // coordinates in meters in the image plane (point 3)
vec_p.push_back(p);
p.set_x(-0.2); p.set_y(-0.15); // coordinates in meters in the image plane (point 4)
vec_p.push_back(p);
vpMomentObject obj(4); // Create an image moment object with 4 as maximum order
obj.setType(vpMomentObject::DISCRETE); // The object is constituted by discrete points
obj.fromVector(vec_p); // Init the dense object with the points
return 0;
}
Examples:
mbot-apriltag-ibvs.cpp, servoBebop2.cpp, servoMomentPoints.cpp, servoMomentPolygon.cpp, and tutorial-count-coins.cpp.

Definition at line 228 of file vpMomentObject.cpp.

References cacheValues(), DENSE_POLYGON, order, type, and values.

◆ get() [1/2]

const std::vector< double > & vpMomentObject::get ( ) const

Returns all basic moment values $m_{ij}$ with $i \in [0:\mbox{order}]$ and $j \in [0:\mbox{order}]$.

Returns
Vector of moment values. To access $m_{ij}$, you have to read vpMomentObject::get()[j*(order+1)+i].

For example, if the maximal order is 3, the following values are provided:

m00 m10 m20 m01 m11 m21 m02 m12 m12 m30 m03

To access for example to the basic moment m12, you should use this kind of code:

// ... initialise the object using fromVector() or fromImage()
std::vector mij = obj.get();
double m12;
m12 = mij[2*(obj.getOrder()+1)+1]; // i=1 and j=2
Examples:
tutorial-count-coins.cpp.

Definition at line 528 of file vpMomentObject.cpp.

References values.

Referenced by vpMomentArea::compute(), vpMomentCentered::compute(), vpMomentGravityCenter::compute(), vpMomentAreaNormalized::compute(), vpMomentCInvariant::compute(), vpFeatureMomentArea::compute_interaction(), vpFeatureMomentBasic::compute_interaction(), vpFeatureMomentCentered::compute_interaction(), vpFeatureMomentAreaNormalized::compute_interaction(), vpFeatureMomentGravityCenter::compute_interaction(), vpFeatureMomentCInvariant::compute_interaction(), vpFeatureMomentGravityCenterNormalized::compute_interaction(), convertTovpMatrix(), vpMomentBasic::get(), vpMomentGravityCenter::printDependencies(), vpMomentAreaNormalized::printDependencies(), printWithIndices(), and vpMomentCInvariant::vpMomentCInvariant().

◆ get() [2/2]

double vpMomentObject::get ( unsigned int  i,
unsigned int  j 
) const

Returns the basic moment value $m_{ij}$ corresponding to i,j indexes

Parameters
i: First moment index, with $i+j \leq order$.
j: Second moment index, with $i+j \leq order$.

Definition at line 536 of file vpMomentObject.cpp.

References vpException::badValue, getOrder(), order, and values.

◆ getOrder()

◆ getType()

◆ init() [1/2]

void vpMomentObject::init ( unsigned int  orderinp)

Does exactly the work of the default constructor as it existed in the very first version of vpMomentObject.

Definition at line 456 of file vpMomentObject.cpp.

References DENSE_FULL_OBJECT, flg_normalize_intensity, order, type, and values.

Referenced by vpMomentObject().

◆ init() [2/2]

void vpMomentObject::init ( const vpMomentObject objin)

Helper to copy constructor.

Definition at line 468 of file vpMomentObject.cpp.

References flg_normalize_intensity, getOrder(), getType(), order, type, and values.

◆ printWithIndices()

void vpMomentObject::printWithIndices ( const vpMomentObject momobj,
std::ostream &  os 
)
static

Outputs raw moments in indexed form like m[1,1] = value of moment m11

Parameters
momobj: A vpMomentObject
os: Output stream.

Outputs the raw moment values $m_{ij}$ in indexed form. The moment values are same as provided by the operator << which outputs x for uncalculated moments.

Definition at line 600 of file vpMomentObject.cpp.

References get(), and getOrder().

Referenced by vpMomentBasic::printDependencies(), and vpMomentCentered::printDependencies().

◆ set()

void vpMomentObject::set ( unsigned int  i,
unsigned int  j,
const double &  value_ij 
)
protected

Sets the basic moment value $m_{ij}$ corresponding to i,j indexes

Parameters
i: First moment index, with $i+j \leq order$.
j: Second moment index, with $i+j \leq order$.
value_ij: Moment value.

Definition at line 553 of file vpMomentObject.cpp.

References vpException::badValue, getOrder(), order, and values.

◆ setType()

void vpMomentObject::setType ( vpObjectType  input_type)
inline

Specifies the type of the input data.

Parameters
input_type: An input type.
Examples:
mbot-apriltag-ibvs.cpp, servoBebop2.cpp, servoMomentImage.cpp, servoMomentPoints.cpp, servoMomentPolygon.cpp, testMomentAlpha.cpp, and tutorial-count-coins.cpp.

Definition at line 297 of file vpMomentObject.h.

Friends And Related Function Documentation

◆ operator<<

VISP_EXPORT std::ostream& operator<< ( std::ostream &  os,
const vpMomentObject v 
)
friend

Outputs the basic moment's values $m_{ij}$ to a stream presented as a matrix. The first line corresponds to $m_{0[0:order]}$, the second one to $m_{1[0:order]}$ Values in table corresponding to a higher order are marked with an "x" and not computed.

For example, if the maximal order is 3, the following values are provided:

m00 m10 m20 m30
m01 m11 m21 x
m02 m12 x x
m03 x x x

Definition at line 577 of file vpMomentObject.cpp.

Member Data Documentation

◆ flg_normalize_intensity

bool vpMomentObject::flg_normalize_intensity

Definition at line 245 of file vpMomentObject.h.

Referenced by fromImage(), and init().

◆ order

unsigned int vpMomentObject::order
protected

Definition at line 306 of file vpMomentObject.h.

Referenced by cacheValues(), convertTovpMatrix(), fromImage(), fromVector(), get(), init(), and set().

◆ type

vpObjectType vpMomentObject::type
protected

Definition at line 307 of file vpMomentObject.h.

Referenced by fromVector(), and init().

◆ values

std::vector<double> vpMomentObject::values
protected

Definition at line 308 of file vpMomentObject.h.

Referenced by fromImage(), fromVector(), get(), init(), and set().