samples/cpp/tutorial_code/ml/introduction_to_pca/introduction_to_pca.cpp#

Check the corresponding tutorial for more details

  1/**
  2 * @file introduction_to_pca.cpp
  3 * @brief This program demonstrates how to use OpenCV PCA to extract the orientation of an object
  4 * @author OpenCV team
  5 */
  6
  7#include "opencv2/core.hpp"
  8#include "opencv2/imgproc.hpp"
  9#include "opencv2/geometry.hpp"
 10#include "opencv2/highgui.hpp"
 11#include <iostream>
 12
 13using namespace std;
 14using namespace cv;
 15
 16// Function declarations
 17void drawAxis(Mat&, Point, Point, Scalar, const float);
 18double getOrientation(const vector<Point> &, Mat&);
 19
 20/**
 21 * @function drawAxis
 22 */
 23void drawAxis(Mat& img, Point p, Point q, Scalar colour, const float scale = 0.2)
 24{
 25    //! [visualization1]
 26    double angle = atan2( (double) p.y - q.y, (double) p.x - q.x ); // angle in radians
 27    double hypotenuse = sqrt( (double) (p.y - q.y) * (p.y - q.y) + (p.x - q.x) * (p.x - q.x));
 28
 29    // Here we lengthen the arrow by a factor of scale
 30    q.x = (int) (p.x - scale * hypotenuse * cos(angle));
 31    q.y = (int) (p.y - scale * hypotenuse * sin(angle));
 32    line(img, p, q, colour, 1, LINE_AA);
 33
 34    // create the arrow hooks
 35    p.x = (int) (q.x + 9 * cos(angle + CV_PI / 4));
 36    p.y = (int) (q.y + 9 * sin(angle + CV_PI / 4));
 37    line(img, p, q, colour, 1, LINE_AA);
 38
 39    p.x = (int) (q.x + 9 * cos(angle - CV_PI / 4));
 40    p.y = (int) (q.y + 9 * sin(angle - CV_PI / 4));
 41    line(img, p, q, colour, 1, LINE_AA);
 42    //! [visualization1]
 43}
 44
 45/**
 46 * @function getOrientation
 47 */
 48double getOrientation(const vector<Point> &pts, Mat &img)
 49{
 50    //! [pca]
 51    //Construct a buffer used by the pca analysis
 52    int sz = static_cast<int>(pts.size());
 53    Mat data_pts = Mat(sz, 2, CV_64F);
 54    for (int i = 0; i < data_pts.rows; i++)
 55    {
 56        data_pts.at<double>(i, 0) = pts[i].x;
 57        data_pts.at<double>(i, 1) = pts[i].y;
 58    }
 59
 60    //Perform PCA analysis
 61    PCA pca_analysis(data_pts, Mat(), PCA::DATA_AS_ROW);
 62
 63    //Store the center of the object
 64    Point cntr = Point(static_cast<int>(pca_analysis.mean.at<double>(0, 0)),
 65                      static_cast<int>(pca_analysis.mean.at<double>(0, 1)));
 66
 67    //Store the eigenvalues and eigenvectors
 68    vector<Point2d> eigen_vecs(2);
 69    vector<double> eigen_val(2);
 70    for (int i = 0; i < 2; i++)
 71    {
 72        eigen_vecs[i] = Point2d(pca_analysis.eigenvectors.at<double>(i, 0),
 73                                pca_analysis.eigenvectors.at<double>(i, 1));
 74
 75        eigen_val[i] = pca_analysis.eigenvalues.at<double>(i);
 76    }
 77    //! [pca]
 78
 79    //! [visualization]
 80    // Draw the principal components
 81    circle(img, cntr, 3, Scalar(255, 0, 255), 2);
 82    Point p1 = cntr + 0.02 * Point(static_cast<int>(eigen_vecs[0].x * eigen_val[0]), static_cast<int>(eigen_vecs[0].y * eigen_val[0]));
 83    Point p2 = cntr - 0.02 * Point(static_cast<int>(eigen_vecs[1].x * eigen_val[1]), static_cast<int>(eigen_vecs[1].y * eigen_val[1]));
 84    drawAxis(img, cntr, p1, Scalar(0, 255, 0), 1);
 85    drawAxis(img, cntr, p2, Scalar(255, 255, 0), 5);
 86
 87    double angle = atan2(eigen_vecs[0].y, eigen_vecs[0].x); // orientation in radians
 88    //! [visualization]
 89
 90    return angle;
 91}
 92
 93/**
 94 * @function main
 95 */
 96int main(int argc, char** argv)
 97{
 98    //! [pre-process]
 99    // Load image
100    CommandLineParser parser(argc, argv, "{@input | pca_test1.jpg | input image}");
101    parser.about( "This program demonstrates how to use OpenCV PCA to extract the orientation of an object.\n" );
102    parser.printMessage();
103
104    Mat src = imread( samples::findFile( parser.get<String>("@input") ) );
105
106    // Check if image is loaded successfully
107    if(src.empty())
108    {
109        cout << "Problem loading image!!!" << endl;
110        return EXIT_FAILURE;
111    }
112
113    imshow("src", src);
114
115    // Convert image to grayscale
116    Mat gray;
117    cvtColor(src, gray, COLOR_BGR2GRAY);
118
119    // Convert image to binary
120    Mat bw;
121    threshold(gray, bw, 50, 255, THRESH_BINARY | THRESH_OTSU);
122    //! [pre-process]
123
124    //! [contours]
125    // Find all the contours in the thresholded image
126    vector<vector<Point> > contours;
127    findContours(bw, contours, RETR_LIST, CHAIN_APPROX_NONE);
128
129    for (size_t i = 0; i < contours.size(); i++)
130    {
131        // Calculate the area of each contour
132        double area = contourArea(contours[i]);
133        // Ignore contours that are too small or too large
134        if (area < 1e2 || 1e5 < area) continue;
135
136        // Draw each contour only for visualisation purposes
137        drawContours(src, contours, static_cast<int>(i), Scalar(0, 0, 255), 2);
138        // Find the orientation of each shape
139        getOrientation(contours[i], src);
140    }
141    //! [contours]
142
143    imshow("output", src);
144
145    waitKey();
146    return EXIT_SUCCESS;
147}