|
#include <vector>
#include <queue>
#include <iostream>
#include <math.h>
#include <opencv/cv.h>
using namespace std;
#ifndef __CANNYEDGE_H__
#define __CANNYEDGE_H__
class AdaptiveCannyEdge {
typedef double Precision;
typedef vector<Precision> dtVect;
typedef vector<dtVect> dtMat;
typedef pair<int,int> Position;
typedef vector< vector<unsigned char> > lgMat;
private:
static const Precision _eps;
static const Precision PercentOfPixelsNotEdges; //Magic number
static const Precision ThresholdRatio; //Magic number
static const Precision gaussKernel[16]; //1-D Gaussian Kernel
static const Precision derivGaussKernel[16]; //1-D Derivative of Gaussian Kernel
public:
// If you do not specify THRESH, CANNYEDGE chooses low and high values automatically.
static bool CannyEdge( const IplImage* src, IplImage* & edgeImg );
static bool CannyGradient( const IplImage* src, IplImage* & gradImg );
private:
// Transform to a double precision intensity image if necessary
static bool im2single( const IplImage* src, dtMat& mat );
// Compute smoothed numerical gradient of image 'src' along x (horizontal) or y (vertical)
// direction.
static bool imfilter( const dtMat& src, const dtVect& filter, dtMat& _mat, bool isHori = true );
// Calculate gradients using a derivative of Gaussian filter
static bool smoothGradient( const dtMat& src, dtMat& dx, dtMat& dy );
// Calculate Magnitude of Gradient
static void hypot( const dtMat& dx, const dtMat& dy, dtMat& magGrad );
// Determine Hysteresis Thresholds
static void selectThresholds( const dtMat& src, double& lowThresh, double& highThresh );
// Perform Non-Maximum Suppression Thining and Hysteresis Thresholding of Edge Strength
static bool cannyFindLocalMaxima( const dtMat& src, const dtMat& dx, const dtMat& dy, const Position& pos );
// This function helps with the non-maximum suppression in the Canny edge detector.
static void thinAndThreshold( const dtMat& src, const dtMat& dx, const dtMat& dy,
const double& lowThresh, const double& highThresh, lgMat& dst );
// Return a binary image 'idxWeak' containing the objects that overlap the pixel (r,c), which
// (r,c) represents each location of 'idxStrongPts'.
static void bwselect( lgMat& idxWeak, lgMat& idxWeakNeg, const vector<Position>& idxStrongPts );
//Transform dtMat to an Image
static void dtMat2Img( const dtMat& mat, IplImage* & img );
// Transform lgMat to an Image
static void lgMat2Img( const lgMat& mat, IplImage* & img );
};
#endif
----------------------------------------------------------------------------------
#include "cannyEdge.h"
const AdaptiveCannyEdge::Precision AdaptiveCannyEdge::_eps = 1.0e-10;
const AdaptiveCannyEdge::Precision AdaptiveCannyEdge::PercentOfPixelsNotEdges = 0.7;
const AdaptiveCannyEdge::Precision AdaptiveCannyEdge::ThresholdRatio = 0.4;
const AdaptiveCannyEdge::Precision AdaptiveCannyEdge::gaussKernel[16] = {
0.000000220358050, 0.000007297256405, 0.000146569312970, 0.001785579770079,
0.013193749090229, 0.059130281094460, 0.160732768610747, 0.265003534507060,
0.265003534507060, 0.160732768610747, 0.059130281094460, 0.013193749090229,
0.001785579770079, 0.000146569312970, 0.000007297256405, 0.000000220358050
};
const AdaptiveCannyEdge::Precision AdaptiveCannyEdge::derivGaussKernel[16] = {
0.000026704586264, 0.000276122963398, 0.003355163265098, 0.024616683775044,
0.108194751875585, 0.278368310241814, 0.388430056419619, 0.196732206873178,
-0.196732206873178, -0.388430056419619, -0.278368310241814, -0.108194751875585,
-0.024616683775044, -0.003355163265098, -0.000276122963398, -0.000026704586264
};
bool AdaptiveCannyEdge::CannyEdge( const IplImage* src, IplImage* & edgeImg )
{
if( 0 == src ) {
return false;
}
dtMat mat;
if( !im2single(src, mat) ) {
return false;
}
dtMat dx, dy;
if( !smoothGradient( mat, dx, dy ) ) {
return false;
}
dtMat mag;
hypot( dx, dy, mag );
double _low=0.0, _high=0.0;
selectThresholds( mag, _low, _high );
lgMat dst;
thinAndThreshold( mag, dx, dy, _low, _high, dst );
lgMat2Img( dst, edgeImg );
return true;
}
bool AdaptiveCannyEdge::CannyGradient( const IplImage* src, IplImage* & gradImg )
{
if( 0 == src ) {
return false;
}
dtMat mat;
if( !im2single(src, mat) ) {
return false;
}
dtMat dx, dy;
if( !smoothGradient( mat, dx, dy ) ) {
return false;
}
dtMat mag;
hypot( dx, dy, mag );
dtMat2Img( mag, gradImg );
return true;
}
bool AdaptiveCannyEdge::im2single( const IplImage* src, dtMat& mat )
{
if( src->nChannels != 1 ) {
return false;
}
Precision max_range = 1.0;
Precision half_range = 0.0;
if( src->depth == IPL_DEPTH_8U ) {
max_range = 255.0;
} else if ( src->depth == IPL_DEPTH_8S ) {
max_range = 255.0;
half_range = 128.0;
} else if( src->depth == IPL_DEPTH_16U ) {
max_range = 65535.0;
} else if (src->depth == IPL_DEPTH_16S ) {
max_range = 65535.0;
half_range = 32768.0;
}
int m = src->height;
int n = src->width;
mat.resize( m, dtVect(n,0.0) );
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
mat[i][j] = (Precision)( ( (uchar*)src->imageData i * src->widthStep )[j] );
mat[i][j] = ( mat[i][j] half_range ) / ( max_range _eps );
}
}
return true;
}
bool AdaptiveCannyEdge::imfilter( const dtMat& src, const dtVect& filter, dtMat& _mat, bool isHori )
{
dtMat mat( src );
int m = mat.size();
int n = mat[0].size();
_mat.resize( m, dtVect(n, 0.0) );
dtVect t_vct;
if( isHori ) {
t_vct.resize( n 16, 0.0 );
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) t_vct[j 8] = mat[i][j];
for( int j=0; j<8; j ) t_vct[j] = mat[i][0];
for( int j=0; j<8; j ) t_vct[j n 8] = mat[i][n-1];
for( int j=0; j<n; j ) {
Precision t_sum = 0.0;
for( int k=0; k<16; k ) {
t_sum = t_vct[j k 1] * filter[15-k];
}
_mat[i][j] = t_sum;
}
}
} else {
t_vct.resize( m 16, 0.0 );
for( int i=0; i<n; i ) {
for( int j=0; j<m; j ) t_vct[j 8] = mat[j][i];
for( int j=0; j<8; j ) t_vct[j] = mat[0][i];
for( int j=0; j<8; j ) t_vct[j m 8] = mat[m-1][i];
for( int j=0; j<m; j ) {
Precision t_sum = 0.0;
for( int k=0; k<16; k ) {
t_sum = t_vct[j k 1] * filter[15-k];
}
_mat[j][i] = t_sum;
}
}
}
return true;
}
bool AdaptiveCannyEdge::smoothGradient( const dtMat& src, dtMat& dx, dtMat& dy )
{
int m = src.size();
if( m <= 0 ) {
return false;
}
int n = src[0].size();
dtMat mat( src );
dtVect filter1( gaussKernel, gaussKernel 16 );
dtVect filter2( derivGaussKernel, derivGaussKernel 16 );
imfilter( mat, filter1, dx, false );
imfilter( dx, filter2, dx, true );
imfilter( mat, filter1, dy, true );
imfilter( dy, filter2, dy, false );
return true;
}
void AdaptiveCannyEdge::hypot( const dtMat& dx, const dtMat& dy, dtMat& magGrad )
{
int m = dx.size();
int n = dx[0].size();
magGrad.resize( m, dtVect(n,0.0) );
Precision maxValue = 0.0;
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
Precision p1 = dx[i][j];
Precision p2 = dy[i][j];
Precision t = sqrt( p1*p1 p2*p2 );
if( maxValue < t ) maxValue = t;
magGrad[i][j] = t;
}
}
if( maxValue > _eps ) {
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
magGrad[i][j] /= maxValue;
}
}
}
}
void AdaptiveCannyEdge::selectThresholds( const dtMat& src, double& lowThresh, double& highThresh )
{
vector<int> counter(64, 0);
int m = src.size();
int n = src[0].size();
double scale = 63.0;
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
int z = floor( scale * src[i][j] 0.5 );
if( z >= 64 ) counter[63];
else if( z<0 ) counter[0];
else counter[ z ];
}
}
highThresh = 0.0;
double pixelsNotEdges = PercentOfPixelsNotEdges * m * n;
int cunsum = 0;
for( int i=0; i<64; i ) {
cunsum = counter[i];
if( cunsum - pixelsNotEdges > _eps ) {
highThresh = double(i 1) / 64.0;
break;
}
}
lowThresh = ThresholdRatio * highThresh;
}
bool AdaptiveCannyEdge::cannyFindLocalMaxima( const dtMat& src, const dtMat& dx, const dtMat& dy, const Position& pos )
{
int ix = pos.first;
int iy = pos.second;
int m = src.size();
int n = src[0].size();
if( ix <= 0 || ix >= m-1 || iy <= 0 || iy >= n-1 ) {
return false;
}
double dx_val = dx[ix][iy];
double dy_val = dy[ix][iy];
double gradmag1, gradmag2, gradmag = src[ix][iy];
bool case1 = false, case2 = false, case3 = false, case4 = false;
if( (dy_val<=0&&dx_val dy_val>0) || (dy_val>=0&&dx_val dy_val<0) )
{ //Case 1
double dam = dx_val >= 0 ? ( dx_val _eps ) : ( dx_val - _eps );
double d = fabs( dy_val / dam );
gradmag1 = src[ix][iy 1] * (1-d) src[ix-1][iy 1] * d;
gradmag2 = src[ix][iy-1] * (1-d) src[ix 1][iy-1] * d;
case1 = (gradmag>=gradmag1&&gradmag>=gradmag2);
}
if( (dx_val>0&&dx_val dy_val<=0) || (dx_val<0&&dx_val dy_val>=0) )
{ //Case 2
double dam = dy_val >= 0 ? ( dy_val _eps ) : ( dy_val - _eps );
double d = fabs( dx_val / dam );
gradmag1 = src[ix-1][iy] * (1-d) src[ix-1][iy 1] * d;
gradmag2 = src[ix 1][iy] * (1-d) src[ix 1][iy-1] * d;
case2 = (gradmag>=gradmag1&&gradmag>=gradmag2);
}
if( (dx_val<=0&&dx_val>dy_val) || (dx_val>=0&&dx_val<dy_val) )
{ //Case 3
double dam = dy_val >= 0 ? ( dy_val _eps ) : ( dy_val - _eps );
double d = fabs( dx_val / dam );
gradmag1 = src[ix-1][iy] * (1-d) src[ix-1][iy-1] * d;
gradmag2 = src[ix 1][iy] * (1-d) src[ix 1][iy 1] * d;
case3 = (gradmag>=gradmag1&&gradmag>=gradmag2);
}
if( (dy_val<0&&dx_val<=dy_val) || (dy_val>0&&dx_val>=dy_val ) )
{ //Case 4
double dam = dx_val >= 0 ? ( dx_val _eps ) : ( dx_val - _eps );
double d = fabs( dy_val / dam );
gradmag1 = src[ix][iy-1] * (1-d) src[ix-1][iy-1] * d;
gradmag2 = src[ix][iy 1] * (1-d) src[ix 1][iy 1] * d;
case4 = (gradmag>=gradmag1&&gradmag>=gradmag2);
}
return ( case1 || case2 || case3 || case4 );
}
void AdaptiveCannyEdge::thinAndThreshold( const dtMat& src, const dtMat& dx, const dtMat& dy,
const double& lowThresh, const double& highThresh, lgMat& dst )
{
const int m = src.size();
const int n = src[0].size();
lgMat idxWeakNeg( m, vector<unsigned char>(n, 1) );
lgMat idxWeak( m, vector<unsigned char>(n, 0) );
vector<Position> idxStrongPts;
idxStrongPts.reserve( ceil(m*n*PercentOfPixelsNotEdges) );
lgMat test( m, vector<unsigned char>(n, 0) );
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
if( cannyFindLocalMaxima( src, dx, dy, Position(i,j) ) ) {
if ( src[i][j] > lowThresh ) {
idxWeak[i][j] = 1;
idxWeakNeg[i][j] = 0;
}
if( src[i][j] > highThresh ) {
test[i][j] = 1;
idxStrongPts.push_back( Position(i,j) );
}
}
}
}
bwselect( idxWeak, idxWeakNeg, idxStrongPts );
dst.swap( idxWeak );
}
void AdaptiveCannyEdge::bwselect( lgMat& idxWeak, lgMat& idxWeakNeg, const vector<Position>& idxStrongPts )
{
const int m = idxWeak.size();
const int n = idxWeak[0].size();
const int ptSize = idxStrongPts.size();
const int d_x[8] = {0, 0,1,1, 1,-1,-1, -1};
const int d_y[8] = {1,-1,0,1,-1, 1, 0, -1};
for( int i=0; i<ptSize; i ) {
int ix = idxStrongPts[i].first;
int iy = idxStrongPts[i].second;
if( idxWeakNeg[ix][iy] == 0 ) {
queue<Position> Q;
Q.push( idxStrongPts[i] );
while( !Q.empty() ) {
Position pos = Q.front();
Q.pop();
ix = pos.first;
iy = pos.second;
idxWeakNeg[ix][iy] = 1;
for( int j=0; j<8; j ) {
if( 0 == idxWeakNeg[ ix d_x[j] ][ iy d_y[j] ] )
{
Position tpos( ix d_x[j], iy d_y[j] );
Q.push( tpos );
idxWeakNeg[ ix d_x[j] ][ iy d_y[j] ] = 1;
}
}
}
}
}
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
idxWeak[i][j] &= idxWeakNeg[i][j];
}
}
}
void AdaptiveCannyEdge::dtMat2Img( const dtMat& mat, IplImage* & img )
{
int m = mat.size();
int n = mat[0].size();
if( img != 0 ) {
try {
cvReleaseImage( &img );
} catch(...) {
//cvReleaseImage Error
}
img = 0;
}
img = cvCreateImage( cvSize(n,m), IPL_DEPTH_8U, 1 );
double maxValue = 1.0;
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
if( maxValue < fabs(mat[i][j]) ) maxValue = fabs(mat[i][j]);
}
}
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
( (uchar*)img->imageData i * img->widthStep )[j] =
(uchar)(int)( 255.0 * fabs( mat[i][j] ) / maxValue );
}
}
}
void AdaptiveCannyEdge::lgMat2Img( const lgMat& mat, IplImage* & img )
{
int m = mat.size();
int n = mat[0].size();
if( img != 0 ) {
try {
cvReleaseImage( &img );
} catch(...) {
//cvReleaseImage Error
}
img = 0;
}
img = cvCreateImage( cvSize(n,m), IPL_DEPTH_8U, 1 );
for( int i=0; i<m; i ) {
for( int j=0; j<n; j ) {
( (uchar*)img->imageData i * img->widthStep )[j] = ( mat[i][j] == 0 ? 0 : 255 );
}
}
}
|
