Skip to content

O
opencv_contrib
Commit 9a52f462 authored by biagio montesano's avatar biagio montesano
Browse files
Options

Created nested classes. Removed some useless headers.

parent e39065bf
Hide whitespace changes
Inline Side-by-side
Showing with 2802 additions and 3254 deletions
modules/line_descriptor/include/opencv2/line_descriptor/array32.hpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __OPENCV_ARRAY32_HPP
#define __OPENCV_ARRAY32_HPP
#ifdef _WIN32
#pragma warning( disable : 4267 )
#endif
#include "types.hpp"
class Array32
{
private:
static double ARRAY_RESIZE_FACTOR;
static double ARRAY_RESIZE_ADD_FACTOR;
public:
/* set ARRAY_RESIZE_FACTOR */
static void setArrayResizeFactor( double arf );
/* constructor */
Array32();
/* destructor */
~Array32();
/* cleaning function used in destructor */
void cleanup();
/* push data */
void push( UINT32 data );
/* insert data at given index */
void insert( UINT32 index, UINT32 data );
/* return data */
UINT32* data();
/* return data size */
UINT32 size();
/* return capacity */
UINT32 capacity();
/* definition of operator = */
void operator=( const Array32& );
/* print data */
void print();
/* initializer */
void init( int size );
/* data */
UINT32 *arr;
};
#endif
modules/line_descriptor/include/opencv2/line_descriptor/bucket_group.hpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __OPENCV_BUCKET_GROUP_HPP
#define __OPENCV_BUCKET_GROUP_HPP
#ifdef _WIN32
#pragma warning( disable : 4267 )
#endif
#include "types.hpp"
#include "array32.hpp"
#include "bitarray.hpp"
class BucketGroup
{
public:
/* constructor */
BucketGroup();
/* destructor */
~BucketGroup();
/* insert data into the bucket */
void insert( int subindex, UINT32 data );
/* perform a query to the bucket */
UINT32* query( int subindex, int *size );
/* data fields */
UINT32 empty;
Array32 *group;
};
#endif
modules/line_descriptor/include/opencv2/line_descriptor/descriptor.hpp
View file @ 9a52f462
...@@ -42,23 +42,36 @@ ...@@ -42,23 +42,36 @@
#ifndef __OPENCV_DESCRIPTOR_HPP__ #ifndef __OPENCV_DESCRIPTOR_HPP__
#define __OPENCV_DESCRIPTOR_HPP__ #define __OPENCV_DESCRIPTOR_HPP__
#include "line_structure.hpp"
#include "array32.hpp"
#include "bitarray.hpp"
#include "bitops.hpp"
#include "bucket_group.hpp"
#include "mihasher.hpp"
#include "sparse_hashtable.hpp"
#include "types.hpp"
#include "ed_line_detector.hpp"
#include <map> #include <map>
#include <vector>
#include <list>
#include <inttypes.h>
#include <stdio.h>
#include <iostream>
#include "opencv2/core/utility.hpp"
#include "opencv2/core/private.hpp"
#include <opencv2/imgproc.hpp>
#include <opencv2/features2d.hpp>
#include <opencv2/highgui.hpp>
#include "opencv2/core.hpp"
/* define data types */
typedef uint64_t UINT64;
typedef uint32_t UINT32;
typedef uint16_t UINT16;
typedef uint8_t UINT8;
/* define constants */
#define UINT64_1 ((UINT64)0x01)
#define UINT32_1 ((UINT32)0x01)
namespace cv namespace cv
{ {
namespace line_descriptor namespace line_descriptor
{ {
//CV_EXPORTS bool initModule_line_descriptor(); CV_EXPORTS bool initModule_line_descriptor();
struct CV_EXPORTS KeyLine struct CV_EXPORTS KeyLine
{ {
...@@ -228,176 +241,906 @@ class CV_EXPORTS BinaryDescriptor : public Algorithm ...@@ -228,176 +241,906 @@ class CV_EXPORTS BinaryDescriptor : public Algorithm
AlgorithmInfo* info() const; AlgorithmInfo* info() const;
private: private:
/* compute Gaussian pyramids */ /* struct to represent lines extracted from an octave */
void computeGaussianPyramid( const Mat& image, const int numOctaves ); struct OctaveLine
{
unsigned int octaveCount; //the octave which this line is detected
unsigned int lineIDInOctave; //the line ID in that octave image
unsigned int lineIDInScaleLineVec; //the line ID in Scale line vector
float lineLength; //the length of line in original image scale
};
// A 2D line (normal equation parameters).
struct SingleLine
{
//note: rho and theta are based on coordinate origin, i.e. the top-left corner of image
double rho; //unit: pixel length
double theta; //unit: rad
double linePointX; // = rho * cos(theta);
double linePointY; // = rho * sin(theta);
//for EndPoints, the coordinate origin is the top-left corner of image.
double startPointX;
double startPointY;
double endPointX;
double endPointY;
//direction of a line, the angle between positive line direction (dark side is in the left) and positive X axis.
double direction;
//mean gradient magnitude
double gradientMagnitude;
//mean gray value of pixels in dark side of line
double darkSideGrayValue;
//mean gray value of pixels in light side of line
double lightSideGrayValue;
//the length of line
double lineLength;
//the width of line;
double width;
//number of pixels
int numOfPixels;
//the decriptor of line
std::vector<double> descriptor;
};
// Specifies a vector of lines.
typedef std::vector<SingleLine> Lines_list;
struct OctaveSingleLine
{
/*endPoints, the coordinate origin is the top-left corner of the original image.
*startPointX = sPointInOctaveX * (factor)^octaveCount; */
float startPointX;
float startPointY;
float endPointX;
float endPointY;
//endPoints, the coordinate origin is the top-left corner of the octave image.
float sPointInOctaveX;
float sPointInOctaveY;
float ePointInOctaveX;
float ePointInOctaveY;
//direction of a line, the angle between positive line direction (dark side is in the left) and positive X axis.
float direction;
//the summation of gradient magnitudes of pixels on lines
float salience;
//the length of line
float lineLength;
//number of pixels
unsigned int numOfPixels;
//the octave which this line is detected
unsigned int octaveCount;
//the decriptor of line
std::vector<float> descriptor;
};
struct Pixel
{
unsigned int x; //X coordinate
unsigned int y; //Y coordinate
};
struct EdgeChains
{
std::vector<unsigned int> xCors; //all the x coordinates of edge points
std::vector<unsigned int> yCors; //all the y coordinates of edge points
std::vector<unsigned int> sId; //the start index of each edge in the coordinate arrays
unsigned int numOfEdges; //the number of edges whose length are larger than minLineLen; numOfEdges < sId.size;
};
struct LineChains
{
std::vector<unsigned int> xCors; //all the x coordinates of line points
std::vector<unsigned int> yCors; //all the y coordinates of line points
std::vector<unsigned int> sId; //the start index of each line in the coordinate arrays
unsigned int numOfLines; //the number of lines whose length are larger than minLineLen; numOfLines < sId.size;
};
typedef std::list<Pixel> PixelChain; //each edge is a pixel chain
struct EDLineParam
{
int ksize;
float sigma;
float gradientThreshold;
float anchorThreshold;
int scanIntervals;
int minLineLen;
double lineFitErrThreshold;
};
/* compute Sobel's derivatives */ #define RELATIVE_ERROR_FACTOR 100.0
void computeSobel( const Mat& image, const int numOctaves ); #define MLN10 2.30258509299404568402
#define log_gamma(x) ((x)>15.0?log_gamma_windschitl(x):log_gamma_lanczos(x))
/* This class is used to detect lines from input image.
* First, edges are extracted from input image following the method presented in Cihan Topal and
* Cuneyt Akinlar's paper:"Edge Drawing: A Heuristic Approach to Robust Real-Time Edge Detection", 2010.
* Then, lines are extracted from the edge image following the method presented in Cuneyt Akinlar and
* Cihan Topal's paper:"EDLines: A real-time line segment detector with a false detection control", 2011
* PS: The linking step of edge detection has a little bit difference with the Edge drawing algorithm
* described in the paper. The edge chain doesn't stop when the pixel direction is changed.
*/
class EDLineDetector
{
public:
EDLineDetector();
EDLineDetector( EDLineParam param );
~EDLineDetector();
/*extract edges from image
*image: In, gray image;
*edges: Out, store the edges, each edge is a pixel chain
*return -1: error happen
*/
int EdgeDrawing( cv::Mat &image, EdgeChains &edgeChains );
/*extract lines from image
*image: In, gray image;
*lines: Out, store the extracted lines,
*return -1: error happen
*/
int EDline( cv::Mat &image, LineChains &lines );
/* extract line from image, and store them */
int EDline( cv::Mat &image );
cv::Mat dxImg_; //store the dxImg;
cv::Mat dyImg_; //store the dyImg;
cv::Mat gImgWO_; //store the gradient image without threshold;
LineChains lines_; //store the detected line chains;
//store the line Equation coefficients, vec3=[w1,w2,w3] for line w1*x + w2*y + w3=0;
std::vector<std::vector<double> > lineEquations_;
//store the line endpoints, [x1,y1,x2,y3]
std::vector<std::vector<float> > lineEndpoints_;
//store the line direction
std::vector<float> lineDirection_;
//store the line salience, which is the summation of gradients of pixels on line
std::vector<float> lineSalience_;
// image sizes
unsigned int imageWidth;
unsigned int imageHeight;
/*The threshold of line fit error;
*If lineFitErr is large than this threshold, then
*the pixel chain is not accepted as a single line segment.*/
double lineFitErrThreshold_;
/*the threshold of pixel gradient magnitude.
*Only those pixel whose gradient magnitude are larger than this threshold will be
*taken as possible edge points. Default value is 36*/
short gradienThreshold_;
/*If the pixel's gradient value is bigger than both of its neighbors by a
*certain threshold (ANCHOR_THRESHOLD), the pixel is marked to be an anchor.
*Default value is 8*/
unsigned char anchorThreshold_;
/* conversion of an LBD descriptor to its binary representation */ /*anchor testing can be performed at different scan intervals, i.e.,
unsigned char binaryConversion( float* f1, float* f2 ); *every row/column, every second row/column etc.
*Default value is 2*/
unsigned int scanIntervals_;
/* compute LBD descriptors using EDLine extractor */ int minLineLen_; //minimal acceptable line length
int computeLBD( ScaleLines &keyLines, bool useDetectionData = false );
/* gathers lines in groups using EDLine extractor. private:
Each group contains the same line, detected in different octaves */ void InitEDLine_();
int OctaveKeyLines( cv::Mat& image, ScaleLines &keyLines );
/* the local gaussian coefficient applied to the orthogonal line direction within each band */ /*For an input edge chain, find the best fit line, the default chain length is minLineLen_
std::vector<double> gaussCoefL_; *xCors: In, pointer to the X coordinates of pixel chain;
*yCors: In, pointer to the Y coordinates of pixel chain;
*offsetS:In, start index of this chain in vector;
*lineEquation: Out, [a,b] which are the coefficient of lines y=ax+b(horizontal) or x=ay+b(vertical);
*return: line fit error; -1:error happens;
*/
double LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, std::vector<double> &lineEquation );
/*For an input pixel chain, find the best fit line. Only do the update based on new points.
*For A*x=v, Least square estimation of x = Inv(A^T * A) * (A^T * v);
*If some new observations are added, i.e, [A; A'] * x = [v; v'],
*then x' = Inv(A^T * A + (A')^T * A') * (A^T * v + (A')^T * v');
*xCors: In, pointer to the X coordinates of pixel chain;
*yCors: In, pointer to the Y coordinates of pixel chain;
*offsetS:In, start index of this chain in vector;
*newOffsetS: In, start index of extended part;
*offsetE:In, end index of this chain in vector;
*lineEquation: Out, [a,b] which are the coefficient of lines y=ax+b(horizontal) or x=ay+b(vertical);
*return: line fit error; -1:error happens;
*/
double LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int newOffsetS, unsigned int offsetE,
std::vector<double> &lineEquation );
/* Validate line based on the Helmholtz principle, which basically states that
* for a structure to be perceptually meaningful, the expectation of this structure
* by chance must be very low.
*/
bool LineValidation_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int offsetE, std::vector<double> &lineEquation,
float &direction );
bool bValidate_; //flag to decide whether line will be validated
int ksize_; //the size of Gaussian kernel: ksize X ksize, default value is 5.
float sigma_; //the sigma of Gaussian kernal, default value is 1.0.
/*For example, there two edges in the image:
*edge1 = [(7,4), (8,5), (9,6),| (10,7)|, (11, 8), (12,9)] and
*edge2 = [(14,9), (15,10), (16,11), (17,12),| (18, 13)|, (19,14)] ; then we store them as following:
*pFirstPartEdgeX_ = [10, 11, 12, 18, 19];//store the first part of each edge[from middle to end]
*pFirstPartEdgeY_ = [7, 8, 9, 13, 14];
*pFirstPartEdgeS_ = [0,3,5];// the index of start point of first part of each edge
*pSecondPartEdgeX_ = [10, 9, 8, 7, 18, 17, 16, 15, 14];//store the second part of each edge[from middle to front]
*pSecondPartEdgeY_ = [7, 6, 5, 4, 13, 12, 11, 10, 9];//anchor points(10, 7) and (18, 13) are stored again
*pSecondPartEdgeS_ = [0, 4, 9];// the index of start point of second part of each edge
*This type of storage order is because of the order of edge detection process.
*For each edge, start from one anchor point, first go right, then go left or first go down, then go up*/
//store the X coordinates of the first part of the pixels for chains
unsigned int *pFirstPartEdgeX_;
//store the Y coordinates of the first part of the pixels for chains
unsigned int *pFirstPartEdgeY_;
//store the start index of every edge chain in the first part arrays
unsigned int *pFirstPartEdgeS_;
//store the X coordinates of the second part of the pixels for chains
unsigned int *pSecondPartEdgeX_;
//store the Y coordinates of the second part of the pixels for chains
unsigned int *pSecondPartEdgeY_;
//store the start index of every edge chain in the second part arrays
unsigned int *pSecondPartEdgeS_;
//store the X coordinates of anchors
unsigned int *pAnchorX_;
//store the Y coordinates of anchors
unsigned int *pAnchorY_;
//edges
cv::Mat edgeImage_;
cv::Mat gImg_; //store the gradient image;
cv::Mat dirImg_; //store the direction image
double logNT_;
cv::Mat_<float> ATA; //the previous matrix of A^T * A;
cv::Mat_<float> ATV; //the previous vector of A^T * V;
cv::Mat_<float> fitMatT; //the matrix used in line fit function;
cv::Mat_<float> fitVec; //the vector used in line fit function;
cv::Mat_<float> tempMatLineFit; //the matrix used in line fit function;
cv::Mat_<float> tempVecLineFit; //the vector used in line fit function;
/** Compare doubles by relative error.
The resulting rounding error after floating point computations
depend on the specific operations done. The same number computed by
different algorithms could present different rounding errors. For a
useful comparison, an estimation of the relative rounding error
should be considered and compared to a factor times EPS. The factor
should be related to the accumulated rounding error in the chain of
computation. Here, as a simplification, a fixed factor is used.
*/
static int double_equal( double a, double b )
{
double abs_diff, aa, bb, abs_max;
/* trivial case */
if( a == b )
return true;
abs_diff = fabs( a - b );
aa = fabs( a );
bb = fabs( b );
abs_max = aa > bb ? aa : bb;
/* DBL_MIN is the smallest normalized number, thus, the smallest
number whose relative error is bounded by DBL_EPSILON. For
smaller numbers, the same quantization steps as for DBL_MIN
are used. Then, for smaller numbers, a meaningful "relative"
error should be computed by dividing the difference by DBL_MIN. */
if( abs_max < DBL_MIN )
abs_max = DBL_MIN;
/* equal if relative error <= factor x eps */
return ( abs_diff / abs_max ) <= ( RELATIVE_ERROR_FACTOR * DBL_EPSILON );
}
/** Computes the natural logarithm of the absolute value of
the gamma function of x using the Lanczos approximation.
See http://www.rskey.org/gamma.htm
The formula used is
@f[
\Gamma(x) = \frac{ \sum_{n=0}^{N} q_n x^n }{ \Pi_{n=0}^{N} (x+n) }
(x+5.5)^{x+0.5} e^{-(x+5.5)}
@f]
so
@f[
\log\Gamma(x) = \log\left( \sum_{n=0}^{N} q_n x^n \right)
+ (x+0.5) \log(x+5.5) - (x+5.5) - \sum_{n=0}^{N} \log(x+n)
@f]
and
q0 = 75122.6331530,
q1 = 80916.6278952,
q2 = 36308.2951477,
q3 = 8687.24529705,
q4 = 1168.92649479,
q5 = 83.8676043424,
q6 = 2.50662827511.
*/
static double log_gamma_lanczos( double x )
{
static double q[7] =
{ 75122.6331530, 80916.6278952, 36308.2951477, 8687.24529705, 1168.92649479, 83.8676043424, 2.50662827511 };
double a = ( x + 0.5 ) * log( x + 5.5 ) - ( x + 5.5 );
double b = 0.0;
int n;
for ( n = 0; n < 7; n++ )
{
a -= log( x + (double) n );
b += q[n] * pow( x, (double) n );
}
return a + log( b );
}
/** Computes the natural logarithm of the absolute value of
the gamma function of x using Windschitl method.
See http://www.rskey.org/gamma.htm
The formula used is
@f[
\Gamma(x) = \sqrt{\frac{2\pi}{x}} \left( \frac{x}{e}
\sqrt{ x\sinh(1/x) + \frac{1}{810x^6} } \right)^x
@f]
so
@f[
\log\Gamma(x) = 0.5\log(2\pi) + (x-0.5)\log(x) - x
+ 0.5x\log\left( x\sinh(1/x) + \frac{1}{810x^6} \right).
@f]
This formula is a good approximation when x > 15.
*/
static double log_gamma_windschitl( double x )
{
return 0.918938533204673 + ( x - 0.5 ) * log( x ) - x + 0.5 * x * log( x * sinh( 1 / x ) + 1 / ( 810.0 * pow( x, 6.0 ) ) );
}
/** Computes -log10(NFA).
NFA stands for Number of False Alarms:
@f[
\mathrm{NFA} = NT \cdot B(n,k,p)
@f]
- NT - number of tests
- B(n,k,p) - tail of binomial distribution with parameters n,k and p:
@f[
B(n,k,p) = \sum_{j=k}^n
\left(\begin{array}{c}n\\j\end{array}\right)
p^{j} (1-p)^{n-j}
@f]
The value -log10(NFA) is equivalent but more intuitive than NFA:
- -1 corresponds to 10 mean false alarms
- 0 corresponds to 1 mean false alarm
- 1 corresponds to 0.1 mean false alarms
- 2 corresponds to 0.01 mean false alarms
- ...
Used this way, the bigger the value, better the detection,
and a logarithmic scale is used.
@param n,k,p binomial parameters.
@param logNT logarithm of Number of Tests
The computation is based in the gamma function by the following
relation:
@f[
\left(\begin{array}{c}n\\k\end{array}\right)
= \frac{ \Gamma(n+1) }{ \Gamma(k+1) \cdot \Gamma(n-k+1) }.
@f]
We use efficient algorithms to compute the logarithm of
the gamma function.
To make the computation faster, not all the sum is computed, part
of the terms are neglected based on a bound to the error obtained
(an error of 10% in the result is accepted).
*/
static double nfa( int n, int k, double p, double logNT )
{
double tolerance = 0.1; /* an error of 10% in the result is accepted */
double log1term, term, bin_term, mult_term, bin_tail, err, p_term;
int i;
/* check parameters */
if( n < 0 || k < 0 || k > n || p <= 0.0 || p >= 1.0 )
{
std::cout << "nfa: wrong n, k or p values." << std::endl;
exit( 0 );
}
/* trivial cases */
if( n == 0 || k == 0 )
return -logNT;
if( n == k )
return -logNT - (double) n * log10( p );
/* probability term */
p_term = p / ( 1.0 - p );
/* compute the first term of the series */
/*
binomial_tail(n,k,p) = sum_{i=k}^n bincoef(n,i) * p^i * (1-p)^{n-i}
where bincoef(n,i) are the binomial coefficients.
But
bincoef(n,k) = gamma(n+1) / ( gamma(k+1) * gamma(n-k+1) ).
We use this to compute the first term. Actually the log of it.
*/
log1term = log_gamma( (double) n + 1.0 )- log_gamma( (double ) k + 1.0 )- log_gamma( (double ) ( n - k ) + 1.0 )
+ (double) k * log( p )
+ (double) ( n - k ) * log( 1.0 - p );
term = exp( log1term );
/* in some cases no more computations are needed */
if( double_equal( term, 0.0 ) )
{ /* the first term is almost zero */
if( (double) k > (double) n * p ) /* at begin or end of the tail? */
return -log1term / MLN10 - logNT; /* end: use just the first term */
else
return -logNT; /* begin: the tail is roughly 1 */
}
/* compute more terms if needed */
bin_tail = term;
for ( i = k + 1; i <= n; i++ )
{
/* As
term_i = bincoef(n,i) * p^i * (1-p)^(n-i)
and
bincoef(n,i)/bincoef(n,i-1) = n-i+1 / i,
then,
term_i / term_i-1 = (n-i+1)/i * p/(1-p)
and
term_i = term_i-1 * (n-i+1)/i * p/(1-p).
p/(1-p) is computed only once and stored in 'p_term'.
*/
bin_term = (double) ( n - i + 1 ) / (double) i;
mult_term = bin_term * p_term;
term *= mult_term;
bin_tail += term;
if( bin_term < 1.0 )
{
/* When bin_term<1 then mult_term_j<mult_term_i for j>i.
Then, the error on the binomial tail when truncated at
the i term can be bounded by a geometric series of form
term_i * sum mult_term_i^j. */
err = term * ( ( 1.0 - pow( mult_term, (double) ( n - i + 1 ) ) ) / ( 1.0 - mult_term ) - 1.0 );
/* One wants an error at most of tolerance*final_result, or:
tolerance * abs(-log10(bin_tail)-logNT).
Now, the error that can be accepted on bin_tail is
given by tolerance*final_result divided by the derivative
of -log10(x) when x=bin_tail. that is:
tolerance * abs(-log10(bin_tail)-logNT) / (1/bin_tail)
Finally, we truncate the tail if the error is less than:
tolerance * abs(-log10(bin_tail)-logNT) * bin_tail */
if( err < tolerance * fabs( -log10( bin_tail ) - logNT ) * bin_tail )
break;
}
}
return -log10( bin_tail ) - logNT;
}
};
/* the global gaussian coefficient applied to each row within line support region */ // Specifies a vector of lines.
std::vector<double> gaussCoefG_; typedef std::vector<OctaveSingleLine> LinesVec;
/* descriptor parameters */ // each element in ScaleLines is a vector of lines
Params params; // which corresponds the same line detected in different octave images.
typedef std::vector<LinesVec> ScaleLines;
/* vector of sizes of downsampled and blurred images */ /* compute Gaussian pyramids */
std::vector<cv::Size> images_sizes; void computeGaussianPyramid( const Mat& image, const int numOctaves );
/*For each octave of image, we define an EDLineDetector, because we can get gradient images (dxImg, dyImg, gImg) /* compute Sobel's derivatives */
*from the EDLineDetector class without extra computation cost. Another reason is that, if we use void computeSobel( const Mat& image, const int numOctaves );
*a single EDLineDetector to detect lines in different octave of images, then we need to allocate and release
*memory for gradient images (dxImg, dyImg, gImg) repeatedly for their varying size*/
std::vector<Ptr<EDLineDetector> > edLineVec_;
/* Sobel's derivatives */ /* conversion of an LBD descriptor to its binary representation */
std::vector<cv::Mat> dxImg_vector, dyImg_vector; unsigned char binaryConversion( float* f1, float* f2 );
/* Gaussian pyramid */ /* compute LBD descriptors using EDLine extractor */
std::vector<cv::Mat> octaveImages; int computeLBD( ScaleLines &keyLines, bool useDetectionData = false );
/* gathers lines in groups using EDLine extractor.
Each group contains the same line, detected in different octaves */
int OctaveKeyLines( cv::Mat& image, ScaleLines &keyLines );
/* the local gaussian coefficient applied to the orthogonal line direction within each band */
std::vector<double> gaussCoefL_;
/* the global gaussian coefficient applied to each row within line support region */
std::vector<double> gaussCoefG_;
/* descriptor parameters */
Params params;
/* vector of sizes of downsampled and blurred images */
std::vector<cv::Size> images_sizes;
/*For each octave of image, we define an EDLineDetector, because we can get gradient images (dxImg, dyImg, gImg)
*from the EDLineDetector class without extra computation cost. Another reason is that, if we use
*a single EDLineDetector to detect lines in different octave of images, then we need to allocate and release
*memory for gradient images (dxImg, dyImg, gImg) repeatedly for their varying size*/
std::vector<Ptr<EDLineDetector> > edLineVec_;
/* Sobel's derivatives */
std::vector<cv::Mat> dxImg_vector, dyImg_vector;
/* Gaussian pyramid */
std::vector<cv::Mat> octaveImages;
}; };
class CV_EXPORTS LSDDetector : public Algorithm class CV_EXPORTS LSDDetector : public Algorithm
{ {
public: public:
/* constructor */ /* constructor */
/*CV_WRAP*/ /*CV_WRAP*/
LSDDetector() LSDDetector()
{ {
} }
; ;
/* constructor with smart pointer */ /* constructor with smart pointer */
/*CV_WRAP*/ /*CV_WRAP*/
static Ptr<LSDDetector> createLSDDetector(); static Ptr<LSDDetector> createLSDDetector();
/* requires line detection (only one image) */ /* requires line detection (only one image) */
/*CV_WRAP*/ /*CV_WRAP*/
void detect( const Mat& image, CV_OUT std::vector<KeyLine>& keypoints, int scale, int numOctaves, const Mat& mask = Mat() ); void detect( const Mat& image, CV_OUT std::vector<KeyLine>& keypoints, int scale, int numOctaves, const Mat& mask = Mat() );
/* requires line detection (more than one image) */ /* requires line detection (more than one image) */
/*CV_WRAP*/ /*CV_WRAP*/
void detect( const std::vector<Mat>& images, std::vector<std::vector<KeyLine> >& keylines, int scale, int numOctaves, void detect( const std::vector<Mat>& images, std::vector<std::vector<KeyLine> >& keylines, int scale, int numOctaves,
const std::vector<Mat>& masks = std::vector<Mat>() ) const; const std::vector<Mat>& masks = std::vector<Mat>() ) const;
private: private:
/* compute Gaussian pyramid of input image */ /* compute Gaussian pyramid of input image */
void computeGaussianPyramid( const Mat& image, int numOctaves, int scale ); void computeGaussianPyramid( const Mat& image, int numOctaves, int scale );
/* implementation of line detection */ /* implementation of line detection */
void detectImpl( const Mat& imageSrc, std::vector<KeyLine>& keylines, int numOctaves, int scale, const Mat& mask ) const; void detectImpl( const Mat& imageSrc, std::vector<KeyLine>& keylines, int numOctaves, int scale, const Mat& mask ) const;
/* matrices for Gaussian pyramids */ /* matrices for Gaussian pyramids */
std::vector<cv::Mat> gaussianPyrs; std::vector<cv::Mat> gaussianPyrs;
protected: protected:
/* function inherited from Algorithm */ /* function inherited from Algorithm */
AlgorithmInfo* info() const; AlgorithmInfo* info() const;
}; };
class CV_EXPORTS BinaryDescriptorMatcher : public Algorithm class CV_EXPORTS BinaryDescriptorMatcher : public Algorithm
{ {
public: public:
/* for every input descriptor, /* for every input descriptor,
find the best matching one (for a pair of images) */ find the best matching one (for a pair of images) */
/*CV_WRAP*/ /*CV_WRAP*/
void match( const Mat& queryDescriptors, const Mat& trainDescriptors, std::vector<DMatch>& matches, const Mat& mask = Mat() ) const; void match( const Mat& queryDescriptors, const Mat& trainDescriptors, std::vector<DMatch>& matches, const Mat& mask = Mat() ) const;
/* for every input descriptor,
find the best matching one (from one image to a set) */
/*CV_WRAP*/
void match( const Mat& queryDescriptors, std::vector<DMatch>& matches, const std::vector<Mat>& masks = std::vector<Mat>() );
/* for every input descriptor,
find the best k matching descriptors (for a pair of images) */
/*CV_WRAP*/
void knnMatch( const Mat& queryDescriptors, const Mat& trainDescriptors, std::vector<std::vector<DMatch> >& matches, int k, const Mat& mask = Mat(),
bool compactResult = false ) const;
/* for every input descriptor,
find the best k matching descriptors (from one image to a set) */
/*CV_WRAP*/
void knnMatch( const Mat& queryDescriptors, std::vector<std::vector<DMatch> >& matches, int k, const std::vector<Mat>& masks = std::vector<Mat>(),
bool compactResult = false );
/* for every input descriptor, find all the ones falling in a
certain matching radius (for a pair of images) */
/*CV_WRAP*/
void radiusMatch( const Mat& queryDescriptors, const Mat& trainDescriptors, std::vector<std::vector<DMatch> >& matches, float maxDistance,
const Mat& mask = Mat(), bool compactResult = false ) const;
/* for every input descriptor, find all the ones falling in a
certain matching radius (from one image to a set) */
/*CV_WRAP*/
void radiusMatch( const Mat& queryDescriptors, std::vector<std::vector<DMatch> >& matches, float maxDistance, const std::vector<Mat>& masks =
std::vector<Mat>(),
bool compactResult = false );
/* store new descriptors to be inserted in dataset */
/*CV_WRAP*/
void add( const std::vector<Mat>& descriptors );
/* store new descriptors into dataset */
/*CV_WRAP*/
void train();
/* constructor with smart pointer */
/*CV_WRAP*/
static Ptr<BinaryDescriptorMatcher> createBinaryDescriptorMatcher();
/* clear dataset and internal data */
/*CV_WRAP*/
void clear();
/* constructor */
/*CV_WRAP*/
BinaryDescriptorMatcher();
/* destructor */
~BinaryDescriptorMatcher()
{
}
/* for every input descriptor, protected:
find the best matching one (from one image to a set) */ /* function inherited from Algorithm */
/*CV_WRAP*/ AlgorithmInfo* info() const;
void match( const Mat& queryDescriptors, std::vector<DMatch>& matches, const std::vector<Mat>& masks = std::vector<Mat>() );
/* for every input descriptor, private:
find the best k matching descriptors (for a pair of images) */ class Array32
/*CV_WRAP*/ {
void knnMatch( const Mat& queryDescriptors, const Mat& trainDescriptors, std::vector<std::vector<DMatch> >& matches, int k, const Mat& mask = Mat(),
bool compactResult = false ) const;
/* for every input descriptor, public:
find the best k matching descriptors (from one image to a set) */ /* set ARRAY_RESIZE_FACTOR */
/*CV_WRAP*/ static void setArrayResizeFactor( double arf );
void knnMatch( const Mat& queryDescriptors, std::vector<std::vector<DMatch> >& matches, int k, const std::vector<Mat>& masks = std::vector<Mat>(),
bool compactResult = false );
/* for every input descriptor, find all the ones falling in a /* constructor */
certain matching radius (for a pair of images) */ Array32();
/*CV_WRAP*/
void radiusMatch( const Mat& queryDescriptors, const Mat& trainDescriptors, std::vector<std::vector<DMatch> >& matches, float maxDistance,
const Mat& mask = Mat(), bool compactResult = false ) const;
/* for every input descriptor, find all the ones falling in a /* destructor */
certain matching radius (from one image to a set) */ ~Array32();
/*CV_WRAP*/
void radiusMatch( const Mat& queryDescriptors, std::vector<std::vector<DMatch> >& matches, float maxDistance, const std::vector<Mat>& masks =
std::vector<Mat>(),
bool compactResult = false );
/* store new descriptors to be inserted in dataset */ /* cleaning function used in destructor */
/*CV_WRAP*/ void cleanup();
void add( const std::vector<Mat>& descriptors );
/* store new descriptors into dataset */ /* push data */
/*CV_WRAP*/ void push( UINT32 data );
void train();
/* constructor with smart pointer */ /* insert data at given index */
/*CV_WRAP*/ void insert( UINT32 index, UINT32 data );
static Ptr<BinaryDescriptorMatcher> createBinaryDescriptorMatcher();
/* clear dataset and internal data */ /* return data */
/*CV_WRAP*/ UINT32* data();
void clear();
/* constructor */ /* return data size */
/*CV_WRAP*/ UINT32 size();
BinaryDescriptorMatcher();
/* destructor */ /* return capacity */
~BinaryDescriptorMatcher() UINT32 capacity();
{
}
protected: /* definition of operator = */
/* function inherited from Algorithm */ void operator=( const Array32& );
AlgorithmInfo* info() const;
private: /* print data */
/* retrieve Hamming distances */ void print();
void checkKDistances( UINT32 * numres, int k, std::vector<int>& k_distances, int row, int string_length ) const;
/* initializer */
void init( int size );
/* matrix to store new descriptors */ /* data */
Mat descriptorsMat; UINT32 *arr;
/* map storing where each bunch of descriptors benins in DS */ };
std::map<int, int> indexesMap;
class BucketGroup
{
/* internal MiHaser representing dataset */ public:
Mihasher* dataset; /* constructor */
BucketGroup();
/* index from which next added descriptors' bunch must begin */ /* destructor */
int nextAddedIndex; ~BucketGroup();
/* number of images whose descriptors are stored in DS */ /* insert data into the bucket */
int numImages; void insert( int subindex, UINT32 data );
/* number of descriptors in dataset */ /* perform a query to the bucket */
int descrInDS; UINT32* query( int subindex, int *size );
/* data fields */
UINT32 empty;
Array32 *group;
};
class SparseHashtable
{
private:
/* Maximum bits per key before folding the table */
static const int MAX_B;
/* Bins (each bin is an Array object for duplicates of the same key) */
BucketGroup *table;
public:
/* constructor */
SparseHashtable();
/* destructor */
~SparseHashtable();
/* initializer */
int init( int _b );
/* insert data */
void insert( UINT64 index, UINT32 data );
/* query data */
UINT32* query( UINT64 index, int* size );
/* Bits per index */
int b;
/* Number of bins */
UINT64 size;
};
/* class defining a sequence of bits */
class bitarray
{
public:
/* pointer to bits sequence and sequence's length */
UINT32 *arr;
UINT32 length;
/* constructor setting default values */
bitarray()
{
arr = NULL;
length = 0;
}
/* constructor setting sequence's length */
bitarray( UINT64 _bits )
{
init( _bits );
}
/* initializer of private fields */
void init( UINT64 _bits )
{
length = (UINT32) ceil( _bits / 32.00 );
arr = new UINT32[length];
erase();
}
/* destructor */
~bitarray()
{
if( arr )
delete[] arr;
}
inline void flip( UINT64 index )
{
arr[index >> 5] ^= ( (UINT32) 0x01 ) << ( index % 32 );
}
inline void set( UINT64 index )
{
arr[index >> 5] |= ( (UINT32) 0x01 ) << ( index % 32 );
}
inline UINT8 get( UINT64 index )
{
return ( arr[index >> 5] & ( ( (UINT32) 0x01 ) << ( index % 32 ) ) ) != 0;
}
/* reserve menory for an UINT32 */
inline void erase()
{
memset( arr, 0, sizeof(UINT32) * length );
}
};
class Mihasher
{
public:
/* Bits per code */
int B;
/* B/8 */
int B_over_8;
/* Bits per chunk (must be less than 64) */
int b;
/* Number of chunks */
int m;
/* Number of chunks with b bits (have 1 bit more than others) */
int mplus;
/* Maximum hamming search radius (we use B/2 by default) */
int D;
/* Maximum hamming search radius per substring */
int d;
/* Maximum results to return */
int K;
/* Number of codes */
UINT64 N;
/* Table of original full-length codes */
cv::Mat codes;
/* Counter for eliminating duplicate results (it is not thread safe) */
bitarray *counter;
/* Array of m hashtables */
SparseHashtable *H;
/* Volume of a b-bit Hamming ball with radius s (for s = 0 to d) */
UINT32 *xornum;
/* Used within generation of binary codes at a certain Hamming distance */
int power[100];
/* constructor */
Mihasher();
/* desctructor */
~Mihasher();
/* constructor 2 */
Mihasher( int B, int m );
/* K setter */
void setK( int K );
/* populate tables */
void populate( cv::Mat & codes, UINT32 N, int dim1codes );
/* execute a batch query */
void batchquery( UINT32 * results, UINT32 *numres/*, qstat *stats*/, const cv::Mat & q, UINT32 numq, int dim1queries );
private:
/* execute a single query */
void query( UINT32 * results, UINT32* numres/*, qstat *stats*/, UINT8 *q, UINT64 * chunks, UINT32 * res );
};
/* retrieve Hamming distances */
void checkKDistances( UINT32 * numres, int k, std::vector<int>& k_distances, int row, int string_length ) const;
/* matrix to store new descriptors */
Mat descriptorsMat;
/* map storing where each bunch of descriptors benins in DS */
std::map<int, int> indexesMap;
/* internal MiHaser representing dataset */
Mihasher* dataset;
/* index from which next added descriptors' bunch must begin */
int nextAddedIndex;
/* number of images whose descriptors are stored in DS */
int numImages;
/* number of descriptors in dataset */
int descrInDS;
}; };
...@@ -408,17 +1151,17 @@ class CV_EXPORTS BinaryDescriptorMatcher : public Algorithm ...@@ -408,17 +1151,17 @@ class CV_EXPORTS BinaryDescriptorMatcher : public Algorithm
/* struct for drawing options */ /* struct for drawing options */
struct CV_EXPORTS DrawLinesMatchesFlags struct CV_EXPORTS DrawLinesMatchesFlags
{ {
enum enum
{ {
DEFAULT = 0, // Output image matrix will be created (Mat::create), DEFAULT = 0, // Output image matrix will be created (Mat::create),
// i.e. existing memory of output image may be reused. // i.e. existing memory of output image may be reused.
// Two source images, matches, and single keylines // Two source images, matches, and single keylines
// will be drawn. // will be drawn.
DRAW_OVER_OUTIMG = 1, // Output image matrix will not be DRAW_OVER_OUTIMG = 1,// Output image matrix will not be
// created (using Mat::create). Matches will be drawn // created (using Mat::create). Matches will be drawn
// on existing content of output image. // on existing content of output image.
NOT_DRAW_SINGLE_LINES = 2 // Single keylines will not be drawn. NOT_DRAW_SINGLE_LINES = 2// Single keylines will not be drawn.
}; };
}; };
/* draw matches between two images */ /* draw matches between two images */
......
modules/line_descriptor/include/opencv2/line_descriptor/ed_line_detector.hpp deleted 100644 → 0
View file @ e39065bf
/*IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
By downloading, copying, installing or using the software you agree to this license.
If you do not agree to this license, do not download, install,
copy or use the software.
License Agreement
For Open Source Computer Vision Library
Copyright (C) 2011-2012, Lilian Zhang, all rights reserved.
Copyright (C) 2013, Manuele Tamburrano, Stefano Fabri, all rights reserved.
Copyright (C) 2014, Biagio Montesano, all rights reserved.
Third party copyrights are property of their respective owners.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright notice,
this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
* The name of the copyright holders may not be used to endorse or promote products
derived from this software without specific prior written permission.
This software is provided by the copyright holders and contributors "as is" and
any express or implied warranties, including, but not limited to, the implied
warranties of merchantability and fitness for a particular purpose are disclaimed.
In no event shall the Intel Corporation or contributors be liable for any direct,
indirect, incidental, special, exemplary, or consequential damages
(including, but not limited to, procurement of substitute goods or services;
loss of use, data, or profits; or business interruption) however caused
and on any theory of liability, whether in contract, strict liability,
or tort (including negligence or otherwise) arising in any way out of
the use of this software, even if advised of the possibility of such damage.
*/
#ifndef __OPENCV_ED_LINE_DETECTOR_HH_
#define __OPENCV_ED_LINE_DETECTOR_HH_
#ifdef _WIN32
#pragma warning( disable : 4267 )
#endif
#include <vector>
#include <iostream>
#include <list>
#include <opencv2/core.hpp>
#include <opencv2/features2d.hpp>
#include <opencv2/imgproc.hpp>
#include <opencv2/highgui.hpp>
#include <opencv2/core/utility.hpp>
struct Pixel
{
unsigned int x; //X coordinate
unsigned int y; //Y coordinate
};
struct EdgeChains
{
std::vector<unsigned int> xCors; //all the x coordinates of edge points
std::vector<unsigned int> yCors; //all the y coordinates of edge points
std::vector<unsigned int> sId; //the start index of each edge in the coordinate arrays
unsigned int numOfEdges; //the number of edges whose length are larger than minLineLen; numOfEdges < sId.size;
};
struct LineChains
{
std::vector<unsigned int> xCors; //all the x coordinates of line points
std::vector<unsigned int> yCors; //all the y coordinates of line points
std::vector<unsigned int> sId; //the start index of each line in the coordinate arrays
unsigned int numOfLines; //the number of lines whose length are larger than minLineLen; numOfLines < sId.size;
};
typedef std::list<Pixel> PixelChain; //each edge is a pixel chain
struct EDLineParam
{
int ksize;
float sigma;
float gradientThreshold;
float anchorThreshold;
int scanIntervals;
int minLineLen;
double lineFitErrThreshold;
};
#define RELATIVE_ERROR_FACTOR 100.0
#define MLN10 2.30258509299404568402
#define log_gamma(x) ((x)>15.0?log_gamma_windschitl(x):log_gamma_lanczos(x))
/* This class is used to detect lines from input image.
* First, edges are extracted from input image following the method presented in Cihan Topal and
* Cuneyt Akinlar's paper:"Edge Drawing: A Heuristic Approach to Robust Real-Time Edge Detection", 2010.
* Then, lines are extracted from the edge image following the method presented in Cuneyt Akinlar and
* Cihan Topal's paper:"EDLines: A real-time line segment detector with a false detection control", 2011
* PS: The linking step of edge detection has a little bit difference with the Edge drawing algorithm
* described in the paper. The edge chain doesn't stop when the pixel direction is changed.
*/
class EDLineDetector
{
public:
EDLineDetector();
EDLineDetector( EDLineParam param );
~EDLineDetector();
/*extract edges from image
*image: In, gray image;
*edges: Out, store the edges, each edge is a pixel chain
*return -1: error happen
*/
int EdgeDrawing( cv::Mat &image, EdgeChains &edgeChains );
/*extract lines from image
*image: In, gray image;
*lines: Out, store the extracted lines,
*return -1: error happen
*/
int EDline( cv::Mat &image, LineChains &lines );
/* extract line from image, and store them */
int EDline( cv::Mat &image );
cv::Mat dxImg_; //store the dxImg;
cv::Mat dyImg_; //store the dyImg;
cv::Mat gImgWO_; //store the gradient image without threshold;
LineChains lines_; //store the detected line chains;
//store the line Equation coefficients, vec3=[w1,w2,w3] for line w1*x + w2*y + w3=0;
std::vector<std::vector<double> > lineEquations_;
//store the line endpoints, [x1,y1,x2,y3]
std::vector<std::vector<float> > lineEndpoints_;
//store the line direction
std::vector<float> lineDirection_;
//store the line salience, which is the summation of gradients of pixels on line
std::vector<float> lineSalience_;
// image sizes
unsigned int imageWidth;
unsigned int imageHeight;
/*The threshold of line fit error;
*If lineFitErr is large than this threshold, then
*the pixel chain is not accepted as a single line segment.*/
double lineFitErrThreshold_;
/*the threshold of pixel gradient magnitude.
*Only those pixel whose gradient magnitude are larger than this threshold will be
*taken as possible edge points. Default value is 36*/
short gradienThreshold_;
/*If the pixel's gradient value is bigger than both of its neighbors by a
*certain threshold (ANCHOR_THRESHOLD), the pixel is marked to be an anchor.
*Default value is 8*/
unsigned char anchorThreshold_;
/*anchor testing can be performed at different scan intervals, i.e.,
*every row/column, every second row/column etc.
*Default value is 2*/
unsigned int scanIntervals_;
int minLineLen_; //minimal acceptable line length
private:
void InitEDLine_();
/*For an input edge chain, find the best fit line, the default chain length is minLineLen_
*xCors: In, pointer to the X coordinates of pixel chain;
*yCors: In, pointer to the Y coordinates of pixel chain;
*offsetS:In, start index of this chain in vector;
*lineEquation: Out, [a,b] which are the coefficient of lines y=ax+b(horizontal) or x=ay+b(vertical);
*return: line fit error; -1:error happens;
*/
double LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, std::vector<double> &lineEquation );
/*For an input pixel chain, find the best fit line. Only do the update based on new points.
*For A*x=v, Least square estimation of x = Inv(A^T * A) * (A^T * v);
*If some new observations are added, i.e, [A; A'] * x = [v; v'],
*then x' = Inv(A^T * A + (A')^T * A') * (A^T * v + (A')^T * v');
*xCors: In, pointer to the X coordinates of pixel chain;
*yCors: In, pointer to the Y coordinates of pixel chain;
*offsetS:In, start index of this chain in vector;
*newOffsetS: In, start index of extended part;
*offsetE:In, end index of this chain in vector;
*lineEquation: Out, [a,b] which are the coefficient of lines y=ax+b(horizontal) or x=ay+b(vertical);
*return: line fit error; -1:error happens;
*/
double LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int newOffsetS, unsigned int offsetE,
std::vector<double> &lineEquation );
/* Validate line based on the Helmholtz principle, which basically states that
* for a structure to be perceptually meaningful, the expectation of this structure
* by chance must be very low.
*/
bool LineValidation_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int offsetE, std::vector<double> &lineEquation,
float &direction );
bool bValidate_; //flag to decide whether line will be validated
int ksize_; //the size of Gaussian kernel: ksize X ksize, default value is 5.
float sigma_; //the sigma of Gaussian kernal, default value is 1.0.
/*For example, there two edges in the image:
*edge1 = [(7,4), (8,5), (9,6),| (10,7)|, (11, 8), (12,9)] and
*edge2 = [(14,9), (15,10), (16,11), (17,12),| (18, 13)|, (19,14)] ; then we store them as following:
*pFirstPartEdgeX_ = [10, 11, 12, 18, 19];//store the first part of each edge[from middle to end]
*pFirstPartEdgeY_ = [7, 8, 9, 13, 14];
*pFirstPartEdgeS_ = [0,3,5];// the index of start point of first part of each edge
*pSecondPartEdgeX_ = [10, 9, 8, 7, 18, 17, 16, 15, 14];//store the second part of each edge[from middle to front]
*pSecondPartEdgeY_ = [7, 6, 5, 4, 13, 12, 11, 10, 9];//anchor points(10, 7) and (18, 13) are stored again
*pSecondPartEdgeS_ = [0, 4, 9];// the index of start point of second part of each edge
*This type of storage order is because of the order of edge detection process.
*For each edge, start from one anchor point, first go right, then go left or first go down, then go up*/
//store the X coordinates of the first part of the pixels for chains
unsigned int *pFirstPartEdgeX_;
//store the Y coordinates of the first part of the pixels for chains
unsigned int *pFirstPartEdgeY_;
//store the start index of every edge chain in the first part arrays
unsigned int *pFirstPartEdgeS_;
//store the X coordinates of the second part of the pixels for chains
unsigned int *pSecondPartEdgeX_;
//store the Y coordinates of the second part of the pixels for chains
unsigned int *pSecondPartEdgeY_;
//store the start index of every edge chain in the second part arrays
unsigned int *pSecondPartEdgeS_;
//store the X coordinates of anchors
unsigned int *pAnchorX_;
//store the Y coordinates of anchors
unsigned int *pAnchorY_;
//edges
cv::Mat edgeImage_;
cv::Mat gImg_; //store the gradient image;
cv::Mat dirImg_; //store the direction image
double logNT_;
cv::Mat_<float> ATA; //the previous matrix of A^T * A;
cv::Mat_<float> ATV; //the previous vector of A^T * V;
cv::Mat_<float> fitMatT; //the matrix used in line fit function;
cv::Mat_<float> fitVec; //the vector used in line fit function;
cv::Mat_<float> tempMatLineFit; //the matrix used in line fit function;
cv::Mat_<float> tempVecLineFit; //the vector used in line fit function;
/** Compare doubles by relative error.
The resulting rounding error after floating point computations
depend on the specific operations done. The same number computed by
different algorithms could present different rounding errors. For a
useful comparison, an estimation of the relative rounding error
should be considered and compared to a factor times EPS. The factor
should be related to the accumulated rounding error in the chain of
computation. Here, as a simplification, a fixed factor is used.
*/
static int double_equal( double a, double b )
{
double abs_diff, aa, bb, abs_max;
/* trivial case */
if( a == b )
return true;
abs_diff = fabs( a - b );
aa = fabs( a );
bb = fabs( b );
abs_max = aa > bb ? aa : bb;
/* DBL_MIN is the smallest normalized number, thus, the smallest
number whose relative error is bounded by DBL_EPSILON. For
smaller numbers, the same quantization steps as for DBL_MIN
are used. Then, for smaller numbers, a meaningful "relative"
error should be computed by dividing the difference by DBL_MIN. */
if( abs_max < DBL_MIN )
abs_max = DBL_MIN;
/* equal if relative error <= factor x eps */
return ( abs_diff / abs_max ) <= ( RELATIVE_ERROR_FACTOR * DBL_EPSILON );
}
/** Computes the natural logarithm of the absolute value of
the gamma function of x using the Lanczos approximation.
See http://www.rskey.org/gamma.htm
The formula used is
@f[
\Gamma(x) = \frac{ \sum_{n=0}^{N} q_n x^n }{ \Pi_{n=0}^{N} (x+n) }
(x+5.5)^{x+0.5} e^{-(x+5.5)}
@f]
so
@f[
\log\Gamma(x) = \log\left( \sum_{n=0}^{N} q_n x^n \right)
+ (x+0.5) \log(x+5.5) - (x+5.5) - \sum_{n=0}^{N} \log(x+n)
@f]
and
q0 = 75122.6331530,
q1 = 80916.6278952,
q2 = 36308.2951477,
q3 = 8687.24529705,
q4 = 1168.92649479,
q5 = 83.8676043424,
q6 = 2.50662827511.
*/
static double log_gamma_lanczos( double x )
{
static double q[7] =
{ 75122.6331530, 80916.6278952, 36308.2951477, 8687.24529705, 1168.92649479, 83.8676043424, 2.50662827511 };
double a = ( x + 0.5 ) * log( x + 5.5 ) - ( x + 5.5 );
double b = 0.0;
int n;
for ( n = 0; n < 7; n++ )
{
a -= log( x + (double) n );
b += q[n] * pow( x, (double) n );
}
return a + log( b );
}
/** Computes the natural logarithm of the absolute value of
the gamma function of x using Windschitl method.
See http://www.rskey.org/gamma.htm
The formula used is
@f[
\Gamma(x) = \sqrt{\frac{2\pi}{x}} \left( \frac{x}{e}
\sqrt{ x\sinh(1/x) + \frac{1}{810x^6} } \right)^x
@f]
so
@f[
\log\Gamma(x) = 0.5\log(2\pi) + (x-0.5)\log(x) - x
+ 0.5x\log\left( x\sinh(1/x) + \frac{1}{810x^6} \right).
@f]
This formula is a good approximation when x > 15.
*/
static double log_gamma_windschitl( double x )
{
return 0.918938533204673 + ( x - 0.5 ) * log( x ) - x + 0.5 * x * log( x * sinh( 1 / x ) + 1 / ( 810.0 * pow( x, 6.0 ) ) );
}
/** Computes -log10(NFA).
NFA stands for Number of False Alarms:
@f[
\mathrm{NFA} = NT \cdot B(n,k,p)
@f]
- NT - number of tests
- B(n,k,p) - tail of binomial distribution with parameters n,k and p:
@f[
B(n,k,p) = \sum_{j=k}^n
\left(\begin{array}{c}n\\j\end{array}\right)
p^{j} (1-p)^{n-j}
@f]
The value -log10(NFA) is equivalent but more intuitive than NFA:
- -1 corresponds to 10 mean false alarms
- 0 corresponds to 1 mean false alarm
- 1 corresponds to 0.1 mean false alarms
- 2 corresponds to 0.01 mean false alarms
- ...
Used this way, the bigger the value, better the detection,
and a logarithmic scale is used.
@param n,k,p binomial parameters.
@param logNT logarithm of Number of Tests
The computation is based in the gamma function by the following
relation:
@f[
\left(\begin{array}{c}n\\k\end{array}\right)
= \frac{ \Gamma(n+1) }{ \Gamma(k+1) \cdot \Gamma(n-k+1) }.
@f]
We use efficient algorithms to compute the logarithm of
the gamma function.
To make the computation faster, not all the sum is computed, part
of the terms are neglected based on a bound to the error obtained
(an error of 10% in the result is accepted).
*/
static double nfa( int n, int k, double p, double logNT )
{
double tolerance = 0.1; /* an error of 10% in the result is accepted */
double log1term, term, bin_term, mult_term, bin_tail, err, p_term;
int i;
/* check parameters */
if( n < 0 || k < 0 || k > n || p <= 0.0 || p >= 1.0 )
{
std::cout << "nfa: wrong n, k or p values." << std::endl;
exit( 0 );
}
/* trivial cases */
if( n == 0 || k == 0 )
return -logNT;
if( n == k )
return -logNT - (double) n * log10( p );
/* probability term */
p_term = p / ( 1.0 - p );
/* compute the first term of the series */
/*
binomial_tail(n,k,p) = sum_{i=k}^n bincoef(n,i) * p^i * (1-p)^{n-i}
where bincoef(n,i) are the binomial coefficients.
But
bincoef(n,k) = gamma(n+1) / ( gamma(k+1) * gamma(n-k+1) ).
We use this to compute the first term. Actually the log of it.
*/
log1term = log_gamma( (double) n + 1.0 ) - log_gamma( (double ) k + 1.0 ) - log_gamma( (double ) ( n - k ) + 1.0 ) + (double) k * log( p )
+ (double) ( n - k ) * log( 1.0 - p );
term = exp( log1term );
/* in some cases no more computations are needed */
if( double_equal( term, 0.0 ) )
{ /* the first term is almost zero */
if( (double) k > (double) n * p ) /* at begin or end of the tail? */
return -log1term / MLN10 - logNT; /* end: use just the first term */
else
return -logNT; /* begin: the tail is roughly 1 */
}
/* compute more terms if needed */
bin_tail = term;
for ( i = k + 1; i <= n; i++ )
{
/* As
term_i = bincoef(n,i) * p^i * (1-p)^(n-i)
and
bincoef(n,i)/bincoef(n,i-1) = n-i+1 / i,
then,
term_i / term_i-1 = (n-i+1)/i * p/(1-p)
and
term_i = term_i-1 * (n-i+1)/i * p/(1-p).
p/(1-p) is computed only once and stored in 'p_term'.
*/
bin_term = (double) ( n - i + 1 ) / (double) i;
mult_term = bin_term * p_term;
term *= mult_term;
bin_tail += term;
if( bin_term < 1.0 )
{
/* When bin_term<1 then mult_term_j<mult_term_i for j>i.
Then, the error on the binomial tail when truncated at
the i term can be bounded by a geometric series of form
term_i * sum mult_term_i^j. */
err = term * ( ( 1.0 - pow( mult_term, (double) ( n - i + 1 ) ) ) / ( 1.0 - mult_term ) - 1.0 );
/* One wants an error at most of tolerance*final_result, or:
tolerance * abs(-log10(bin_tail)-logNT).
Now, the error that can be accepted on bin_tail is
given by tolerance*final_result divided by the derivative
of -log10(x) when x=bin_tail. that is:
tolerance * abs(-log10(bin_tail)-logNT) / (1/bin_tail)
Finally, we truncate the tail if the error is less than:
tolerance * abs(-log10(bin_tail)-logNT) * bin_tail */
if( err < tolerance * fabs( -log10( bin_tail ) - logNT ) * bin_tail )
break;
}
}
return -log10( bin_tail ) - logNT;
}
};
#endif /* EDLINEDETECTOR_HH_ */
modules/line_descriptor/include/opencv2/line_descriptor/line_structure.hpp deleted 100644 → 0
View file @ e39065bf
/*IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
By downloading, copying, installing or using the software you agree to this license.
If you do not agree to this license, do not download, install,
copy or use the software.
License Agreement
For Open Source Computer Vision Library
Copyright (C) 2011-2012, Lilian Zhang, all rights reserved.
Copyright (C) 2013, Manuele Tamburrano, Stefano Fabri, all rights reserved.
Copyright (C) 2014, Biagio Montesano, all rights reserved.
Third party copyrights are property of their respective owners.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright notice,
this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
* The name of the copyright holders may not be used to endorse or promote products
derived from this software without specific prior written permission.
This software is provided by the copyright holders and contributors "as is" and
any express or implied warranties, including, but not limited to, the implied
warranties of merchantability and fitness for a particular purpose are disclaimed.
In no event shall the Intel Corporation or contributors be liable for any direct,
indirect, incidental, special, exemplary, or consequential damages
(including, but not limited to, procurement of substitute goods or services;
loss of use, data, or profits; or business interruption) however caused
and on any theory of liability, whether in contract, strict liability,
or tort (including negligence or otherwise) arising in any way out of
the use of this software, even if advised of the possibility of such damage.
*/
#ifndef __OPENCV_LINE_STRUCTURE_HH_
#define __OPENCV_LINE_STRUCTURE_HH_
#include <vector>
/* struct to represent lines extracted from an octave */
struct OctaveLine
{
unsigned int octaveCount; //the octave which this line is detected
unsigned int lineIDInOctave; //the line ID in that octave image
unsigned int lineIDInScaleLineVec; //the line ID in Scale line vector
float lineLength; //the length of line in original image scale
};
// A 2D line (normal equation parameters).
struct SingleLine
{
//note: rho and theta are based on coordinate origin, i.e. the top-left corner of image
double rho; //unit: pixel length
double theta; //unit: rad
double linePointX; // = rho * cos(theta);
double linePointY; // = rho * sin(theta);
//for EndPoints, the coordinate origin is the top-left corner of image.
double startPointX;
double startPointY;
double endPointX;
double endPointY;
//direction of a line, the angle between positive line direction (dark side is in the left) and positive X axis.
double direction;
//mean gradient magnitude
double gradientMagnitude;
//mean gray value of pixels in dark side of line
double darkSideGrayValue;
//mean gray value of pixels in light side of line
double lightSideGrayValue;
//the length of line
double lineLength;
//the width of line;
double width;
//number of pixels
int numOfPixels;
//the decriptor of line
std::vector<double> descriptor;
};
// Specifies a vector of lines.
typedef std::vector<SingleLine> Lines_list;
struct OctaveSingleLine
{
/*endPoints, the coordinate origin is the top-left corner of the original image.
*startPointX = sPointInOctaveX * (factor)^octaveCount; */
float startPointX;
float startPointY;
float endPointX;
float endPointY;
//endPoints, the coordinate origin is the top-left corner of the octave image.
float sPointInOctaveX;
float sPointInOctaveY;
float ePointInOctaveX;
float ePointInOctaveY;
//direction of a line, the angle between positive line direction (dark side is in the left) and positive X axis.
float direction;
//the summation of gradient magnitudes of pixels on lines
float salience;
//the length of line
float lineLength;
//number of pixels
unsigned int numOfPixels;
//the octave which this line is detected
unsigned int octaveCount;
//the decriptor of line
std::vector<float> descriptor;
};
// Specifies a vector of lines.
typedef std::vector<OctaveSingleLine> LinesVec;
// each element in ScaleLines is a vector of lines
// which corresponds the same line detected in different octave images.
typedef std::vector<LinesVec> ScaleLines;
#endif
modules/line_descriptor/include/opencv2/line_descriptor/mihasher.hpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __OPENCV_MIHASHER_HPP
#define __OPENCV_MIHASHER_HPP
#ifdef _WIN32
#pragma warning( disable : 4267 )
#endif
#include "types.hpp"
#include "bitops.hpp"
#include "sparse_hashtable.hpp"
#include "bitarray.hpp"
#include <opencv2/core.hpp>
#include <opencv2/imgproc.hpp>
class Mihasher
{
private:
/* Bits per code */
int B;
/* B/8 */
int B_over_8;
/* Bits per chunk (must be less than 64) */
int b;
/* Number of chunks */
int m;
/* Number of chunks with b bits (have 1 bit more than others) */
int mplus;
/* Maximum hamming search radius (we use B/2 by default) */
int D;
/* Maximum hamming search radius per substring */
int d;
/* Maximum results to return */
int K;
/* Number of codes */
UINT64 N;
/* Table of original full-length codes */
cv::Mat codes;
/* Counter for eliminating duplicate results (it is not thread safe) */
bitarray *counter;
/* Array of m hashtables */
SparseHashtable *H;
/* Volume of a b-bit Hamming ball with radius s (for s = 0 to d) */
UINT32 *xornum;
/* Used within generation of binary codes at a certain Hamming distance */
int power[100];
public:
/* constructor */
Mihasher();
/* desctructor */
~Mihasher();
/* constructor 2 */
Mihasher( int B, int m );
/* K setter */
void setK( int K );
/* populate tables */
void populate( cv::Mat & codes, UINT32 N, int dim1codes );
/* execute a batch query */
void batchquery( UINT32 * results, UINT32 *numres/*, qstat *stats*/, const cv::Mat & q, UINT32 numq, int dim1queries );
private:
/* execute a single query */
void query( UINT32 * results, UINT32* numres/*, qstat *stats*/, UINT8 *q, UINT64 * chunks, UINT32 * res );
};
#endif
modules/line_descriptor/include/opencv2/line_descriptor/sparse_hashtable.hpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __OPENCV_SPARSE_HASHTABLE_HPP
#define __OPENCV_SPARSE_HASHTABLE_HPP
#ifdef _WIN32
#pragma warning( disable : 4267 )
#endif
#include "types.hpp"
#include "bucket_group.hpp"
class SparseHashtable
{
private:
/* Maximum bits per key before folding the table */
static const int MAX_B;
/* Bins (each bin is an Array object for duplicates of the same key) */
BucketGroup *table;
public:
/* constructor */
SparseHashtable();
/* destructor */
~SparseHashtable();
/* initializer */
int init( int _b );
/* insert data */
void insert( UINT64 index, UINT32 data );
/* query data */
UINT32* query( UINT64 index, int* size );
/* Bits per index */
int b;
/* Number of bins */
UINT64 size;
};
#endif
modules/line_descriptor/samples/radius_matching.cpp
View file @ 9a52f462
...@@ -62,22 +62,22 @@ static const char* keys = ...@@ -62,22 +62,22 @@ static const char* keys =
static void help() static void help()
{ {
std::cout << "\nThis example shows the functionalities of radius matching " << "Please, run this sample using a command in the form\n" std::cout << "\nThis example shows the functionalities of radius matching " << "Please, run this sample using a command in the form\n"
<< "./example_line_descriptor_radius_matching <path_to_input_images>/" << std::endl; << "./example_line_descriptor_radius_matching <path_to_input_images>/" << std::endl;
} }
int main( int argc, char** argv ) int main( int argc, char** argv )
{ {
/* get parameters from comand line */ /* get parameters from comand line */
CommandLineParser parser( argc, argv, keys ); CommandLineParser parser( argc, argv, keys );
String pathToImages = parser.get<String>( 0 ); String pathToImages = parser.get < String > ( 0 );
/* create structures for hosting KeyLines and descriptors */ /* create structures for hosting KeyLines and descriptors */
int num_elements = sizeof ( images ) / sizeof ( images[0] ); int num_elements = sizeof ( images ) / sizeof ( images[0] );
std::vector<Mat> descriptorsMat; std::vector < Mat > descriptorsMat;
std::vector<std::vector<KeyLine> > linesMat; std::vector < std::vector<KeyLine> > linesMat;
/*create a pointer to a BinaryDescriptor object */ /*create a pointer to a BinaryDescriptor object */
Ptr<BinaryDescriptor> bd = BinaryDescriptor::createBinaryDescriptor(); Ptr < BinaryDescriptor > bd = BinaryDescriptor::createBinaryDescriptor();
/* compute lines and descriptors */ /* compute lines and descriptors */
for ( int i = 0; i < num_elements; i++ ) for ( int i = 0; i < num_elements; i++ )
...@@ -97,7 +97,7 @@ int main( int argc, char** argv ) ...@@ -97,7 +97,7 @@ int main( int argc, char** argv )
} }
/* compute lines and descriptors */ /* compute lines and descriptors */
std::vector<KeyLine> lines; std::vector < KeyLine > lines;
Mat computedDescr; Mat computedDescr;
bd->detect( loadedImage, lines ); bd->detect( loadedImage, lines );
bd->compute( loadedImage, lines, computedDescr ); bd->compute( loadedImage, lines, computedDescr );
...@@ -121,19 +121,19 @@ int main( int argc, char** argv ) ...@@ -121,19 +121,19 @@ int main( int argc, char** argv )
std::cout << "It has been generated a matrix of " << queries.rows << " descriptors" << std::endl; std::cout << "It has been generated a matrix of " << queries.rows << " descriptors" << std::endl;
/* create a BinaryDescriptorMatcher object */ /* create a BinaryDescriptorMatcher object */
Ptr<BinaryDescriptorMatcher> bdm = BinaryDescriptorMatcher::createBinaryDescriptorMatcher(); Ptr < BinaryDescriptorMatcher > bdm = BinaryDescriptorMatcher::createBinaryDescriptorMatcher();
/* populate matcher */ /* populate matcher */
bdm->add( descriptorsMat ); bdm->add( descriptorsMat );
/* compute matches */ /* compute matches */
std::vector<std::vector<DMatch> > matches; std::vector < std::vector<DMatch> > matches;
bdm->radiusMatch( queries, matches, 30 ); bdm->radiusMatch( queries, matches, 30 );
std::cout << "size matches sample " << matches.size() << std::endl; std::cout << "size matches sample " << matches.size() << std::endl;
for ( int i = 0; i < matches.size(); i++ ) for ( int i = 0; i < (int) matches.size(); i++ )
{ {
for ( int j = 0; j < matches[i].size(); j++ ) for ( int j = 0; j < (int) matches[i].size(); j++ )
{ {
std::cout << "match: " << matches[i][j].queryIdx << " " << matches[i][j].trainIdx << " " << matches[i][j].distance << std::endl; std::cout << "match: " << matches[i][j].queryIdx << " " << matches[i][j].trainIdx << " " << matches[i][j].distance << std::endl;
} }
......
modules/line_descriptor/src/array32.cpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#include "precomp.hpp"
/* no need for the static keyword in the definition */
double Array32::ARRAY_RESIZE_FACTOR = 1.1; // minimum is 1.0
double Array32::ARRAY_RESIZE_ADD_FACTOR = 4; // minimum is 1
/* set ARRAY_RESIZE_FACTOR */
void Array32::setArrayResizeFactor( double arf )
{
ARRAY_RESIZE_FACTOR = arf;
}
/* constructor */
Array32::Array32()
{
arr = NULL;
}
/* definition of operator =
Array32& Array32::operator = (const Array32 &rhs) */
void Array32::operator =( const Array32 &rhs )
{
if( &rhs != this )
this->arr = rhs.arr;
}
/* destructor */
Array32::~Array32()
{
cleanup();
}
/* cleaning function used in destructor */
void Array32::cleanup()
{
free( arr );
}
/* push data */
void Array32::push( UINT32 Data )
{
if( arr )
{
if( arr[0] == arr[1] )
{
arr[1] = (UINT32) std::max( ceil( arr[1] * ARRAY_RESIZE_FACTOR ), arr[1] + ARRAY_RESIZE_ADD_FACTOR );
UINT32* new_Data = static_cast<UINT32*>( realloc( arr, sizeof(UINT32) * ( 2 + arr[1] ) ) );
if( new_Data == NULL )
{
/* could not realloc, but orig still valid */
std::cout << "ALERT!!!! Not enough memory, operation aborted!" << std::endl;
exit( 0 );
}
else
{
arr = new_Data;
}
}
arr[2 + arr[0]] = Data;
arr[0]++;
}
else
{
arr = (UINT32*) malloc( (size_t) ( 2 + ARRAY_RESIZE_ADD_FACTOR ) * sizeof(UINT32) );
arr[0] = 1;
arr[1] = 1;
arr[2] = Data;
}
}
/* insert data at given index */
void Array32::insert( UINT32 index, UINT32 Data )
{
if( arr )
{
if( arr[0] == arr[1] )
{
arr[1] = (UINT32) ceil( arr[0] * 1.1 );
UINT32* new_data = static_cast<UINT32*>( realloc( arr, sizeof(UINT32) * ( 2 + arr[1] ) ) );
if( new_data == NULL )
{
// could not realloc, but orig still valid
std::cout << "ALERT!!!! Not enough memory, operation aborted!" << std::endl;
exit( 0 );
}
else
{
arr = new_data;
}
}
memmove( arr + ( 2 + index ) + 1, arr + ( 2 + index ), ( arr[0] - index ) * sizeof ( *arr ) );
arr[2 + index] = Data;
arr[0]++;
}
else
{
arr = (UINT32*) malloc( 3 * sizeof(UINT32) );
arr[0] = 1;
arr[1] = 1;
arr[2] = Data;
}
}
/* return data */
UINT32* Array32::data()
{
return arr ? arr + 2 : NULL;
}
/* return data size */
UINT32 Array32::size()
{
return arr ? arr[0] : 0;
}
/* return capacity */
UINT32 Array32::capacity()
{
return arr ? arr[1] : 0;
}
/* print data */
void Array32::print()
{
for ( int i = 0; i < (int) size(); i++ )
printf( "%d, ", arr[i + 2] );
printf( "\n" );
}
/* initializer */
void Array32::init( int Size )
{
if( arr == NULL )
{
arr = (UINT32*) malloc( ( 2 + Size ) * sizeof(UINT32) );
arr[0] = 0;
arr[1] = Size;
}
}
modules/line_descriptor/src/binary_descriptor.cpp
View file @ 9a52f462
...@@ -42,6 +42,14 @@ ...@@ -42,6 +42,14 @@
#include "precomp.hpp" #include "precomp.hpp"
#define NUM_OF_BANDS 9 #define NUM_OF_BANDS 9
#define Horizontal 255//if |dx|<|dy|;
#define Vertical 0//if |dy|<=|dx|;
#define UpDir 1
#define RightDir 2
#define DownDir 3
#define LeftDir 4
#define TryTime 6
#define SkipEdgePoint 2
//using namespace cv; //using namespace cv;
namespace cv namespace cv
...@@ -119,7 +127,7 @@ void BinaryDescriptor::setWidthOfBand( int width ) ...@@ -119,7 +127,7 @@ void BinaryDescriptor::setWidthOfBand( int width )
images_sizes.resize( params.numOfOctave_ ); images_sizes.resize( params.numOfOctave_ );
for ( int i = 0; i < params.numOfOctave_; i++ ) for ( int i = 0; i < params.numOfOctave_; i++ )
edLineVec_[i] = Ptr<EDLineDetector>( new EDLineDetector() ); edLineVec_[i] = Ptr < EDLineDetector > ( new EDLineDetector() );
/* prepare a vector to host local weights F_l*/ /* prepare a vector to host local weights F_l*/
gaussCoefL_.resize( params.widthOfBand_ * 3 ); gaussCoefL_.resize( params.widthOfBand_ * 3 );
...@@ -184,12 +192,12 @@ void BinaryDescriptor::Params::write( cv::FileStorage& fs ) const ...@@ -184,12 +192,12 @@ void BinaryDescriptor::Params::write( cv::FileStorage& fs ) const
Ptr<BinaryDescriptor> BinaryDescriptor::createBinaryDescriptor() Ptr<BinaryDescriptor> BinaryDescriptor::createBinaryDescriptor()
{ {
return Ptr<BinaryDescriptor>( new BinaryDescriptor() ); return Ptr < BinaryDescriptor > ( new BinaryDescriptor() );
} }
Ptr<BinaryDescriptor> BinaryDescriptor::createBinaryDescriptor( Params parameters ) Ptr<BinaryDescriptor> BinaryDescriptor::createBinaryDescriptor( Params parameters )
{ {
return Ptr<BinaryDescriptor>( new BinaryDescriptor( parameters ) ); return Ptr < BinaryDescriptor > ( new BinaryDescriptor( parameters ) );
} }
/* construct a BinaryDescrptor object and compute external private parameters */ /* construct a BinaryDescrptor object and compute external private parameters */
...@@ -201,7 +209,7 @@ BinaryDescriptor::BinaryDescriptor( const BinaryDescriptor::Params &parameters ) ...@@ -201,7 +209,7 @@ BinaryDescriptor::BinaryDescriptor( const BinaryDescriptor::Params &parameters )
images_sizes.resize( params.numOfOctave_ ); images_sizes.resize( params.numOfOctave_ );
for ( int i = 0; i < params.numOfOctave_; i++ ) for ( int i = 0; i < params.numOfOctave_; i++ )
edLineVec_[i] = Ptr<EDLineDetector>( new EDLineDetector() ); edLineVec_[i] = Ptr < EDLineDetector > ( new EDLineDetector() );
/* prepare a vector to host local weights F_l*/ /* prepare a vector to host local weights F_l*/
gaussCoefL_.resize( params.widthOfBand_ * 3 ); gaussCoefL_.resize( params.widthOfBand_ * 3 );
...@@ -239,7 +247,7 @@ BinaryDescriptor::BinaryDescriptor( const BinaryDescriptor::Params &parameters ) ...@@ -239,7 +247,7 @@ BinaryDescriptor::BinaryDescriptor( const BinaryDescriptor::Params &parameters )
/* definition of operator () */ /* definition of operator () */
void BinaryDescriptor::operator()( InputArray image, InputArray mask, CV_OUT std::vector<KeyLine>& keylines, OutputArray descriptors, void BinaryDescriptor::operator()( InputArray image, InputArray mask, CV_OUT std::vector<KeyLine>& keylines, OutputArray descriptors,
bool useProvidedKeyLines, bool returnFloatDescr ) const bool useProvidedKeyLines, bool returnFloatDescr ) const
{ {
/* create some matrix objects */ /* create some matrix objects */
...@@ -258,7 +266,7 @@ void BinaryDescriptor::operator()( InputArray image, InputArray mask, CV_OUT std ...@@ -258,7 +266,7 @@ void BinaryDescriptor::operator()( InputArray image, InputArray mask, CV_OUT std
bn->edLineVec_.resize( params.numOfOctave_ ); bn->edLineVec_.resize( params.numOfOctave_ );
for ( int i = 0; i < params.numOfOctave_; i++ ) for ( int i = 0; i < params.numOfOctave_; i++ )
bn->edLineVec_[i] = Ptr<EDLineDetector>( new EDLineDetector() ); bn->edLineVec_[i] = Ptr<EDLineDetector>( new EDLineDetector() );
detectImpl( imageMat, keylines, maskMat ); detectImpl( imageMat, keylines, maskMat );
...@@ -270,10 +278,10 @@ void BinaryDescriptor::operator()( InputArray image, InputArray mask, CV_OUT std ...@@ -270,10 +278,10 @@ void BinaryDescriptor::operator()( InputArray image, InputArray mask, CV_OUT std
//descrMat = descriptors.getMat(); //descrMat = descriptors.getMat();
/* compute descriptors */ /* compute descriptors */
if( !useProvidedKeyLines ) if( !useProvidedKeyLines )
computeImpl( imageMat, keylines, descrMat, returnFloatDescr, true ); computeImpl( imageMat, keylines, descrMat, returnFloatDescr, true );
else else
computeImpl( imageMat, keylines, descrMat, returnFloatDescr, false ); computeImpl( imageMat, keylines, descrMat, returnFloatDescr, false );
descrMat.copyTo( descriptors ); descrMat.copyTo( descriptors );
} }
...@@ -400,10 +408,10 @@ void BinaryDescriptor::detect( const Mat& image, CV_OUT std::vector<KeyLine>& ke ...@@ -400,10 +408,10 @@ void BinaryDescriptor::detect( const Mat& image, CV_OUT std::vector<KeyLine>& ke
} }
if( mask.data != NULL && ( mask.size() != image.size() || mask.type() != CV_8UC1 ) ) if( mask.data != NULL && ( mask.size() != image.size() || mask.type() != CV_8UC1 ) )
throw std::runtime_error( "Mask error while detecting lines: please check its dimensions and that data type is CV_8UC1" ); throw std::runtime_error( "Mask error while detecting lines: please check its dimensions and that data type is CV_8UC1" );
else else
detectImpl( image, keylines, mask ); detectImpl( image, keylines, mask );
} }
/* requires line detection (more than one image) */ /* requires line detection (more than one image) */
...@@ -491,7 +499,7 @@ void BinaryDescriptor::detectImpl( const Mat& imageSrc, std::vector<KeyLine>& ke ...@@ -491,7 +499,7 @@ void BinaryDescriptor::detectImpl( const Mat& imageSrc, std::vector<KeyLine>& ke
for ( size_t keyCounter = 0; keyCounter < keylines.size(); keyCounter++ ) for ( size_t keyCounter = 0; keyCounter < keylines.size(); keyCounter++ )
{ {
KeyLine kl = keylines[keyCounter]; KeyLine kl = keylines[keyCounter];
if( mask.at<uchar>( (int) kl.startPointY, (int) kl.startPointX ) == 0 && mask.at<uchar>( (int) kl.endPointY, (int) kl.endPointX ) == 0 ) if( mask.at < uchar > ( (int) kl.startPointY, (int) kl.startPointX ) == 0 && mask.at < uchar > ( (int) kl.endPointY, (int) kl.endPointX ) == 0 )
keylines.erase( keylines.begin() + keyCounter ); keylines.erase( keylines.begin() + keyCounter );
} }
} }
...@@ -500,7 +508,7 @@ void BinaryDescriptor::detectImpl( const Mat& imageSrc, std::vector<KeyLine>& ke ...@@ -500,7 +508,7 @@ void BinaryDescriptor::detectImpl( const Mat& imageSrc, std::vector<KeyLine>& ke
/* requires descriptors computation (only one image) */ /* requires descriptors computation (only one image) */
void BinaryDescriptor::compute( const Mat& image, CV_OUT CV_IN_OUT std::vector<KeyLine>& keylines, CV_OUT Mat& descriptors, void BinaryDescriptor::compute( const Mat& image, CV_OUT CV_IN_OUT std::vector<KeyLine>& keylines, CV_OUT Mat& descriptors,
bool returnFloatDescr ) const bool returnFloatDescr ) const
{ {
computeImpl( image, keylines, descriptors, returnFloatDescr, false ); computeImpl( image, keylines, descriptors, returnFloatDescr, false );
} }
...@@ -703,7 +711,7 @@ int BinaryDescriptor::OctaveKeyLines( cv::Mat& image, ScaleLines &keyLines ) ...@@ -703,7 +711,7 @@ int BinaryDescriptor::OctaveKeyLines( cv::Mat& image, ScaleLines &keyLines )
} /* end of loop over number of octaves */ } /* end of loop over number of octaves */
/* prepare a vector to store octave information associated to extracted lines */ /* prepare a vector to store octave information associated to extracted lines */
std::vector<OctaveLine> octaveLines( numOfFinalLine ); std::vector < OctaveLine > octaveLines( numOfFinalLine );
/* set lines' counter to 0 for reuse */ /* set lines' counter to 0 for reuse */
numOfFinalLine = 0; numOfFinalLine = 0;
...@@ -1346,5 +1354,1385 @@ int BinaryDescriptor::computeLBD( ScaleLines &keyLines, bool useDetectionData ) ...@@ -1346,5 +1354,1385 @@ int BinaryDescriptor::computeLBD( ScaleLines &keyLines, bool useDetectionData )
return 1; return 1;
} }
BinaryDescriptor::EDLineDetector::EDLineDetector()
{
//set parameters for line segment detection
ksize_ = 15; //15
sigma_ = 30.0; //30
gradienThreshold_ = 80; // ***** ORIGINAL WAS 25
anchorThreshold_ = 8; //8
scanIntervals_ = 2; //2
minLineLen_ = 15; //15
lineFitErrThreshold_ = 1.6; //1.4
InitEDLine_();
}
BinaryDescriptor::EDLineDetector::EDLineDetector( EDLineParam param )
{
//set parameters for line segment detection
ksize_ = param.ksize;
sigma_ = param.sigma;
gradienThreshold_ = (short) param.gradientThreshold;
anchorThreshold_ = (unsigned char) param.anchorThreshold;
scanIntervals_ = param.scanIntervals;
minLineLen_ = param.minLineLen;
lineFitErrThreshold_ = param.lineFitErrThreshold;
InitEDLine_();
}
void BinaryDescriptor::EDLineDetector::InitEDLine_()
{
bValidate_ = true;
ATA = cv::Mat_<int>( 2, 2 );
ATV = cv::Mat_<int>( 1, 2 );
tempMatLineFit = cv::Mat_<int>( 2, 2 );
tempVecLineFit = cv::Mat_<int>( 1, 2 );
fitMatT = cv::Mat_<int>( 2, minLineLen_ );
fitVec = cv::Mat_<int>( 1, minLineLen_ );
for ( int i = 0; i < minLineLen_; i++ )
{
fitMatT[1][i] = 1;
}
dxImg_.create( 1, 1, CV_16SC1 );
dyImg_.create( 1, 1, CV_16SC1 );
gImgWO_.create( 1, 1, CV_8SC1 );
pFirstPartEdgeX_ = NULL;
pFirstPartEdgeY_ = NULL;
pFirstPartEdgeS_ = NULL;
pSecondPartEdgeX_ = NULL;
pSecondPartEdgeY_ = NULL;
pSecondPartEdgeS_ = NULL;
pAnchorX_ = NULL;
pAnchorY_ = NULL;
}
BinaryDescriptor::EDLineDetector::~EDLineDetector()
{
if( pFirstPartEdgeX_ != NULL )
{
delete[] pFirstPartEdgeX_;
delete[] pFirstPartEdgeY_;
delete[] pSecondPartEdgeX_;
delete[] pSecondPartEdgeY_;
delete[] pAnchorX_;
delete[] pAnchorY_;
}
if( pFirstPartEdgeS_ != NULL )
{
delete[] pFirstPartEdgeS_;
delete[] pSecondPartEdgeS_;
}
}
int BinaryDescriptor::EDLineDetector::EdgeDrawing( cv::Mat &image, EdgeChains &edgeChains )
{
imageWidth = image.cols;
imageHeight = image.rows;
unsigned int pixelNum = imageWidth * imageHeight;
unsigned int edgePixelArraySize = pixelNum / 5;
unsigned int maxNumOfEdge = edgePixelArraySize / 20;
//compute dx, dy images
if( gImg_.cols != (int) imageWidth || gImg_.rows != (int) imageHeight )
{
if( pFirstPartEdgeX_ != NULL )
{
delete[] pFirstPartEdgeX_;
delete[] pFirstPartEdgeY_;
delete[] pSecondPartEdgeX_;
delete[] pSecondPartEdgeY_;
delete[] pFirstPartEdgeS_;
delete[] pSecondPartEdgeS_;
delete[] pAnchorX_;
delete[] pAnchorY_;
}
dxImg_.create( imageHeight, imageWidth, CV_16SC1 );
dyImg_.create( imageHeight, imageWidth, CV_16SC1 );
gImgWO_.create( imageHeight, imageWidth, CV_16SC1 );
gImg_.create( imageHeight, imageWidth, CV_16SC1 );
dirImg_.create( imageHeight, imageWidth, CV_8UC1 );
edgeImage_.create( imageHeight, imageWidth, CV_8UC1 );
pFirstPartEdgeX_ = new unsigned int[edgePixelArraySize];
pFirstPartEdgeY_ = new unsigned int[edgePixelArraySize];
pSecondPartEdgeX_ = new unsigned int[edgePixelArraySize];
pSecondPartEdgeY_ = new unsigned int[edgePixelArraySize];
pFirstPartEdgeS_ = new unsigned int[maxNumOfEdge];
pSecondPartEdgeS_ = new unsigned int[maxNumOfEdge];
pAnchorX_ = new unsigned int[edgePixelArraySize];
pAnchorY_ = new unsigned int[edgePixelArraySize];
}
cv::Sobel( image, dxImg_, CV_16SC1, 1, 0, 3 );
cv::Sobel( image, dyImg_, CV_16SC1, 0, 1, 3 );
//compute gradient and direction images
cv::Mat dxABS_m = cv::abs( dxImg_ );
cv::Mat dyABS_m = cv::abs( dyImg_ );
cv::Mat sumDxDy;
cv::add( dyABS_m, dxABS_m, sumDxDy );
cv::threshold( sumDxDy, gImg_, gradienThreshold_ + 1, 255, cv::THRESH_TOZERO );
gImg_ = gImg_ / 4;
gImgWO_ = sumDxDy / 4;
cv::compare( dxABS_m, dyABS_m, dirImg_, cv::CMP_LT );
short *pgImg = gImg_.ptr<short>();
unsigned char *pdirImg = dirImg_.ptr();
//extract the anchors in the gradient image, store into a vector
memset( pAnchorX_, 0, edgePixelArraySize * sizeof(unsigned int) ); //initialization
memset( pAnchorY_, 0, edgePixelArraySize * sizeof(unsigned int) );
unsigned int anchorsSize = 0;
int indexInArray;
unsigned char gValue1, gValue2, gValue3;
for ( unsigned int w = 1; w < imageWidth - 1; w = w + scanIntervals_ )
{
for ( unsigned int h = 1; h < imageHeight - 1; h = h + scanIntervals_ )
{
indexInArray = h * imageWidth + w;
//gValue1 = pdirImg[indexInArray];
if( pdirImg[indexInArray] == Horizontal )
{ //if the direction of pixel is horizontal, then compare with up and down
//gValue2 = pgImg[indexInArray];
if( pgImg[indexInArray] >= pgImg[indexInArray - imageWidth] + anchorThreshold_
&& pgImg[indexInArray] >= pgImg[indexInArray + imageWidth] + anchorThreshold_ )
{ // (w,h) is accepted as an anchor
pAnchorX_[anchorsSize] = w;
pAnchorY_[anchorsSize++] = h;
}
}
else
{ // if(pdirImg[indexInArray]==Vertical){//it is vertical edge, should be compared with left and right
//gValue2 = pgImg[indexInArray];
if( pgImg[indexInArray] >= pgImg[indexInArray - 1] + anchorThreshold_ && pgImg[indexInArray] >= pgImg[indexInArray + 1] + anchorThreshold_ )
{ // (w,h) is accepted as an anchor
pAnchorX_[anchorsSize] = w;
pAnchorY_[anchorsSize++] = h;
}
}
}
}
if( anchorsSize > edgePixelArraySize )
{
std::cout << "anchor size is larger than its maximal size. anchorsSize=" << anchorsSize << ", maximal size = " << edgePixelArraySize << std::endl;
return -1;
}
//link the anchors by smart routing
edgeImage_.setTo( 0 );
unsigned char *pEdgeImg = edgeImage_.data;
memset( pFirstPartEdgeX_, 0, edgePixelArraySize * sizeof(unsigned int) ); //initialization
memset( pFirstPartEdgeY_, 0, edgePixelArraySize * sizeof(unsigned int) );
memset( pSecondPartEdgeX_, 0, edgePixelArraySize * sizeof(unsigned int) );
memset( pSecondPartEdgeY_, 0, edgePixelArraySize * sizeof(unsigned int) );
memset( pFirstPartEdgeS_, 0, maxNumOfEdge * sizeof(unsigned int) );
memset( pSecondPartEdgeS_, 0, maxNumOfEdge * sizeof(unsigned int) );
unsigned int offsetPFirst = 0, offsetPSecond = 0;
unsigned int offsetPS = 0;
unsigned int x, y;
unsigned int lastX = 0;
unsigned int lastY = 0;
unsigned char lastDirection; //up = 1, right = 2, down = 3, left = 4;
unsigned char shouldGoDirection; //up = 1, right = 2, down = 3, left = 4;
int edgeLenFirst, edgeLenSecond;
for ( unsigned int i = 0; i < anchorsSize; i++ )
{
x = pAnchorX_[i];
y = pAnchorY_[i];
indexInArray = y * imageWidth + x;
if( pEdgeImg[indexInArray] )
{ //if anchor i is already been an edge pixel.
continue;
}
/*The walk stops under 3 conditions:
* 1. We move out of the edge areas, i.e., the thresholded gradient value
* of the current pixel is 0.
* 2. The current direction of the edge changes, i.e., from horizontal
* to vertical or vice versa.?? (This is turned out not correct. From the online edge draw demo
* we can figure out that authors don't implement this rule either because their extracted edge
* chain could be a circle which means pixel directions would definitely be different
* in somewhere on the chain.)
* 3. We encounter a previously detected edge pixel. */
pFirstPartEdgeS_[offsetPS] = offsetPFirst;
if( pdirImg[indexInArray] == Horizontal )
{ //if the direction of this pixel is horizontal, then go left and right.
//fist go right, pixel direction may be different during linking.
lastDirection = RightDir;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pFirstPartEdgeX_[offsetPFirst] = x;
pFirstPartEdgeY_[offsetPFirst++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go right
//then go left, pixel direction may be different during linking.
x = pAnchorX_[i];
y = pAnchorY_[i];
indexInArray = y * imageWidth + x;
pEdgeImg[indexInArray] = 0; //mark the anchor point be a non-edge pixel and
lastDirection = LeftDir;
pSecondPartEdgeS_[offsetPS] = offsetPSecond;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pSecondPartEdgeX_[offsetPSecond] = x;
pSecondPartEdgeY_[offsetPSecond++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go left
//end anchor is Horizontal
}
else
{ //the direction of this pixel is vertical, go up and down
//fist go down, pixel direction may be different during linking.
lastDirection = DownDir;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pFirstPartEdgeX_[offsetPFirst] = x;
pFirstPartEdgeY_[offsetPFirst++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go down
//then go up, pixel direction may be different during linking.
lastDirection = UpDir;
x = pAnchorX_[i];
y = pAnchorY_[i];
indexInArray = y * imageWidth + x;
pEdgeImg[indexInArray] = 0; //mark the anchor point be a non-edge pixel and
pSecondPartEdgeS_[offsetPS] = offsetPSecond;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pSecondPartEdgeX_[offsetPSecond] = x;
pSecondPartEdgeY_[offsetPSecond++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go up
} //end anchor is Vertical
//only keep the edge chains whose length is larger than the minLineLen_;
edgeLenFirst = offsetPFirst - pFirstPartEdgeS_[offsetPS];
edgeLenSecond = offsetPSecond - pSecondPartEdgeS_[offsetPS];
if( edgeLenFirst + edgeLenSecond < minLineLen_ + 1 )
{ //short edge, drop it
offsetPFirst = pFirstPartEdgeS_[offsetPS];
offsetPSecond = pSecondPartEdgeS_[offsetPS];
}
else
{
offsetPS++;
}
}
//store the last index
pFirstPartEdgeS_[offsetPS] = offsetPFirst;
pSecondPartEdgeS_[offsetPS] = offsetPSecond;
if( offsetPS > maxNumOfEdge )
{
std::cout << "Edge drawing Error: The total number of edges is larger than MaxNumOfEdge, "
"numofedge = " << offsetPS << ", MaxNumOfEdge=" << maxNumOfEdge << std::endl;
return -1;
}
if( offsetPFirst > edgePixelArraySize || offsetPSecond > edgePixelArraySize )
{
std::cout << "Edge drawing Error: The total number of edge pixels is larger than MaxNumOfEdgePixels, "
"numofedgePixel1 = " << offsetPFirst << ", numofedgePixel2 = " << offsetPSecond << ", MaxNumOfEdgePixel=" << edgePixelArraySize << std::endl;
return -1;
}
/*now all the edge information are stored in pFirstPartEdgeX_, pFirstPartEdgeY_,
*pFirstPartEdgeS_, pSecondPartEdgeX_, pSecondPartEdgeY_, pSecondPartEdgeS_;
*we should reorganize them into edgeChains for easily using. */
int tempID;
edgeChains.xCors.resize( offsetPFirst + offsetPSecond );
edgeChains.yCors.resize( offsetPFirst + offsetPSecond );
edgeChains.sId.resize( offsetPS + 1 );
unsigned int *pxCors = edgeChains.xCors.data();
unsigned int *pyCors = edgeChains.yCors.data();
unsigned int *psId = edgeChains.sId.data();
offsetPFirst = 0;
offsetPSecond = 0;
unsigned int indexInCors = 0;
unsigned int numOfEdges = 0;
for ( unsigned int edgeId = 0; edgeId < offsetPS; edgeId++ )
{
//step1, put the first and second parts edge coordinates together from edge start to edge end
psId[numOfEdges++] = indexInCors;
indexInArray = pFirstPartEdgeS_[edgeId];
offsetPFirst = pFirstPartEdgeS_[edgeId + 1];
for ( tempID = offsetPFirst - 1; tempID >= indexInArray; tempID-- )
{ //add first part edge
pxCors[indexInCors] = pFirstPartEdgeX_[tempID];
pyCors[indexInCors++] = pFirstPartEdgeY_[tempID];
}
indexInArray = pSecondPartEdgeS_[edgeId];
offsetPSecond = pSecondPartEdgeS_[edgeId + 1];
for ( tempID = indexInArray + 1; tempID < (int) offsetPSecond; tempID++ )
{ //add second part edge
pxCors[indexInCors] = pSecondPartEdgeX_[tempID];
pyCors[indexInCors++] = pSecondPartEdgeY_[tempID];
}
}
psId[numOfEdges] = indexInCors; //the end index of the last edge
edgeChains.numOfEdges = numOfEdges;
return 1;
}
int BinaryDescriptor::EDLineDetector::EDline( cv::Mat &image, LineChains &lines )
{
//first, call EdgeDrawing function to extract edges
EdgeChains edges;
if( ( EdgeDrawing( image, edges ) ) != 1 )
{
std::cout << "Line Detection not finished" << std::endl;
return -1;
}
//detect lines
unsigned int linePixelID = edges.sId[edges.numOfEdges];
lines.xCors.resize( linePixelID );
lines.yCors.resize( linePixelID );
lines.sId.resize( 5 * edges.numOfEdges );
unsigned int *pEdgeXCors = edges.xCors.data();
unsigned int *pEdgeYCors = edges.yCors.data();
unsigned int *pEdgeSID = edges.sId.data();
unsigned int *pLineXCors = lines.xCors.data();
unsigned int *pLineYCors = lines.yCors.data();
unsigned int *pLineSID = lines.sId.data();
logNT_ = 2.0 * ( log10( (double) imageWidth ) + log10( (double) imageHeight ) );
double lineFitErr = 0; //the line fit error;
std::vector<double> lineEquation( 2, 0 );
lineEquations_.clear();
lineEndpoints_.clear();
lineDirection_.clear();
unsigned char *pdirImg = dirImg_.data;
unsigned int numOfLines = 0;
unsigned int newOffsetS = 0;
unsigned int offsetInEdgeArrayS, offsetInEdgeArrayE; //start index and end index
unsigned int offsetInLineArray = 0;
float direction; //line direction
for ( unsigned int edgeID = 0; edgeID < edges.numOfEdges; edgeID++ )
{
offsetInEdgeArrayS = pEdgeSID[edgeID];
offsetInEdgeArrayE = pEdgeSID[edgeID + 1];
while ( offsetInEdgeArrayE > offsetInEdgeArrayS + minLineLen_ )
{ //extract line segments from an edge, may find more than one segments
//find an initial line segment
while ( offsetInEdgeArrayE > offsetInEdgeArrayS + minLineLen_ )
{
lineFitErr = LeastSquaresLineFit_( pEdgeXCors, pEdgeYCors, offsetInEdgeArrayS, lineEquation );
if( lineFitErr <= lineFitErrThreshold_ )
break; //ok, an initial line segment detected
offsetInEdgeArrayS += SkipEdgePoint; //skip the first two pixel in the chain and try with the remaining pixels
}
if( lineFitErr > lineFitErrThreshold_ )
break; //no line is detected
//An initial line segment is detected. Try to extend this line segment
pLineSID[numOfLines] = offsetInLineArray;
double coef1 = 0; //for a line ax+by+c=0, coef1 = 1/sqrt(a^2+b^2);
double pointToLineDis; //for a line ax+by+c=0 and a point(xi, yi), pointToLineDis = coef1*|a*xi+b*yi+c|
bool bExtended = true;
bool bFirstTry = true;
int numOfOutlier; //to against noise, we accept a few outlier of a line.
int tryTimes = 0;
if( pdirImg[pEdgeYCors[offsetInEdgeArrayS] * imageWidth + pEdgeXCors[offsetInEdgeArrayS]] == Horizontal )
{ //y=ax+b, i.e. ax-y+b=0
while ( bExtended )
{
tryTimes++;
if( bFirstTry )
{
bFirstTry = false;
for ( int i = 0; i < minLineLen_; i++ )
{ //First add the initial line segment to the line array
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
}
}
else
{ //after each try, line is extended, line equation should be re-estimated
//adjust the line equation
lineFitErr = LeastSquaresLineFit_( pLineXCors, pLineYCors, pLineSID[numOfLines], newOffsetS, offsetInLineArray, lineEquation );
}
coef1 = 1 / sqrt( lineEquation[0] * lineEquation[0] + 1 );
numOfOutlier = 0;
newOffsetS = offsetInLineArray;
while ( offsetInEdgeArrayE > offsetInEdgeArrayS )
{
pointToLineDis = fabs( lineEquation[0] * pEdgeXCors[offsetInEdgeArrayS] - pEdgeYCors[offsetInEdgeArrayS] + lineEquation[1] ) * coef1;
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
if( pointToLineDis > lineFitErrThreshold_ )
{
numOfOutlier++;
if( numOfOutlier > 3 )
break;
}
else
{ //we count number of connective outliers.
numOfOutlier = 0;
}
}
//pop back the last few outliers from lines and return them to edge chain
offsetInLineArray -= numOfOutlier;
offsetInEdgeArrayS -= numOfOutlier;
if( offsetInLineArray - newOffsetS > 0 && tryTimes < TryTime )
{ //some new pixels are added to the line
}
else
{
bExtended = false; //no new pixels are added.
}
}
//the line equation coefficients,for line w1x+w2y+w3 =0, we normalize it to make w1^2+w2^2 = 1.
std::vector<double> lineEqu( 3, 0 );
lineEqu[0] = lineEquation[0] * coef1;
lineEqu[1] = -1 * coef1;
lineEqu[2] = lineEquation[1] * coef1;
if( LineValidation_( pLineXCors, pLineYCors, pLineSID[numOfLines], offsetInLineArray, lineEqu, direction ) )
{ //check the line
//store the line equation coefficients
lineEquations_.push_back( lineEqu );
/*At last, compute the line endpoints and store them.
*we project the first and last pixels in the pixelChain onto the best fit line
*to get the line endpoints.
*xp= (w2^2*x0-w1*w2*y0-w3*w1)/(w1^2+w2^2)
*yp= (w1^2*y0-w1*w2*x0-w3*w2)/(w1^2+w2^2) */
std::vector<float> lineEndP( 4, 0 ); //line endpoints
double a1 = lineEqu[1] * lineEqu[1];
double a2 = lineEqu[0] * lineEqu[0];
double a3 = lineEqu[0] * lineEqu[1];
double a4 = lineEqu[2] * lineEqu[0];
double a5 = lineEqu[2] * lineEqu[1];
unsigned int Px = pLineXCors[pLineSID[numOfLines]]; //first pixel
unsigned int Py = pLineYCors[pLineSID[numOfLines]];
lineEndP[0] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[1] = (float) ( a2 * Py - a3 * Px - a5 ); //y
Px = pLineXCors[offsetInLineArray - 1]; //last pixel
Py = pLineYCors[offsetInLineArray - 1];
lineEndP[2] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[3] = (float) ( a2 * Py - a3 * Px - a5 ); //y
lineEndpoints_.push_back( lineEndP );
lineDirection_.push_back( direction );
numOfLines++;
}
else
{
offsetInLineArray = pLineSID[numOfLines]; // line was not accepted, the offset is set back
}
}
else
{ //x=ay+b, i.e. x-ay-b=0
while ( bExtended )
{
tryTimes++;
if( bFirstTry )
{
bFirstTry = false;
for ( int i = 0; i < minLineLen_; i++ )
{ //First add the initial line segment to the line array
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
}
}
else
{ //after each try, line is extended, line equation should be re-estimated
//adjust the line equation
lineFitErr = LeastSquaresLineFit_( pLineXCors, pLineYCors, pLineSID[numOfLines], newOffsetS, offsetInLineArray, lineEquation );
}
coef1 = 1 / sqrt( 1 + lineEquation[0] * lineEquation[0] );
numOfOutlier = 0;
newOffsetS = offsetInLineArray;
while ( offsetInEdgeArrayE > offsetInEdgeArrayS )
{
pointToLineDis = fabs( pEdgeXCors[offsetInEdgeArrayS] - lineEquation[0] * pEdgeYCors[offsetInEdgeArrayS] - lineEquation[1] ) * coef1;
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
if( pointToLineDis > lineFitErrThreshold_ )
{
numOfOutlier++;
if( numOfOutlier > 3 )
break;
}
else
{ //we count number of connective outliers.
numOfOutlier = 0;
}
}
//pop back the last few outliers from lines and return them to edge chain
offsetInLineArray -= numOfOutlier;
offsetInEdgeArrayS -= numOfOutlier;
if( offsetInLineArray - newOffsetS > 0 && tryTimes < TryTime )
{ //some new pixels are added to the line
}
else
{
bExtended = false; //no new pixels are added.
}
}
//the line equation coefficients,for line w1x+w2y+w3 =0, we normalize it to make w1^2+w2^2 = 1.
std::vector<double> lineEqu( 3, 0 );
lineEqu[0] = 1 * coef1;
lineEqu[1] = -lineEquation[0] * coef1;
lineEqu[2] = -lineEquation[1] * coef1;
if( LineValidation_( pLineXCors, pLineYCors, pLineSID[numOfLines], offsetInLineArray, lineEqu, direction ) )
{ //check the line
//store the line equation coefficients
lineEquations_.push_back( lineEqu );
/*At last, compute the line endpoints and store them.
*we project the first and last pixels in the pixelChain onto the best fit line
*to get the line endpoints.
*xp= (w2^2*x0-w1*w2*y0-w3*w1)/(w1^2+w2^2)
*yp= (w1^2*y0-w1*w2*x0-w3*w2)/(w1^2+w2^2) */
std::vector<float> lineEndP( 4, 0 ); //line endpoints
double a1 = lineEqu[1] * lineEqu[1];
double a2 = lineEqu[0] * lineEqu[0];
double a3 = lineEqu[0] * lineEqu[1];
double a4 = lineEqu[2] * lineEqu[0];
double a5 = lineEqu[2] * lineEqu[1];
unsigned int Px = pLineXCors[pLineSID[numOfLines]]; //first pixel
unsigned int Py = pLineYCors[pLineSID[numOfLines]];
lineEndP[0] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[1] = (float) ( a2 * Py - a3 * Px - a5 ); //y
Px = pLineXCors[offsetInLineArray - 1]; //last pixel
Py = pLineYCors[offsetInLineArray - 1];
lineEndP[2] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[3] = (float) ( a2 * Py - a3 * Px - a5 ); //y
lineEndpoints_.push_back( lineEndP );
lineDirection_.push_back( direction );
numOfLines++;
}
else
{
offsetInLineArray = pLineSID[numOfLines]; // line was not accepted, the offset is set back
}
}
//Extract line segments from the remaining pixel; Current chain has been shortened already.
}
} //end for(unsigned int edgeID=0; edgeID<edges.numOfEdges; edgeID++)
pLineSID[numOfLines] = offsetInLineArray;
lines.numOfLines = numOfLines;
return 1;
}
double BinaryDescriptor::EDLineDetector::LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS,
std::vector<double> &lineEquation )
{
float * pMatT;
float * pATA;
double fitError = 0;
double coef;
unsigned char *pdirImg = dirImg_.data;
unsigned int offset = offsetS;
/*If the first pixel in this chain is horizontal,
*then we try to find a horizontal line, y=ax+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Horizontal )
{
/*Build the system,and solve it using least square regression: mat * [a,b]^T = vec
* [x0,1] [y0]
* [x1,1] [a] [y1]
* . [b] = .
* [xn,1] [yn]*/
pMatT = fitMatT.ptr<float>(); //fitMatT = [x0, x1, ... xn; 1,1,...,1];
for ( int i = 0; i < minLineLen_; i++ )
{
//*(pMatT+minLineLen_) = 1; //the value are not changed;
* ( pMatT++ ) = (float) xCors[offsetS];
fitVec[0][i] = (float) yCors[offsetS++];
}
ATA = fitMatT * fitMatT.t();
ATV = fitMatT * fitVec.t();
/* [a,b]^T = Inv(mat^T * mat) * mat^T * vec */
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
// lineEquation = svd.Invert(ATA) * matT * vec;
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
/*compute line fit error */
for ( int i = 0; i < minLineLen_; i++ )
{
//coef = double(yCors[offset]) - double(xCors[offset++]) * lineEquation[0] - lineEquation[1];
coef = double( yCors[offset] ) - double( xCors[offset] ) * lineEquation[0] - lineEquation[1];
offset++;
fitError += coef * coef;
}
return sqrt( fitError );
}
/*If the first pixel in this chain is vertical,
*then we try to find a vertical line, x=ay+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Vertical )
{
/*Build the system,and solve it using least square regression: mat * [a,b]^T = vec
* [y0,1] [x0]
* [y1,1] [a] [x1]
* . [b] = .
* [yn,1] [xn]*/
pMatT = fitMatT.ptr<float>(); //fitMatT = [y0, y1, ... yn; 1,1,...,1];
for ( int i = 0; i < minLineLen_; i++ )
{
//*(pMatT+minLineLen_) = 1;//the value are not changed;
* ( pMatT++ ) = (float) yCors[offsetS];
fitVec[0][i] = (float) xCors[offsetS++];
}
ATA = fitMatT * ( fitMatT.t() );
ATV = fitMatT * fitVec.t();
/* [a,b]^T = Inv(mat^T * mat) * mat^T * vec */
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
// lineEquation = svd.Invert(ATA) * matT * vec;
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
/*compute line fit error */
for ( int i = 0; i < minLineLen_; i++ )
{
//coef = double(xCors[offset]) - double(yCors[offset++]) * lineEquation[0] - lineEquation[1];
coef = double( xCors[offset] ) - double( yCors[offset] ) * lineEquation[0] - lineEquation[1];
offset++;
fitError += coef * coef;
}
return sqrt( fitError );
}
return 0;
}
double BinaryDescriptor::EDLineDetector::LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS,
unsigned int newOffsetS, unsigned int offsetE, std::vector<double> &lineEquation )
{
int length = offsetE - offsetS;
int newLength = offsetE - newOffsetS;
if( length <= 0 || newLength <= 0 )
{
std::cout << "EDLineDetector::LeastSquaresLineFit_ Error:"
" the expected line index is wrong...offsetE = " << offsetE << ", offsetS=" << offsetS << ", newOffsetS=" << newOffsetS << std::endl;
return -1;
}
if( lineEquation.size() != 2 )
{
std::cout << "SHOULD NOT BE != 2" << std::endl;
}
cv::Mat_<float> matT( 2, newLength );
cv::Mat_<float> vec( newLength, 1 );
float * pMatT;
float * pATA;
double coef;
unsigned char *pdirImg = dirImg_.data;
/*If the first pixel in this chain is horizontal,
*then we try to find a horizontal line, y=ax+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Horizontal )
{
/*Build the new system,and solve it using least square regression: mat * [a,b]^T = vec
* [x0',1] [y0']
* [x1',1] [a] [y1']
* . [b] = .
* [xn',1] [yn']*/
pMatT = matT.ptr<float>(); //matT = [x0', x1', ... xn'; 1,1,...,1]
for ( int i = 0; i < newLength; i++ )
{
* ( pMatT + newLength ) = 1;
* ( pMatT++ ) = (float) xCors[newOffsetS];
vec[0][i] = (float) yCors[newOffsetS++];
}
/* [a,b]^T = Inv(ATA + mat^T * mat) * (ATV + mat^T * vec) */
tempMatLineFit = matT * matT.t();
tempVecLineFit = matT * vec;
ATA = ATA + tempMatLineFit;
ATV = ATV + tempVecLineFit;
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
return 0;
}
/*If the first pixel in this chain is vertical,
*then we try to find a vertical line, x=ay+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Vertical )
{
/*Build the system,and solve it using least square regression: mat * [a,b]^T = vec
* [y0',1] [x0']
* [y1',1] [a] [x1']
* . [b] = .
* [yn',1] [xn']*/
pMatT = matT.ptr<float>(); //matT = [y0', y1', ... yn'; 1,1,...,1]
for ( int i = 0; i < newLength; i++ )
{
* ( pMatT + newLength ) = 1;
* ( pMatT++ ) = (float) yCors[newOffsetS];
vec[0][i] = (float) xCors[newOffsetS++];
}
/* [a,b]^T = Inv(ATA + mat^T * mat) * (ATV + mat^T * vec) */
// matT.MultiplyWithTransposeOf(matT, tempMatLineFit);
tempMatLineFit = matT * matT.t();
tempVecLineFit = matT * vec;
ATA = ATA + tempMatLineFit;
ATV = ATV + tempVecLineFit;
// pATA = ATA.GetData();
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
}
return 0;
}
bool BinaryDescriptor::EDLineDetector::LineValidation_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int offsetE,
std::vector<double> &lineEquation, float &direction )
{
if( bValidate_ )
{
int n = offsetE - offsetS;
/*first compute the direction of line, make sure that the dark side always be the
*left side of a line.*/
int meanGradientX = 0, meanGradientY = 0;
short *pdxImg = dxImg_.ptr<short>();
short *pdyImg = dyImg_.ptr<short>();
double dx, dy;
std::vector<double> pointDirection;
int index;
for ( int i = 0; i < n; i++ )
{
index = yCors[offsetS] * imageWidth + xCors[offsetS];
offsetS++;
meanGradientX += pdxImg[index];
meanGradientY += pdyImg[index];
dx = (double) pdxImg[index];
dy = (double) pdyImg[index];
pointDirection.push_back( atan2( -dx, dy ) );
}
dx = fabs( lineEquation[1] );
dy = fabs( lineEquation[0] );
if( meanGradientX == 0 && meanGradientY == 0 )
{ //not possible, if happens, it must be a wrong line,
return false;
}
if( meanGradientX > 0 && meanGradientY >= 0 )
{ //first quadrant, and positive direction of X axis.
direction = (float) atan2( -dy, dx ); //line direction is in fourth quadrant
}
if( meanGradientX <= 0 && meanGradientY > 0 )
{ //second quadrant, and positive direction of Y axis.
direction = (float) atan2( dy, dx ); //line direction is in first quadrant
}
if( meanGradientX < 0 && meanGradientY <= 0 )
{ //third quadrant, and negative direction of X axis.
direction = (float) atan2( dy, -dx ); //line direction is in second quadrant
}
if( meanGradientX >= 0 && meanGradientY < 0 )
{ //fourth quadrant, and negative direction of Y axis.
direction = (float) atan2( -dy, -dx ); //line direction is in third quadrant
}
/*then check whether the line is on the border of the image. We don't keep the border line.*/
if( fabs( direction ) < 0.15 || M_PI - fabs( direction ) < 0.15 )
{ //Horizontal line
if( fabs( lineEquation[2] ) < 10 || fabs( imageHeight - fabs( lineEquation[2] ) ) < 10 )
{ //upper border or lower border
return false;
}
}
if( fabs( fabs( direction ) - M_PI * 0.5 ) < 0.15 )
{ //Vertical line
if( fabs( lineEquation[2] ) < 10 || fabs( imageWidth - fabs( lineEquation[2] ) ) < 10 )
{ //left border or right border
return false;
}
}
//count the aligned points on the line which have the same direction as the line.
double disDirection;
int k = 0;
for ( int i = 0; i < n; i++ )
{
disDirection = fabs( direction - pointDirection[i] );
if( fabs( 2 * M_PI - disDirection ) < 0.392699 || disDirection < 0.392699 )
{ //same direction, pi/8 = 0.392699081698724
k++;
}
}
//now compute NFA(Number of False Alarms)
double ret = nfa( n, k, 0.125, logNT_ );
return ( ret > 0 ); //0 corresponds to 1 mean false alarm
}
else
{
return true;
}
}
int BinaryDescriptor::EDLineDetector::EDline( cv::Mat &image )
{
if( ( EDline( image, lines_/*, smoothed*/) ) != 1 )
{
return -1;
}
lineSalience_.clear();
lineSalience_.resize( lines_.numOfLines );
unsigned char *pgImg = gImgWO_.ptr();
unsigned int indexInLineArray;
unsigned int *pXCor = lines_.xCors.data();
unsigned int *pYCor = lines_.yCors.data();
unsigned int *pSID = lines_.sId.data();
for ( unsigned int i = 0; i < lineSalience_.size(); i++ )
{
int salience = 0;
for ( indexInLineArray = pSID[i]; indexInLineArray < pSID[i + 1]; indexInLineArray++ )
{
salience += pgImg[pYCor[indexInLineArray] * imageWidth + pXCor[indexInLineArray]];
}
lineSalience_[i] = (float) salience;
}
return 1;
}
} }
} }
modules/line_descriptor/src/binary_descriptor_matcher.cpp
View file @ 9a52f462
...@@ -41,6 +41,10 @@ ...@@ -41,6 +41,10 @@
#include "precomp.hpp" #include "precomp.hpp"
#define MAX_B 37
double ARRAY_RESIZE_FACTOR = 1.1; // minimum is 1.0
double ARRAY_RESIZE_ADD_FACTOR = 4; // minimum is 1
//using namespace cv; //using namespace cv;
namespace cv namespace cv
{ {
...@@ -59,7 +63,7 @@ BinaryDescriptorMatcher::BinaryDescriptorMatcher() ...@@ -59,7 +63,7 @@ BinaryDescriptorMatcher::BinaryDescriptorMatcher()
/* constructor with smart pointer */ /* constructor with smart pointer */
Ptr<BinaryDescriptorMatcher> BinaryDescriptorMatcher::createBinaryDescriptorMatcher() Ptr<BinaryDescriptorMatcher> BinaryDescriptorMatcher::createBinaryDescriptorMatcher()
{ {
return Ptr<BinaryDescriptorMatcher>( new BinaryDescriptorMatcher() ); return Ptr < BinaryDescriptorMatcher > ( new BinaryDescriptorMatcher() );
} }
/* store new descriptors to be inserted in dataset */ /* store new descriptors to be inserted in dataset */
...@@ -132,7 +136,7 @@ void BinaryDescriptorMatcher::match( const Mat& queryDescriptors, std::vector<DM ...@@ -132,7 +136,7 @@ void BinaryDescriptorMatcher::match( const Mat& queryDescriptors, std::vector<DM
if( masks.size() != 0 && (int) masks.size() != numImages ) if( masks.size() != 0 && (int) masks.size() != numImages )
{ {
std::cout << "Error: the number of images in dataset is " << numImages << " but match function received " << masks.size() std::cout << "Error: the number of images in dataset is " << numImages << " but match function received " << masks.size()
<< " masks. Program will be terminated" << std::endl; << " masks. Program will be terminated" << std::endl;
return; return;
} }
...@@ -163,12 +167,12 @@ void BinaryDescriptorMatcher::match( const Mat& queryDescriptors, std::vector<DM ...@@ -163,12 +167,12 @@ void BinaryDescriptorMatcher::match( const Mat& queryDescriptors, std::vector<DM
{ {
std::stringstream ss; std::stringstream ss;
ss << "Error: mask " << itup->second << " in knnMatch function " << "should have " << queryDescriptors.rows << " and " ss << "Error: mask " << itup->second << " in knnMatch function " << "should have " << queryDescriptors.rows << " and "
<< "1 column. Program will be terminated"; << "1 column. Program will be terminated";
//throw std::runtime_error( ss.str() ); //throw std::runtime_error( ss.str() );
} }
/* create a DMatch object if required by mask or if there is /* create a DMatch object if required by mask or if there is
no mask at all */ no mask at all */
else if( masks.empty() || masks[itup->second].at<uchar>( counter ) != 0 ) else if( masks.empty() || masks[itup->second].at < uchar > ( counter ) != 0 )
{ {
std::vector<int> k_distances; std::vector<int> k_distances;
checkKDistances( numres, 1, k_distances, counter, 256 ); checkKDistances( numres, 1, k_distances, counter, 256 );
...@@ -227,7 +231,7 @@ void BinaryDescriptorMatcher::match( const Mat& queryDescriptors, const Mat& tra ...@@ -227,7 +231,7 @@ void BinaryDescriptorMatcher::match( const Mat& queryDescriptors, const Mat& tra
{ {
/* create a DMatch object if required by mask or if there is /* create a DMatch object if required by mask or if there is
no mask at all */ no mask at all */
if( mask.empty() || ( !mask.empty() && mask.at<uchar>( counter ) != 0 ) ) if( mask.empty() || ( !mask.empty() && mask.at < uchar > ( counter ) != 0 ) )
{ {
std::vector<int> k_distances; std::vector<int> k_distances;
checkKDistances( numres, 1, k_distances, counter, 256 ); checkKDistances( numres, 1, k_distances, counter, 256 );
...@@ -291,10 +295,10 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, const Mat& ...@@ -291,10 +295,10 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, const Mat&
for ( int counter = 0; counter < queryDescriptors.rows; counter++ ) for ( int counter = 0; counter < queryDescriptors.rows; counter++ )
{ {
/* initialize a vector of matches */ /* initialize a vector of matches */
std::vector<DMatch> tempVec; std::vector < DMatch > tempVec;
/* chech whether query should be ignored */ /* chech whether query should be ignored */
if( !mask.empty() && mask.at<uchar>( counter ) == 0 ) if( !mask.empty() && mask.at < uchar > ( counter ) == 0 )
{ {
/* if compact result is not requested, add an empty vector */ /* if compact result is not requested, add an empty vector */
if( !compactResult ) if( !compactResult )
...@@ -346,7 +350,7 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, std::vector ...@@ -346,7 +350,7 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, std::vector
if( masks.size() != 0 && (int) masks.size() != numImages ) if( masks.size() != 0 && (int) masks.size() != numImages )
{ {
std::cout << "Error: the number of images in dataset is " << numImages << " but knnMatch function received " << masks.size() std::cout << "Error: the number of images in dataset is " << numImages << " but knnMatch function received " << masks.size()
<< " masks. Program will be terminated" << std::endl; << " masks. Program will be terminated" << std::endl;
return; return;
} }
...@@ -369,7 +373,7 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, std::vector ...@@ -369,7 +373,7 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, std::vector
for ( int counter = 0; counter < queryDescriptors.rows; counter++ ) for ( int counter = 0; counter < queryDescriptors.rows; counter++ )
{ {
/* create a void vector of matches */ /* create a void vector of matches */
std::vector<DMatch> tempVector; std::vector < DMatch > tempVector;
/* loop over k results returned for every query */ /* loop over k results returned for every query */
for ( int j = index; j < index + k; j++ ) for ( int j = index; j < index + k; j++ )
...@@ -384,14 +388,14 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, std::vector ...@@ -384,14 +388,14 @@ void BinaryDescriptorMatcher::knnMatch( const Mat& queryDescriptors, std::vector
if( !masks.empty() && ( masks[itup->second].rows != queryDescriptors.rows || masks[itup->second].cols != 1 ) ) if( !masks.empty() && ( masks[itup->second].rows != queryDescriptors.rows || masks[itup->second].cols != 1 ) )
{ {
std::cout << "Error: mask " << itup->second << " in knnMatch function " << "should have " << queryDescriptors.rows << " and " std::cout << "Error: mask " << itup->second << " in knnMatch function " << "should have " << queryDescriptors.rows << " and "
<< "1 column. Program will be terminated" << std::endl; << "1 column. Program will be terminated" << std::endl;
return; return;
} }
/* decide if, according to relative mask, returned match should be /* decide if, according to relative mask, returned match should be
considered */ considered */
else if( masks.size() == 0 || masks[itup->second].at<uchar>( counter ) != 0 ) else if( masks.size() == 0 || masks[itup->second].at < uchar > ( counter ) != 0 )
{ {
std::vector<int> k_distances; std::vector<int> k_distances;
checkKDistances( numres, k, k_distances, counter, 256 ); checkKDistances( numres, k, k_distances, counter, 256 );
...@@ -465,13 +469,13 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, const Ma ...@@ -465,13 +469,13 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, const Ma
std::vector<int> k_distances; std::vector<int> k_distances;
checkKDistances( numres, trainDescriptors.rows, k_distances, i, 256 ); checkKDistances( numres, trainDescriptors.rows, k_distances, i, 256 );
std::vector<DMatch> tempVector; std::vector < DMatch > tempVector;
for ( int j = index; j < index + trainDescriptors.rows; j++ ) for ( int j = index; j < index + trainDescriptors.rows; j++ )
{ {
// if( numres[j] <= maxDistance ) // if( numres[j] <= maxDistance )
if( k_distances[j - index] <= maxDistance ) if( k_distances[j - index] <= maxDistance )
{ {
if( mask.empty() || mask.at<uchar>( i ) != 0 ) if( mask.empty() || mask.at < uchar > ( i ) != 0 )
{ {
DMatch dm; DMatch dm;
dm.queryIdx = i; dm.queryIdx = i;
...@@ -515,7 +519,7 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec ...@@ -515,7 +519,7 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec
if( masks.size() != 0 && (int) masks.size() != numImages ) if( masks.size() != 0 && (int) masks.size() != numImages )
{ {
std::cout << "Error: the number of images in dataset is " << numImages << " but radiusMatch function received " << masks.size() std::cout << "Error: the number of images in dataset is " << numImages << " but radiusMatch function received " << masks.size()
<< " masks. Program will be terminated" << std::endl; << " masks. Program will be terminated" << std::endl;
return; return;
} }
...@@ -537,7 +541,7 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec ...@@ -537,7 +541,7 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec
int index = 0; int index = 0;
for ( int counter = 0; counter < queryDescriptors.rows; counter++ ) for ( int counter = 0; counter < queryDescriptors.rows; counter++ )
{ {
std::vector<DMatch> tempVector; std::vector < DMatch > tempVector;
for ( int j = index; j < index + descrInDS; j++ ) for ( int j = index; j < index + descrInDS; j++ )
{ {
std::vector<int> k_distances; std::vector<int> k_distances;
...@@ -554,13 +558,13 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec ...@@ -554,13 +558,13 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec
if( !masks.empty() && ( masks[itup->second].rows != queryDescriptors.rows || masks[itup->second].cols != 1 ) ) if( !masks.empty() && ( masks[itup->second].rows != queryDescriptors.rows || masks[itup->second].cols != 1 ) )
{ {
std::cout << "Error: mask " << itup->second << " in radiusMatch function " << "should have " << queryDescriptors.rows << " and " std::cout << "Error: mask " << itup->second << " in radiusMatch function " << "should have " << queryDescriptors.rows << " and "
<< "1 column. Program will be terminated" << std::endl; << "1 column. Program will be terminated" << std::endl;
return; return;
} }
/* add match if necessary */ /* add match if necessary */
else if( masks.empty() || masks[itup->second].at<uchar>( counter ) != 0 ) else if( masks.empty() || masks[itup->second].at < uchar > ( counter ) != 0 )
{ {
DMatch dm; DMatch dm;
...@@ -587,5 +591,477 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec ...@@ -587,5 +591,477 @@ void BinaryDescriptorMatcher::radiusMatch( const Mat& queryDescriptors, std::vec
delete numres; delete numres;
} }
/* execute a batch query */
void BinaryDescriptorMatcher::Mihasher::batchquery( UINT32 * results, UINT32 *numres, const cv::Mat & queries, UINT32 numq, int dim1queries )
{
/* create and initialize a bitarray */
counter = new bitarray;
counter->init( N );
UINT32 *res = new UINT32[K * ( D + 1 )];
UINT64 *chunks = new UINT64[m];
UINT32 * presults = results;
UINT32 *pnumres = numres;
/* make a copy of input queries */
cv::Mat queries_clone = queries.clone();
/* set a pointer to first query (row) */
UINT8 *pq = queries_clone.ptr();
/* loop over number of descriptors */
for ( size_t i = 0; i < numq; i++ )
{
/* for every descriptor, query database */
query( presults, pnumres, pq, chunks, res );
/* move pointer to write next K indeces */
presults += K;
pnumres += B + 1;
/* move forward pointer to current row in descriptors matrix */
pq += dim1queries;
}
delete[] res;
delete[] chunks;
delete counter;
}
/* execute a single query */
void BinaryDescriptorMatcher::Mihasher::query( UINT32* results, UINT32* numres, UINT8 * Query, UINT64 *chunks, UINT32 *res )
{
/* if K == 0 that means we want everything to be processed.
So maxres = N in that case. Otherwise K limits the results processed */
UINT32 maxres = K ? K : (UINT32) N;
/* number of results so far obtained (up to a distance of s per chunk) */
UINT32 n = 0;
/* number of candidates tested with full codes (not counting duplicates) */
UINT32 nc = 0;
/* counting everything retrieved (duplicates are counted multiple times)
number of lookups (and xors) */
UINT32 nl = 0;
UINT32 nd = 0;
UINT32 *arr;
int size = 0;
UINT32 index;
int hammd;
counter->erase();
memset( numres, 0, ( B + 1 ) * sizeof ( *numres ) );
split( chunks, Query, m, mplus, b );
/* the growing search radius per substring */
int s;
/* current b: for the first mplus substrings it is b, for the rest it is (b-1) */
int curb = b;
for ( s = 0; s <= d && n < maxres; s++ )
{
for ( int k = 0; k < m; k++ )
{
if( k < mplus )
curb = b;
else
curb = b - 1;
UINT64 chunksk = chunks[k];
/* number of bit-strings with s number of 1s */
nl += xornum[s + 1] - xornum[s];
/* the bit-string with s number of 1s */
UINT64 bitstr = 0;
for ( int i = 0; i < s; i++ )
/* power[i] stores the location of the i'th 1 */
power[i] = i;
/* used for stopping criterion (location of (s+1)th 1) */
power[s] = curb + 1;
/* bit determines the 1 that should be moving to the left */
int bit = s - 1;
/* start from the left-most 1, and move it to the left until
it touches another one */
/* the loop for changing bitstr */
bool infiniteWhile = true;
while ( infiniteWhile )
{
if( bit != -1 )
{
bitstr ^= ( power[bit] == bit ) ? (UINT64) 1 << power[bit] : (UINT64) 3 << ( power[bit] - 1 );
power[bit]++;
bit--;
}
else
{ /* bit == -1 */
/* the binary code bitstr is available for processing */
arr = H[k].query( chunksk ^ bitstr, &size ); // lookup
if( size )
{ /* the corresponding bucket is not empty */
nd += size;
for ( int c = 0; c < size; c++ )
{
index = arr[c];
if( !counter->get( index ) )
{ /* if it is not a duplicate */
counter->set( index );
hammd = cv::line_descriptor::match( codes.ptr() + (UINT64) index * ( B_over_8 ), Query, B_over_8 );
nc++;
if( hammd <= D && numres[hammd] < maxres )
res[hammd * K + numres[hammd]] = index + 1;
numres[hammd]++;
}
}
}
/* end of processing */
while ( ++bit < s && power[bit] == power[bit + 1] - 1 )
{
bitstr ^= (UINT64) 1 << ( power[bit] - 1 );
power[bit] = bit;
}
if( bit == s )
break;
}
}
n = n + numres[s * m + k];
if( n >= maxres )
break;
}
}
n = 0;
for ( s = 0; s <= D && (int) n < K; s++ )
{
for ( int c = 0; c < (int) numres[s] && (int) n < K; c++ )
results[n++] = res[s * K + c];
}
}
/* constructor 2 */
BinaryDescriptorMatcher::Mihasher::Mihasher( int _B, int _m )
{
B = _B;
B_over_8 = B / 8;
m = _m;
b = (int) ceil( (double) B / m );
/* assuming that B/2 is large enough radius to include
all of the k nearest neighbors */
D = (int) ceil( B / 2.0 );
d = (int) ceil( (double) D / m );
/* mplus is the number of chunks with b bits
(m-mplus) is the number of chunks with (b-1) bits */
mplus = B - m * ( b - 1 );
xornum = new UINT32[d + 2];
xornum[0] = 0;
for ( int i = 0; i <= d; i++ )
xornum[i + 1] = xornum[i] + (UINT32) choose( b, i );
H = new SparseHashtable[m];
/* H[i].init might fail */
for ( int i = 0; i < mplus; i++ )
H[i].init( b );
for ( int i = mplus; i < m; i++ )
H[i].init( b - 1 );
}
/* K setter */
void BinaryDescriptorMatcher::Mihasher::setK( int _K )
{
K = _K;
}
/* desctructor */
BinaryDescriptorMatcher::Mihasher::~Mihasher()
{
delete[] xornum;
delete[] H;
}
/* populate tables */
void BinaryDescriptorMatcher::Mihasher::populate( cv::Mat & _codes, UINT32 _N, int dim1codes )
{
N = _N;
codes = _codes;
UINT64 * chunks = new UINT64[m];
UINT8 * pcodes = codes.ptr();
for ( UINT64 i = 0; i < N; i++, pcodes += dim1codes )
{
split( chunks, pcodes, m, mplus, b );
for ( int k = 0; k < m; k++ )
H[k].insert( chunks[k], (UINT32) i );
if( i % (int) ceil( N / 1000.0 ) == 0 )
fflush (stdout);
}
delete[] chunks;
}
/* constructor */
BinaryDescriptorMatcher::SparseHashtable::SparseHashtable()
{
table = NULL;
size = 0;
b = 0;
}
/* initializer */
int BinaryDescriptorMatcher::SparseHashtable::init( int _b )
{
b = _b;
if( b < 5 || b > MAX_B || b > (int) ( sizeof(UINT64) * 8 ) )
return 1;
size = UINT64_1 << ( b - 5 ); // size = 2 ^ b
table = (BucketGroup*) calloc( size, sizeof(BucketGroup) );
return 0;
}
/* destructor */
BinaryDescriptorMatcher::SparseHashtable::~SparseHashtable()
{
free( table );
}
/* insert data */
void BinaryDescriptorMatcher::SparseHashtable::insert( UINT64 index, UINT32 data )
{
table[index >> 5].insert( (int) ( index % 32 ), data );
}
/* query data */
UINT32* BinaryDescriptorMatcher::SparseHashtable::query( UINT64 index, int *Size )
{
return table[index >> 5].query( (int) ( index % 32 ), Size );
}
/* constructor */
BinaryDescriptorMatcher::BucketGroup::BucketGroup()
{
empty = 0;
group = NULL;
}
/* destructor */
BinaryDescriptorMatcher::BucketGroup::~BucketGroup()
{
if( group != NULL )
delete group;
}
/* insert data into the bucket */
void BinaryDescriptorMatcher::BucketGroup::insert( int subindex, UINT32 data )
{
if( group == NULL )
{
group = new Array32();
group->push( 0 );
}
UINT32 lowerbits = ( (UINT32) 1 << subindex ) - 1;
int end = popcnt( empty & lowerbits );
if( ! ( empty & ( (UINT32) 1 << subindex ) ) )
{
group->insert( end, group->arr[end + 2] );
empty |= (UINT32) 1 << subindex;
}
int totones = popcnt( empty );
group->insert( totones + 1 + group->arr[2 + end + 1], data );
for ( int i = end + 1; i < totones + 1; i++ )
group->arr[2 + i]++;
}
/* perform a query to the bucket */
UINT32* BinaryDescriptorMatcher::BucketGroup::query( int subindex, int *size )
{
if( empty & ( (UINT32) 1 << subindex ) )
{
UINT32 lowerbits = ( (UINT32) 1 << subindex ) - 1;
int end = popcnt( empty & lowerbits );
int totones = popcnt( empty );
*size = group->arr[2 + end + 1] - group->arr[2 + end];
return group->arr + 2 + totones + 1 + group->arr[2 + end];
}
else
{
*size = 0;
return NULL;
}
}
/* set ARRAY_RESIZE_FACTOR */
void BinaryDescriptorMatcher::Array32::setArrayResizeFactor( double arf )
{
ARRAY_RESIZE_FACTOR = arf;
}
/* constructor */
BinaryDescriptorMatcher::Array32::Array32()
{
arr = NULL;
ARRAY_RESIZE_FACTOR = 1.1; // minimum is 1.0
ARRAY_RESIZE_ADD_FACTOR = 4; // minimum is 1
}
/* definition of operator =
Array32& Array32::operator = (const Array32 &rhs) */
void BinaryDescriptorMatcher::Array32::operator =( const Array32 &rhs )
{
if( &rhs != this )
this->arr = rhs.arr;
}
/* destructor */
BinaryDescriptorMatcher::Array32::~Array32()
{
cleanup();
}
/* cleaning function used in destructor */
void BinaryDescriptorMatcher::Array32::cleanup()
{
free( arr );
}
/* push data */
void BinaryDescriptorMatcher::Array32::push( UINT32 Data )
{
if( arr )
{
if( arr[0] == arr[1] )
{
arr[1] = (UINT32) std::max( ceil( arr[1] * ARRAY_RESIZE_FACTOR ), arr[1] + ARRAY_RESIZE_ADD_FACTOR );
UINT32* new_Data = static_cast<UINT32*>( realloc( arr, sizeof(UINT32) * ( 2 + arr[1] ) ) );
if( new_Data == NULL )
{
/* could not realloc, but orig still valid */
std::cout << "ALERT!!!! Not enough memory, operation aborted!" << std::endl;
exit( 0 );
}
else
{
arr = new_Data;
}
}
arr[2 + arr[0]] = Data;
arr[0]++;
}
else
{
arr = (UINT32*) malloc( ( size_t )( 2 + ARRAY_RESIZE_ADD_FACTOR ) * sizeof(UINT32) );
arr[0] = 1;
arr[1] = 1;
arr[2] = Data;
}
}
/* insert data at given index */
void BinaryDescriptorMatcher::Array32::insert( UINT32 index, UINT32 Data )
{
if( arr )
{
if( arr[0] == arr[1] )
{
arr[1] = (UINT32) ceil( arr[0] * 1.1 );
UINT32* new_data = static_cast<UINT32*>( realloc( arr, sizeof(UINT32) * ( 2 + arr[1] ) ) );
if( new_data == NULL )
{
// could not realloc, but orig still valid
std::cout << "ALERT!!!! Not enough memory, operation aborted!" << std::endl;
exit( 0 );
}
else
{
arr = new_data;
}
}
memmove( arr + ( 2 + index ) + 1, arr + ( 2 + index ), ( arr[0] - index ) * sizeof ( *arr ) );
arr[2 + index] = Data;
arr[0]++;
}
else
{
arr = (UINT32*) malloc( 3 * sizeof(UINT32) );
arr[0] = 1;
arr[1] = 1;
arr[2] = Data;
}
}
/* return data */
UINT32* BinaryDescriptorMatcher::Array32::data()
{
return arr ? arr + 2 : NULL;
}
/* return data size */
UINT32 BinaryDescriptorMatcher::Array32::size()
{
return arr ? arr[0] : 0;
}
/* return capacity */
UINT32 BinaryDescriptorMatcher::Array32::capacity()
{
return arr ? arr[1] : 0;
}
/* print data */
void BinaryDescriptorMatcher::Array32::print()
{
for ( int i = 0; i < (int) size(); i++ )
printf( "%d, ", arr[i + 2] );
printf( "\n" );
}
/* initializer */
void BinaryDescriptorMatcher::Array32::init( int Size )
{
if( arr == NULL )
{
arr = (UINT32*) malloc( ( 2 + Size ) * sizeof(UINT32) );
arr[0] = 0;
arr[1] = Size;
}
} }
} }
}
modules/line_descriptor/include/opencv2/line_descriptor/bitops.hpp → modules/line_descriptor/src/bitops.hpp
View file @ 9a52f462
...@@ -43,6 +43,8 @@ ...@@ -43,6 +43,8 @@
#ifndef __OPENCV_BITOPTS_HPP #ifndef __OPENCV_BITOPTS_HPP
#define __OPENCV_BITOPTS_HPP #define __OPENCV_BITOPTS_HPP
#include "precomp.hpp"
#ifdef _WIN32 #ifdef _WIN32
# include <intrin.h> # include <intrin.h>
# define popcnt __popcnt # define popcnt __popcnt
...@@ -63,6 +65,10 @@ const int lookup[] = ...@@ -63,6 +65,10 @@ const int lookup[] =
3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 4, 5, 5, 6, 5, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 4, 5, 5, 6, 5,
6, 6, 7, 5, 6, 6, 7, 6, 7, 7, 8 }; 6, 6, 7, 5, 6, 6, 7, 6, 7, 7, 8 };
namespace cv
{
namespace line_descriptor
{
/*matching function */ /*matching function */
inline int match( UINT8*P, UINT8*Q, int codelb ) inline int match( UINT8*P, UINT8*Q, int codelb )
{ {
...@@ -163,5 +169,7 @@ inline UINT64 choose( int n, int r ) ...@@ -163,5 +169,7 @@ inline UINT64 choose( int n, int r )
return nchooser; return nchooser;
} }
}
}
#endif #endif
modules/line_descriptor/src/bucket_group.cpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#include "precomp.hpp"
/* constructor */
BucketGroup::BucketGroup()
{
empty = 0;
group = NULL;
}
/* destructor */
BucketGroup::~BucketGroup()
{
if( group != NULL )
delete group;
}
/* insert data into the bucket */
void BucketGroup::insert( int subindex, UINT32 data )
{
if( group == NULL )
{
group = new Array32();
group->push( 0 );
}
UINT32 lowerbits = ( (UINT32) 1 << subindex ) - 1;
int end = popcnt( empty & lowerbits );
if( ! ( empty & ( (UINT32) 1 << subindex ) ) )
{
group->insert( end, group->arr[end + 2] );
empty |= (UINT32) 1 << subindex;
}
int totones = popcnt( empty );
group->insert( totones + 1 + group->arr[2 + end + 1], data );
for ( int i = end + 1; i < totones + 1; i++ )
group->arr[2 + i]++;
}
/* perform a query to the bucket */
UINT32* BucketGroup::query( int subindex, int *size )
{
if( empty & ( (UINT32) 1 << subindex ) )
{
UINT32 lowerbits = ( (UINT32) 1 << subindex ) - 1;
int end = popcnt( empty & lowerbits );
int totones = popcnt( empty );
*size = group->arr[2 + end + 1] - group->arr[2 + end];
return group->arr + 2 + totones + 1 + group->arr[2 + end];
}
else
{
*size = 0;
return NULL;
}
}
modules/line_descriptor/src/ed_line_detector.cpp deleted 100644 → 0
View file @ e39065bf
/*IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
By downloading, copying, installing or using the software you agree to this license.
If you do not agree to this license, do not download, install,
copy or use the software.
License Agreement
For Open Source Computer Vision Library
Copyright (C) 2011-2012, Lilian Zhang, all rights reserved.
Copyright (C) 2013, Manuele Tamburrano, Stefano Fabri, all rights reserved.
Third party copyrights are property of their respective owners.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright notice,
this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
* The name of the copyright holders may not be used to endorse or promote products
derived from this software without specific prior written permission.
This software is provided by the copyright holders and contributors "as is" and
any express or implied warranties, including, but not limited to, the implied
warranties of merchantability and fitness for a particular purpose are disclaimed.
In no event shall the Intel Corporation or contributors be liable for any direct,
indirect, incidental, special, exemplary, or consequential damages
(including, but not limited to, procurement of substitute goods or services;
loss of use, data, or profits; or business interruption) however caused
and on any theory of liability, whether in contract, strict liability,
or tort (including negligence or otherwise) arising in any way out of
the use of this software, even if advised of the possibility of such damage.
*/
#include "precomp.hpp"
#define Horizontal 255//if |dx|<|dy|;
#define Vertical 0//if |dy|<=|dx|;
#define UpDir 1
#define RightDir 2
#define DownDir 3
#define LeftDir 4
#define TryTime 6
#define SkipEdgePoint 2
using namespace std;
EDLineDetector::EDLineDetector()
{
//set parameters for line segment detection
ksize_ = 15; //15
sigma_ = 30.0; //30
gradienThreshold_ = 80; // ***** ORIGINAL WAS 25
anchorThreshold_ = 8; //8
scanIntervals_ = 2; //2
minLineLen_ = 15; //15
lineFitErrThreshold_ = 1.6; //1.4
InitEDLine_();
}
EDLineDetector::EDLineDetector( EDLineParam param )
{
//set parameters for line segment detection
ksize_ = param.ksize;
sigma_ = param.sigma;
gradienThreshold_ = (short) param.gradientThreshold;
anchorThreshold_ = (unsigned char) param.anchorThreshold;
scanIntervals_ = param.scanIntervals;
minLineLen_ = param.minLineLen;
lineFitErrThreshold_ = param.lineFitErrThreshold;
InitEDLine_();
}
void EDLineDetector::InitEDLine_()
{
bValidate_ = true;
ATA = cv::Mat_<int>( 2, 2 );
ATV = cv::Mat_<int>( 1, 2 );
tempMatLineFit = cv::Mat_<int>( 2, 2 );
tempVecLineFit = cv::Mat_<int>( 1, 2 );
fitMatT = cv::Mat_<int>( 2, minLineLen_ );
fitVec = cv::Mat_<int>( 1, minLineLen_ );
for ( int i = 0; i < minLineLen_; i++ )
{
fitMatT[1][i] = 1;
}
dxImg_.create( 1, 1, CV_16SC1 );
dyImg_.create( 1, 1, CV_16SC1 );
gImgWO_.create( 1, 1, CV_8SC1 );
pFirstPartEdgeX_ = NULL;
pFirstPartEdgeY_ = NULL;
pFirstPartEdgeS_ = NULL;
pSecondPartEdgeX_ = NULL;
pSecondPartEdgeY_ = NULL;
pSecondPartEdgeS_ = NULL;
pAnchorX_ = NULL;
pAnchorY_ = NULL;
}
EDLineDetector::~EDLineDetector()
{
if( pFirstPartEdgeX_ != NULL )
{
delete[] pFirstPartEdgeX_;
delete[] pFirstPartEdgeY_;
delete[] pSecondPartEdgeX_;
delete[] pSecondPartEdgeY_;
delete[] pAnchorX_;
delete[] pAnchorY_;
}
if( pFirstPartEdgeS_ != NULL )
{
delete[] pFirstPartEdgeS_;
delete[] pSecondPartEdgeS_;
}
}
int EDLineDetector::EdgeDrawing( cv::Mat &image, EdgeChains &edgeChains )
{
imageWidth = image.cols;
imageHeight = image.rows;
unsigned int pixelNum = imageWidth * imageHeight;
unsigned int edgePixelArraySize = pixelNum / 5;
unsigned int maxNumOfEdge = edgePixelArraySize / 20;
//compute dx, dy images
if( gImg_.cols != (int) imageWidth || gImg_.rows != (int) imageHeight )
{
if( pFirstPartEdgeX_ != NULL )
{
delete[] pFirstPartEdgeX_;
delete[] pFirstPartEdgeY_;
delete[] pSecondPartEdgeX_;
delete[] pSecondPartEdgeY_;
delete[] pFirstPartEdgeS_;
delete[] pSecondPartEdgeS_;
delete[] pAnchorX_;
delete[] pAnchorY_;
}
dxImg_.create( imageHeight, imageWidth, CV_16SC1 );
dyImg_.create( imageHeight, imageWidth, CV_16SC1 );
gImgWO_.create( imageHeight, imageWidth, CV_16SC1 );
gImg_.create( imageHeight, imageWidth, CV_16SC1 );
dirImg_.create( imageHeight, imageWidth, CV_8UC1 );
edgeImage_.create( imageHeight, imageWidth, CV_8UC1 );
pFirstPartEdgeX_ = new unsigned int[edgePixelArraySize];
pFirstPartEdgeY_ = new unsigned int[edgePixelArraySize];
pSecondPartEdgeX_ = new unsigned int[edgePixelArraySize];
pSecondPartEdgeY_ = new unsigned int[edgePixelArraySize];
pFirstPartEdgeS_ = new unsigned int[maxNumOfEdge];
pSecondPartEdgeS_ = new unsigned int[maxNumOfEdge];
pAnchorX_ = new unsigned int[edgePixelArraySize];
pAnchorY_ = new unsigned int[edgePixelArraySize];
}
cv::Sobel( image, dxImg_, CV_16SC1, 1, 0, 3 );
cv::Sobel( image, dyImg_, CV_16SC1, 0, 1, 3 );
//compute gradient and direction images
cv::Mat dxABS_m = cv::abs( dxImg_ );
cv::Mat dyABS_m = cv::abs( dyImg_ );
cv::Mat sumDxDy;
cv::add( dyABS_m, dxABS_m, sumDxDy );
cv::threshold( sumDxDy, gImg_, gradienThreshold_ + 1, 255, cv::THRESH_TOZERO );
gImg_ = gImg_ / 4;
gImgWO_ = sumDxDy / 4;
cv::compare( dxABS_m, dyABS_m, dirImg_, cv::CMP_LT );
short *pgImg = gImg_.ptr<short>();
unsigned char *pdirImg = dirImg_.ptr();
//extract the anchors in the gradient image, store into a vector
memset( pAnchorX_, 0, edgePixelArraySize * sizeof(unsigned int) ); //initialization
memset( pAnchorY_, 0, edgePixelArraySize * sizeof(unsigned int) );
unsigned int anchorsSize = 0;
int indexInArray;
unsigned char gValue1, gValue2, gValue3;
for ( unsigned int w = 1; w < imageWidth - 1; w = w + scanIntervals_ )
{
for ( unsigned int h = 1; h < imageHeight - 1; h = h + scanIntervals_ )
{
indexInArray = h * imageWidth + w;
//gValue1 = pdirImg[indexInArray];
if( pdirImg[indexInArray] == Horizontal )
{ //if the direction of pixel is horizontal, then compare with up and down
//gValue2 = pgImg[indexInArray];
if( pgImg[indexInArray] >= pgImg[indexInArray - imageWidth] + anchorThreshold_
&& pgImg[indexInArray] >= pgImg[indexInArray + imageWidth] + anchorThreshold_ )
{ // (w,h) is accepted as an anchor
pAnchorX_[anchorsSize] = w;
pAnchorY_[anchorsSize++] = h;
}
}
else
{ // if(pdirImg[indexInArray]==Vertical){//it is vertical edge, should be compared with left and right
//gValue2 = pgImg[indexInArray];
if( pgImg[indexInArray] >= pgImg[indexInArray - 1] + anchorThreshold_ && pgImg[indexInArray] >= pgImg[indexInArray + 1] + anchorThreshold_ )
{ // (w,h) is accepted as an anchor
pAnchorX_[anchorsSize] = w;
pAnchorY_[anchorsSize++] = h;
}
}
}
}
if( anchorsSize > edgePixelArraySize )
{
cout << "anchor size is larger than its maximal size. anchorsSize=" << anchorsSize << ", maximal size = " << edgePixelArraySize << endl;
return -1;
}
//link the anchors by smart routing
edgeImage_.setTo( 0 );
unsigned char *pEdgeImg = edgeImage_.data;
memset( pFirstPartEdgeX_, 0, edgePixelArraySize * sizeof(unsigned int) ); //initialization
memset( pFirstPartEdgeY_, 0, edgePixelArraySize * sizeof(unsigned int) );
memset( pSecondPartEdgeX_, 0, edgePixelArraySize * sizeof(unsigned int) );
memset( pSecondPartEdgeY_, 0, edgePixelArraySize * sizeof(unsigned int) );
memset( pFirstPartEdgeS_, 0, maxNumOfEdge * sizeof(unsigned int) );
memset( pSecondPartEdgeS_, 0, maxNumOfEdge * sizeof(unsigned int) );
unsigned int offsetPFirst = 0, offsetPSecond = 0;
unsigned int offsetPS = 0;
unsigned int x, y;
unsigned int lastX = 0;
unsigned int lastY = 0;
unsigned char lastDirection; //up = 1, right = 2, down = 3, left = 4;
unsigned char shouldGoDirection; //up = 1, right = 2, down = 3, left = 4;
int edgeLenFirst, edgeLenSecond;
for ( unsigned int i = 0; i < anchorsSize; i++ )
{
x = pAnchorX_[i];
y = pAnchorY_[i];
indexInArray = y * imageWidth + x;
if( pEdgeImg[indexInArray] )
{ //if anchor i is already been an edge pixel.
continue;
}
/*The walk stops under 3 conditions:
* 1. We move out of the edge areas, i.e., the thresholded gradient value
* of the current pixel is 0.
* 2. The current direction of the edge changes, i.e., from horizontal
* to vertical or vice versa.?? (This is turned out not correct. From the online edge draw demo
* we can figure out that authors don't implement this rule either because their extracted edge
* chain could be a circle which means pixel directions would definitely be different
* in somewhere on the chain.)
* 3. We encounter a previously detected edge pixel. */
pFirstPartEdgeS_[offsetPS] = offsetPFirst;
if( pdirImg[indexInArray] == Horizontal )
{ //if the direction of this pixel is horizontal, then go left and right.
//fist go right, pixel direction may be different during linking.
lastDirection = RightDir;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pFirstPartEdgeX_[offsetPFirst] = x;
pFirstPartEdgeY_[offsetPFirst++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go right
//then go left, pixel direction may be different during linking.
x = pAnchorX_[i];
y = pAnchorY_[i];
indexInArray = y * imageWidth + x;
pEdgeImg[indexInArray] = 0; //mark the anchor point be a non-edge pixel and
lastDirection = LeftDir;
pSecondPartEdgeS_[offsetPS] = offsetPSecond;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pSecondPartEdgeX_[offsetPSecond] = x;
pSecondPartEdgeY_[offsetPSecond++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go left
//end anchor is Horizontal
}
else
{ //the direction of this pixel is vertical, go up and down
//fist go down, pixel direction may be different during linking.
lastDirection = DownDir;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pFirstPartEdgeX_[offsetPFirst] = x;
pFirstPartEdgeY_[offsetPFirst++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go down
//then go up, pixel direction may be different during linking.
lastDirection = UpDir;
x = pAnchorX_[i];
y = pAnchorY_[i];
indexInArray = y * imageWidth + x;
pEdgeImg[indexInArray] = 0; //mark the anchor point be a non-edge pixel and
pSecondPartEdgeS_[offsetPS] = offsetPSecond;
while ( pgImg[indexInArray] > 0 && !pEdgeImg[indexInArray] )
{
pEdgeImg[indexInArray] = 1; // Mark this pixel as an edge pixel
pSecondPartEdgeX_[offsetPSecond] = x;
pSecondPartEdgeY_[offsetPSecond++] = y;
shouldGoDirection = 0; //unknown
if( pdirImg[indexInArray] == Horizontal )
{ //should go left or right
if( lastDirection == UpDir || lastDirection == DownDir )
{ //change the pixel direction now
if( x > lastX )
{ //should go right
shouldGoDirection = RightDir;
}
else
{ //should go left
shouldGoDirection = LeftDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == RightDir || shouldGoDirection == RightDir )
{ //go right
if( x == imageWidth - 1 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the right and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-right
x = x + 1;
y = y + 1;
}
else
{ //straight-right
x = x + 1;
}
lastDirection = RightDir;
}
else if( lastDirection == LeftDir || shouldGoDirection == LeftDir )
{ //go left
if( x == 0 || y == 0 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the left and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
gValue2 = (unsigned char) pgImg[indexInArray - 1];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-left
x = x - 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-left
x = x - 1;
}
lastDirection = LeftDir;
}
}
else
{ //should go up or down.
if( lastDirection == RightDir || lastDirection == LeftDir )
{ //change the pixel direction now
if( y > lastY )
{ //should go down
shouldGoDirection = DownDir;
}
else
{ //should go up
shouldGoDirection = UpDir;
}
}
lastX = x;
lastY = y;
if( lastDirection == DownDir || shouldGoDirection == DownDir )
{ //go down
if( x == 0 || x == imageWidth - 1 || y == imageHeight - 1 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the down and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray + imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray + imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray + imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //down-right
x = x + 1;
y = y + 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //down-left
x = x - 1;
y = y + 1;
}
else
{ //straight-down
y = y + 1;
}
lastDirection = DownDir;
}
else if( lastDirection == UpDir || shouldGoDirection == UpDir )
{ //go up
if( x == 0 || x == imageWidth - 1 || y == 0 )
{ //reach the image border
break;
}
// Look at 3 neighbors to the up and pick the one with the max. gradient value
gValue1 = (unsigned char) pgImg[indexInArray - imageWidth + 1];
gValue2 = (unsigned char) pgImg[indexInArray - imageWidth];
gValue3 = (unsigned char) pgImg[indexInArray - imageWidth - 1];
if( gValue1 >= gValue2 && gValue1 >= gValue3 )
{ //up-right
x = x + 1;
y = y - 1;
}
else if( gValue3 >= gValue2 && gValue3 >= gValue1 )
{ //up-left
x = x - 1;
y = y - 1;
}
else
{ //straight-up
y = y - 1;
}
lastDirection = UpDir;
}
}
indexInArray = y * imageWidth + x;
} //end while go up
} //end anchor is Vertical
//only keep the edge chains whose length is larger than the minLineLen_;
edgeLenFirst = offsetPFirst - pFirstPartEdgeS_[offsetPS];
edgeLenSecond = offsetPSecond - pSecondPartEdgeS_[offsetPS];
if( edgeLenFirst + edgeLenSecond < minLineLen_ + 1 )
{ //short edge, drop it
offsetPFirst = pFirstPartEdgeS_[offsetPS];
offsetPSecond = pSecondPartEdgeS_[offsetPS];
}
else
{
offsetPS++;
}
}
//store the last index
pFirstPartEdgeS_[offsetPS] = offsetPFirst;
pSecondPartEdgeS_[offsetPS] = offsetPSecond;
if( offsetPS > maxNumOfEdge )
{
cout << "Edge drawing Error: The total number of edges is larger than MaxNumOfEdge, "
"numofedge = "
<< offsetPS << ", MaxNumOfEdge=" << maxNumOfEdge << endl;
return -1;
}
if( offsetPFirst > edgePixelArraySize || offsetPSecond > edgePixelArraySize )
{
cout << "Edge drawing Error: The total number of edge pixels is larger than MaxNumOfEdgePixels, "
"numofedgePixel1 = "
<< offsetPFirst << ", numofedgePixel2 = " << offsetPSecond << ", MaxNumOfEdgePixel=" << edgePixelArraySize << endl;
return -1;
}
/*now all the edge information are stored in pFirstPartEdgeX_, pFirstPartEdgeY_,
*pFirstPartEdgeS_, pSecondPartEdgeX_, pSecondPartEdgeY_, pSecondPartEdgeS_;
*we should reorganize them into edgeChains for easily using. */
int tempID;
edgeChains.xCors.resize( offsetPFirst + offsetPSecond );
edgeChains.yCors.resize( offsetPFirst + offsetPSecond );
edgeChains.sId.resize( offsetPS + 1 );
unsigned int *pxCors = edgeChains.xCors.data();
unsigned int *pyCors = edgeChains.yCors.data();
unsigned int *psId = edgeChains.sId.data();
offsetPFirst = 0;
offsetPSecond = 0;
unsigned int indexInCors = 0;
unsigned int numOfEdges = 0;
for ( unsigned int edgeId = 0; edgeId < offsetPS; edgeId++ )
{
//step1, put the first and second parts edge coordinates together from edge start to edge end
psId[numOfEdges++] = indexInCors;
indexInArray = pFirstPartEdgeS_[edgeId];
offsetPFirst = pFirstPartEdgeS_[edgeId + 1];
for ( tempID = offsetPFirst - 1; tempID >= indexInArray; tempID-- )
{ //add first part edge
pxCors[indexInCors] = pFirstPartEdgeX_[tempID];
pyCors[indexInCors++] = pFirstPartEdgeY_[tempID];
}
indexInArray = pSecondPartEdgeS_[edgeId];
offsetPSecond = pSecondPartEdgeS_[edgeId + 1];
for ( tempID = indexInArray + 1; tempID < (int) offsetPSecond; tempID++ )
{ //add second part edge
pxCors[indexInCors] = pSecondPartEdgeX_[tempID];
pyCors[indexInCors++] = pSecondPartEdgeY_[tempID];
}
}
psId[numOfEdges] = indexInCors; //the end index of the last edge
edgeChains.numOfEdges = numOfEdges;
return 1;
}
int EDLineDetector::EDline( cv::Mat &image, LineChains &lines )
{
//first, call EdgeDrawing function to extract edges
EdgeChains edges;
if( ( EdgeDrawing( image, edges ) ) != 1 )
{
cout << "Line Detection not finished" << endl;
return -1;
}
//detect lines
unsigned int linePixelID = edges.sId[edges.numOfEdges];
lines.xCors.resize( linePixelID );
lines.yCors.resize( linePixelID );
lines.sId.resize( 5 * edges.numOfEdges );
unsigned int *pEdgeXCors = edges.xCors.data();
unsigned int *pEdgeYCors = edges.yCors.data();
unsigned int *pEdgeSID = edges.sId.data();
unsigned int *pLineXCors = lines.xCors.data();
unsigned int *pLineYCors = lines.yCors.data();
unsigned int *pLineSID = lines.sId.data();
logNT_ = 2.0 * ( log10( (double) imageWidth ) + log10( (double) imageHeight ) );
double lineFitErr = 0; //the line fit error;
std::vector<double> lineEquation( 2, 0 );
lineEquations_.clear();
lineEndpoints_.clear();
lineDirection_.clear();
unsigned char *pdirImg = dirImg_.data;
unsigned int numOfLines = 0;
unsigned int newOffsetS = 0;
unsigned int offsetInEdgeArrayS, offsetInEdgeArrayE; //start index and end index
unsigned int offsetInLineArray = 0;
float direction; //line direction
for ( unsigned int edgeID = 0; edgeID < edges.numOfEdges; edgeID++ )
{
offsetInEdgeArrayS = pEdgeSID[edgeID];
offsetInEdgeArrayE = pEdgeSID[edgeID + 1];
while ( offsetInEdgeArrayE > offsetInEdgeArrayS + minLineLen_ )
{ //extract line segments from an edge, may find more than one segments
//find an initial line segment
while ( offsetInEdgeArrayE > offsetInEdgeArrayS + minLineLen_ )
{
lineFitErr = LeastSquaresLineFit_( pEdgeXCors, pEdgeYCors, offsetInEdgeArrayS, lineEquation );
if( lineFitErr <= lineFitErrThreshold_ )
break; //ok, an initial line segment detected
offsetInEdgeArrayS += SkipEdgePoint; //skip the first two pixel in the chain and try with the remaining pixels
}
if( lineFitErr > lineFitErrThreshold_ )
break; //no line is detected
//An initial line segment is detected. Try to extend this line segment
pLineSID[numOfLines] = offsetInLineArray;
double coef1 = 0; //for a line ax+by+c=0, coef1 = 1/sqrt(a^2+b^2);
double pointToLineDis; //for a line ax+by+c=0 and a point(xi, yi), pointToLineDis = coef1*|a*xi+b*yi+c|
bool bExtended = true;
bool bFirstTry = true;
int numOfOutlier; //to against noise, we accept a few outlier of a line.
int tryTimes = 0;
if( pdirImg[pEdgeYCors[offsetInEdgeArrayS] * imageWidth + pEdgeXCors[offsetInEdgeArrayS]] == Horizontal )
{ //y=ax+b, i.e. ax-y+b=0
while ( bExtended )
{
tryTimes++;
if( bFirstTry )
{
bFirstTry = false;
for ( int i = 0; i < minLineLen_; i++ )
{ //First add the initial line segment to the line array
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
}
}
else
{ //after each try, line is extended, line equation should be re-estimated
//adjust the line equation
lineFitErr = LeastSquaresLineFit_( pLineXCors, pLineYCors, pLineSID[numOfLines], newOffsetS, offsetInLineArray, lineEquation );
}
coef1 = 1 / sqrt( lineEquation[0] * lineEquation[0] + 1 );
numOfOutlier = 0;
newOffsetS = offsetInLineArray;
while ( offsetInEdgeArrayE > offsetInEdgeArrayS )
{
pointToLineDis = fabs( lineEquation[0] * pEdgeXCors[offsetInEdgeArrayS] - pEdgeYCors[offsetInEdgeArrayS] + lineEquation[1] ) * coef1;
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
if( pointToLineDis > lineFitErrThreshold_ )
{
numOfOutlier++;
if( numOfOutlier > 3 )
break;
}
else
{ //we count number of connective outliers.
numOfOutlier = 0;
}
}
//pop back the last few outliers from lines and return them to edge chain
offsetInLineArray -= numOfOutlier;
offsetInEdgeArrayS -= numOfOutlier;
if( offsetInLineArray - newOffsetS > 0 && tryTimes < TryTime )
{ //some new pixels are added to the line
}
else
{
bExtended = false; //no new pixels are added.
}
}
//the line equation coefficients,for line w1x+w2y+w3 =0, we normalize it to make w1^2+w2^2 = 1.
std::vector<double> lineEqu( 3, 0 );
lineEqu[0] = lineEquation[0] * coef1;
lineEqu[1] = -1 * coef1;
lineEqu[2] = lineEquation[1] * coef1;
if( LineValidation_( pLineXCors, pLineYCors, pLineSID[numOfLines], offsetInLineArray, lineEqu, direction ) )
{ //check the line
//store the line equation coefficients
lineEquations_.push_back( lineEqu );
/*At last, compute the line endpoints and store them.
*we project the first and last pixels in the pixelChain onto the best fit line
*to get the line endpoints.
*xp= (w2^2*x0-w1*w2*y0-w3*w1)/(w1^2+w2^2)
*yp= (w1^2*y0-w1*w2*x0-w3*w2)/(w1^2+w2^2) */
std::vector<float> lineEndP( 4, 0 ); //line endpoints
double a1 = lineEqu[1] * lineEqu[1];
double a2 = lineEqu[0] * lineEqu[0];
double a3 = lineEqu[0] * lineEqu[1];
double a4 = lineEqu[2] * lineEqu[0];
double a5 = lineEqu[2] * lineEqu[1];
unsigned int Px = pLineXCors[pLineSID[numOfLines]]; //first pixel
unsigned int Py = pLineYCors[pLineSID[numOfLines]];
lineEndP[0] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[1] = (float) ( a2 * Py - a3 * Px - a5 ); //y
Px = pLineXCors[offsetInLineArray - 1]; //last pixel
Py = pLineYCors[offsetInLineArray - 1];
lineEndP[2] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[3] = (float) ( a2 * Py - a3 * Px - a5 ); //y
lineEndpoints_.push_back( lineEndP );
lineDirection_.push_back( direction );
numOfLines++;
}
else
{
offsetInLineArray = pLineSID[numOfLines]; // line was not accepted, the offset is set back
}
}
else
{ //x=ay+b, i.e. x-ay-b=0
while ( bExtended )
{
tryTimes++;
if( bFirstTry )
{
bFirstTry = false;
for ( int i = 0; i < minLineLen_; i++ )
{ //First add the initial line segment to the line array
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
}
}
else
{ //after each try, line is extended, line equation should be re-estimated
//adjust the line equation
lineFitErr = LeastSquaresLineFit_( pLineXCors, pLineYCors, pLineSID[numOfLines], newOffsetS, offsetInLineArray, lineEquation );
}
coef1 = 1 / sqrt( 1 + lineEquation[0] * lineEquation[0] );
numOfOutlier = 0;
newOffsetS = offsetInLineArray;
while ( offsetInEdgeArrayE > offsetInEdgeArrayS )
{
pointToLineDis = fabs( pEdgeXCors[offsetInEdgeArrayS] - lineEquation[0] * pEdgeYCors[offsetInEdgeArrayS] - lineEquation[1] ) * coef1;
pLineXCors[offsetInLineArray] = pEdgeXCors[offsetInEdgeArrayS];
pLineYCors[offsetInLineArray++] = pEdgeYCors[offsetInEdgeArrayS++];
if( pointToLineDis > lineFitErrThreshold_ )
{
numOfOutlier++;
if( numOfOutlier > 3 )
break;
}
else
{ //we count number of connective outliers.
numOfOutlier = 0;
}
}
//pop back the last few outliers from lines and return them to edge chain
offsetInLineArray -= numOfOutlier;
offsetInEdgeArrayS -= numOfOutlier;
if( offsetInLineArray - newOffsetS > 0 && tryTimes < TryTime )
{ //some new pixels are added to the line
}
else
{
bExtended = false; //no new pixels are added.
}
}
//the line equation coefficients,for line w1x+w2y+w3 =0, we normalize it to make w1^2+w2^2 = 1.
std::vector<double> lineEqu( 3, 0 );
lineEqu[0] = 1 * coef1;
lineEqu[1] = -lineEquation[0] * coef1;
lineEqu[2] = -lineEquation[1] * coef1;
if( LineValidation_( pLineXCors, pLineYCors, pLineSID[numOfLines], offsetInLineArray, lineEqu, direction ) )
{ //check the line
//store the line equation coefficients
lineEquations_.push_back( lineEqu );
/*At last, compute the line endpoints and store them.
*we project the first and last pixels in the pixelChain onto the best fit line
*to get the line endpoints.
*xp= (w2^2*x0-w1*w2*y0-w3*w1)/(w1^2+w2^2)
*yp= (w1^2*y0-w1*w2*x0-w3*w2)/(w1^2+w2^2) */
std::vector<float> lineEndP( 4, 0 ); //line endpoints
double a1 = lineEqu[1] * lineEqu[1];
double a2 = lineEqu[0] * lineEqu[0];
double a3 = lineEqu[0] * lineEqu[1];
double a4 = lineEqu[2] * lineEqu[0];
double a5 = lineEqu[2] * lineEqu[1];
unsigned int Px = pLineXCors[pLineSID[numOfLines]]; //first pixel
unsigned int Py = pLineYCors[pLineSID[numOfLines]];
lineEndP[0] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[1] = (float) ( a2 * Py - a3 * Px - a5 ); //y
Px = pLineXCors[offsetInLineArray - 1]; //last pixel
Py = pLineYCors[offsetInLineArray - 1];
lineEndP[2] = (float) ( a1 * Px - a3 * Py - a4 ); //x
lineEndP[3] = (float) ( a2 * Py - a3 * Px - a5 ); //y
lineEndpoints_.push_back( lineEndP );
lineDirection_.push_back( direction );
numOfLines++;
}
else
{
offsetInLineArray = pLineSID[numOfLines]; // line was not accepted, the offset is set back
}
}
//Extract line segments from the remaining pixel; Current chain has been shortened already.
}
} //end for(unsigned int edgeID=0; edgeID<edges.numOfEdges; edgeID++)
pLineSID[numOfLines] = offsetInLineArray;
lines.numOfLines = numOfLines;
return 1;
}
double EDLineDetector::LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, std::vector<double> &lineEquation )
{
float * pMatT;
float * pATA;
double fitError = 0;
double coef;
unsigned char *pdirImg = dirImg_.data;
unsigned int offset = offsetS;
/*If the first pixel in this chain is horizontal,
*then we try to find a horizontal line, y=ax+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Horizontal )
{
/*Build the system,and solve it using least square regression: mat * [a,b]^T = vec
* [x0,1] [y0]
* [x1,1] [a] [y1]
* . [b] = .
* [xn,1] [yn]*/
pMatT = fitMatT.ptr<float>(); //fitMatT = [x0, x1, ... xn; 1,1,...,1];
for ( int i = 0; i < minLineLen_; i++ )
{
//*(pMatT+minLineLen_) = 1; //the value are not changed;
* ( pMatT++ ) = (float) xCors[offsetS];
fitVec[0][i] = (float) yCors[offsetS++];
}
ATA = fitMatT * fitMatT.t();
ATV = fitMatT * fitVec.t();
/* [a,b]^T = Inv(mat^T * mat) * mat^T * vec */
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
// lineEquation = svd.Invert(ATA) * matT * vec;
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
/*compute line fit error */
for ( int i = 0; i < minLineLen_; i++ )
{
//coef = double(yCors[offset]) - double(xCors[offset++]) * lineEquation[0] - lineEquation[1];
coef = double( yCors[offset] ) - double( xCors[offset] ) * lineEquation[0] - lineEquation[1];
offset++;
fitError += coef * coef;
}
return sqrt( fitError );
}
/*If the first pixel in this chain is vertical,
*then we try to find a vertical line, x=ay+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Vertical )
{
/*Build the system,and solve it using least square regression: mat * [a,b]^T = vec
* [y0,1] [x0]
* [y1,1] [a] [x1]
* . [b] = .
* [yn,1] [xn]*/
pMatT = fitMatT.ptr<float>(); //fitMatT = [y0, y1, ... yn; 1,1,...,1];
for ( int i = 0; i < minLineLen_; i++ )
{
//*(pMatT+minLineLen_) = 1;//the value are not changed;
* ( pMatT++ ) = (float) yCors[offsetS];
fitVec[0][i] = (float) xCors[offsetS++];
}
ATA = fitMatT * ( fitMatT.t() );
ATV = fitMatT * fitVec.t();
/* [a,b]^T = Inv(mat^T * mat) * mat^T * vec */
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
// lineEquation = svd.Invert(ATA) * matT * vec;
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
/*compute line fit error */
for ( int i = 0; i < minLineLen_; i++ )
{
//coef = double(xCors[offset]) - double(yCors[offset++]) * lineEquation[0] - lineEquation[1];
coef = double( xCors[offset] ) - double( yCors[offset] ) * lineEquation[0] - lineEquation[1];
offset++;
fitError += coef * coef;
}
return sqrt( fitError );
}
return 0;
}
double EDLineDetector::LeastSquaresLineFit_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int newOffsetS,
unsigned int offsetE, std::vector<double> &lineEquation )
{
int length = offsetE - offsetS;
int newLength = offsetE - newOffsetS;
if( length <= 0 || newLength <= 0 )
{
cout << "EDLineDetector::LeastSquaresLineFit_ Error:"
" the expected line index is wrong...offsetE = "
<< offsetE << ", offsetS=" << offsetS << ", newOffsetS=" << newOffsetS << endl;
return -1;
}
if( lineEquation.size() != 2 )
{
std::cout << "SHOULD NOT BE != 2" << std::endl;
}
cv::Mat_<float> matT( 2, newLength );
cv::Mat_<float> vec( newLength, 1 );
float * pMatT;
float * pATA;
double coef;
unsigned char *pdirImg = dirImg_.data;
/*If the first pixel in this chain is horizontal,
*then we try to find a horizontal line, y=ax+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Horizontal )
{
/*Build the new system,and solve it using least square regression: mat * [a,b]^T = vec
* [x0',1] [y0']
* [x1',1] [a] [y1']
* . [b] = .
* [xn',1] [yn']*/
pMatT = matT.ptr<float>(); //matT = [x0', x1', ... xn'; 1,1,...,1]
for ( int i = 0; i < newLength; i++ )
{
* ( pMatT + newLength ) = 1;
* ( pMatT++ ) = (float) xCors[newOffsetS];
vec[0][i] = (float) yCors[newOffsetS++];
}
/* [a,b]^T = Inv(ATA + mat^T * mat) * (ATV + mat^T * vec) */
tempMatLineFit = matT * matT.t();
tempVecLineFit = matT * vec;
ATA = ATA + tempMatLineFit;
ATV = ATV + tempVecLineFit;
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
return 0;
}
/*If the first pixel in this chain is vertical,
*then we try to find a vertical line, x=ay+b;*/
if( pdirImg[yCors[offsetS] * imageWidth + xCors[offsetS]] == Vertical )
{
/*Build the system,and solve it using least square regression: mat * [a,b]^T = vec
* [y0',1] [x0']
* [y1',1] [a] [x1']
* . [b] = .
* [yn',1] [xn']*/
pMatT = matT.ptr<float>(); //matT = [y0', y1', ... yn'; 1,1,...,1]
for ( int i = 0; i < newLength; i++ )
{
* ( pMatT + newLength ) = 1;
* ( pMatT++ ) = (float) yCors[newOffsetS];
vec[0][i] = (float) xCors[newOffsetS++];
}
/* [a,b]^T = Inv(ATA + mat^T * mat) * (ATV + mat^T * vec) */
// matT.MultiplyWithTransposeOf(matT, tempMatLineFit);
tempMatLineFit = matT * matT.t();
tempVecLineFit = matT * vec;
ATA = ATA + tempMatLineFit;
ATV = ATV + tempVecLineFit;
// pATA = ATA.GetData();
pATA = ATA.ptr<float>();
coef = 1.0 / ( double( pATA[0] ) * double( pATA[3] ) - double( pATA[1] ) * double( pATA[2] ) );
lineEquation[0] = coef * ( double( pATA[3] ) * double( ATV[0][0] ) - double( pATA[1] ) * double( ATV[0][1] ) );
lineEquation[1] = coef * ( double( pATA[0] ) * double( ATV[0][1] ) - double( pATA[2] ) * double( ATV[0][0] ) );
}
return 0;
}
bool EDLineDetector::LineValidation_( unsigned int *xCors, unsigned int *yCors, unsigned int offsetS, unsigned int offsetE,
std::vector<double> &lineEquation, float &direction )
{
if( bValidate_ )
{
int n = offsetE - offsetS;
/*first compute the direction of line, make sure that the dark side always be the
*left side of a line.*/
int meanGradientX = 0, meanGradientY = 0;
short *pdxImg = dxImg_.ptr<short>();
short *pdyImg = dyImg_.ptr<short>();
double dx, dy;
std::vector<double> pointDirection;
int index;
for ( int i = 0; i < n; i++ )
{
index = yCors[offsetS] * imageWidth + xCors[offsetS];
offsetS++;
meanGradientX += pdxImg[index];
meanGradientY += pdyImg[index];
dx = (double) pdxImg[index];
dy = (double) pdyImg[index];
pointDirection.push_back( atan2( -dx, dy ) );
}
dx = fabs( lineEquation[1] );
dy = fabs( lineEquation[0] );
if( meanGradientX == 0 && meanGradientY == 0 )
{ //not possible, if happens, it must be a wrong line,
return false;
}
if( meanGradientX > 0 && meanGradientY >= 0 )
{ //first quadrant, and positive direction of X axis.
direction = (float) atan2( -dy, dx ); //line direction is in fourth quadrant
}
if( meanGradientX <= 0 && meanGradientY > 0 )
{ //second quadrant, and positive direction of Y axis.
direction = (float) atan2( dy, dx ); //line direction is in first quadrant
}
if( meanGradientX < 0 && meanGradientY <= 0 )
{ //third quadrant, and negative direction of X axis.
direction = (float) atan2( dy, -dx ); //line direction is in second quadrant
}
if( meanGradientX >= 0 && meanGradientY < 0 )
{ //fourth quadrant, and negative direction of Y axis.
direction = (float) atan2( -dy, -dx ); //line direction is in third quadrant
}
/*then check whether the line is on the border of the image. We don't keep the border line.*/
if( fabs( direction ) < 0.15 || M_PI - fabs( direction ) < 0.15 )
{ //Horizontal line
if( fabs( lineEquation[2] ) < 10 || fabs( imageHeight - fabs( lineEquation[2] ) ) < 10 )
{ //upper border or lower border
return false;
}
}
if( fabs( fabs( direction ) - M_PI * 0.5 ) < 0.15 )
{ //Vertical line
if( fabs( lineEquation[2] ) < 10 || fabs( imageWidth - fabs( lineEquation[2] ) ) < 10 )
{ //left border or right border
return false;
}
}
//count the aligned points on the line which have the same direction as the line.
double disDirection;
int k = 0;
for ( int i = 0; i < n; i++ )
{
disDirection = fabs( direction - pointDirection[i] );
if( fabs( 2 * M_PI - disDirection ) < 0.392699 || disDirection < 0.392699 )
{ //same direction, pi/8 = 0.392699081698724
k++;
}
}
//now compute NFA(Number of False Alarms)
double ret = nfa( n, k, 0.125, logNT_ );
return ( ret > 0 ); //0 corresponds to 1 mean false alarm
}
else
{
return true;
}
}
int EDLineDetector::EDline( cv::Mat &image )
{
if( ( EDline( image, lines_/*, smoothed*/ ) ) != 1 )
{
return -1;
}
lineSalience_.clear();
lineSalience_.resize( lines_.numOfLines );
unsigned char *pgImg = gImgWO_.ptr();
unsigned int indexInLineArray;
unsigned int *pXCor = lines_.xCors.data();
unsigned int *pYCor = lines_.yCors.data();
unsigned int *pSID = lines_.sId.data();
for ( unsigned int i = 0; i < lineSalience_.size(); i++ )
{
int salience = 0;
for ( indexInLineArray = pSID[i]; indexInLineArray < pSID[i + 1]; indexInLineArray++ )
{
salience += pgImg[pYCor[indexInLineArray] * imageWidth + pXCor[indexInLineArray]];
}
lineSalience_[i] = (float) salience;
}
return 1;
}
modules/line_descriptor/src/mihasher.cpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Mohammad Norouzi, Ali Punjani, David J. Fleet,
// all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#include "precomp.hpp"
/* execute a batch query */
void Mihasher::batchquery( UINT32 * results, UINT32 *numres, const cv::Mat & queries, UINT32 numq, int dim1queries )
{
/* create and initialize a bitarray */
counter = new bitarray;
counter->init( N );
UINT32 *res = new UINT32[K * ( D + 1 )];
UINT64 *chunks = new UINT64[m];
UINT32 * presults = results;
UINT32 *pnumres = numres;
/* make a copy of input queries */
cv::Mat queries_clone = queries.clone();
/* set a pointer to first query (row) */
UINT8 *pq = queries_clone.ptr();
/* loop over number of descriptors */
for ( size_t i = 0; i < numq; i++ )
{
/* for every descriptor, query database */
query( presults, pnumres, pq, chunks, res );
/* move pointer to write next K indeces */
presults += K;
pnumres += B + 1;
/* move forward pointer to current row in descriptors matrix */
pq += dim1queries;
}
delete[] res;
delete[] chunks;
delete counter;
}
/* execute a single query */
void Mihasher::query( UINT32* results, UINT32* numres, UINT8 * Query, UINT64 *chunks, UINT32 *res )
{
/* if K == 0 that means we want everything to be processed.
So maxres = N in that case. Otherwise K limits the results processed */
UINT32 maxres = K ? K : (UINT32) N;
/* number of results so far obtained (up to a distance of s per chunk) */
UINT32 n = 0;
/* number of candidates tested with full codes (not counting duplicates) */
UINT32 nc = 0;
/* counting everything retrieved (duplicates are counted multiple times)
number of lookups (and xors) */
UINT32 nl = 0;
UINT32 nd = 0;
UINT32 *arr;
int size = 0;
UINT32 index;
int hammd;
counter->erase();
memset( numres, 0, ( B + 1 ) * sizeof ( *numres ) );
split( chunks, Query, m, mplus, b );
/* the growing search radius per substring */
int s;
/* current b: for the first mplus substrings it is b, for the rest it is (b-1) */
int curb = b;
for ( s = 0; s <= d && n < maxres; s++ )
{
for ( int k = 0; k < m; k++ )
{
if( k < mplus )
curb = b;
else
curb = b - 1;
UINT64 chunksk = chunks[k];
/* number of bit-strings with s number of 1s */
nl += xornum[s + 1] - xornum[s];
/* the bit-string with s number of 1s */
UINT64 bitstr = 0;
for ( int i = 0; i < s; i++ )
/* power[i] stores the location of the i'th 1 */
power[i] = i;
/* used for stopping criterion (location of (s+1)th 1) */
power[s] = curb + 1;
/* bit determines the 1 that should be moving to the left */
int bit = s - 1;
/* start from the left-most 1, and move it to the left until
it touches another one */
/* the loop for changing bitstr */
bool infiniteWhile = true;
while ( infiniteWhile )
{
if( bit != -1 )
{
bitstr ^= ( power[bit] == bit ) ? (UINT64) 1 << power[bit] : (UINT64) 3 << ( power[bit] - 1 );
power[bit]++;
bit--;
}
else
{ /* bit == -1 */
/* the binary code bitstr is available for processing */
arr = H[k].query( chunksk ^ bitstr, &size ); // lookup
if( size )
{ /* the corresponding bucket is not empty */
nd += size;
for ( int c = 0; c < size; c++ )
{
index = arr[c];
if( !counter->get( index ) )
{ /* if it is not a duplicate */
counter->set( index );
hammd = match( codes.ptr() + (UINT64) index * ( B_over_8 ), Query, B_over_8 );
nc++;
if( hammd <= D && numres[hammd] < maxres )
res[hammd * K + numres[hammd]] = index + 1;
numres[hammd]++;
}
}
}
/* end of processing */
while ( ++bit < s && power[bit] == power[bit + 1] - 1 )
{
bitstr ^= (UINT64) 1 << ( power[bit] - 1 );
power[bit] = bit;
}
if( bit == s )
break;
}
}
n = n + numres[s * m + k];
if( n >= maxres )
break;
}
}
n = 0;
for ( s = 0; s <= D && (int) n < K; s++ )
{
for ( int c = 0; c < (int) numres[s] && (int) n < K; c++ )
results[n++] = res[s * K + c];
}
}
/* constructor 2 */
Mihasher::Mihasher( int _B, int _m )
{
B = _B;
B_over_8 = B / 8;
m = _m;
b = (int) ceil( (double) B / m );
/* assuming that B/2 is large enough radius to include
all of the k nearest neighbors */
D = (int) ceil( B / 2.0 );
d = (int) ceil( (double) D / m );
/* mplus is the number of chunks with b bits
(m-mplus) is the number of chunks with (b-1) bits */
mplus = B - m * ( b - 1 );
xornum = new UINT32[d + 2];
xornum[0] = 0;
for ( int i = 0; i <= d; i++ )
xornum[i + 1] = xornum[i] + (UINT32) choose( b, i );
H = new SparseHashtable[m];
/* H[i].init might fail */
for ( int i = 0; i < mplus; i++ )
H[i].init( b );
for ( int i = mplus; i < m; i++ )
H[i].init( b - 1 );
}
/* K setter */
void Mihasher::setK( int _K )
{
K = _K;
}
/* desctructor */
Mihasher::~Mihasher()
{
delete[] xornum;
delete[] H;
}
/* populate tables */
void Mihasher::populate( cv::Mat & _codes, UINT32 _N, int dim1codes )
{
N = _N;
codes = _codes;
UINT64 * chunks = new UINT64[m];
UINT8 * pcodes = codes.ptr();
for ( UINT64 i = 0; i < N; i++, pcodes += dim1codes )
{
split( chunks, pcodes, m, mplus, b );
for ( int k = 0; k < m; k++ )
H[k].insert( chunks[k], (UINT32) i );
if( i % (int) ceil( N / 1000.0 ) == 0 )
fflush( stdout );
}
delete[] chunks;
}
modules/line_descriptor/src/precomp.hpp
View file @ 9a52f462
...@@ -68,6 +68,10 @@ ...@@ -68,6 +68,10 @@
#include <sstream> #include <sstream>
#include <vector> #include <vector>
#include "bitarray.hpp"
#include "bitops.hpp"
#include "types.hpp"
#include "opencv2/line_descriptor.hpp" #include "opencv2/line_descriptor.hpp"
#endif #endif
modules/line_descriptor/src/sparse_hashtable.cpp deleted 100644 → 0
View file @ e39065bf
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2014, Biagio Montesano, all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#include "precomp.hpp"
const int SparseHashtable::MAX_B = 37;
/* constructor */
SparseHashtable::SparseHashtable()
{
table = NULL;
size = 0;
b = 0;
}
/* initializer */
int SparseHashtable::init( int _b )
{
b = _b;
if( b < 5 || b > MAX_B || b > (int) ( sizeof(UINT64) * 8 ) )
return 1;
size = UINT64_1 << ( b - 5 ); // size = 2 ^ b
table = (BucketGroup*) calloc( size, sizeof(BucketGroup) );
return 0;
}
/* destructor */
SparseHashtable::~SparseHashtable()
{
free( table );
}
/* insert data */
void SparseHashtable::insert( UINT64 index, UINT32 data )
{
table[index >> 5].insert( (int) ( index % 32 ), data );
}
/* query data */
UINT32* SparseHashtable::query( UINT64 index, int *Size )
{
return table[index >> 5].query( (int) ( index % 32 ), Size );
}
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment