## Tuesday, March 25, 2014

### Edge Detection with 1D Inplace Fast Haar Wavelet Transform Applied to Image Rows

Problem Statement

In this post, we will use the 1D inplace fast Haar wavelet transform to design a simple technique to detect edges in images. We will use OpenCV to do the coding.  We will start with an implementation of the 1D inplace fast Haar wavelet transform.

const double loge_of_2 = log(2);
const double log10_of_2 = log10(2);
double log2(double x) { return log10(x)/log10_of_2; }

void inPlaceFastHaarWaveletTransformWithNumSweeps(double* row_data, int size, int num_sweeps) {
int I = 1;
int GAP_SIZE = 2;
int NUM_SAMPLE_VALS = size;
if ( num_sweeps < 1 || num_sweeps > size ) return;
double a = 0;
double c = 0;
for(int SWEEP_NUM = 1; SWEEP_NUM <= num_sweeps; SWEEP_NUM++) {
NUM_SAMPLE_VALS /= 2;
for(int K = 0; K < NUM_SAMPLE_VALS; K++) {
a = (row_data[GAP_SIZE * K] + row_data[GAP_SIZE * K + I])/2;
c = (row_data[GAP_SIZE * K] - row_data[GAP_SIZE * K + I])/2;
row_data[GAP_SIZE * K] = a;
row_data[GAP_SIZE * K + I] = c;
}
I = GAP_SIZE;
GAP_SIZE *= 2;
}
}

void inPlaceFastHaarWaveletTransform(double* row_data, int n) {
if ( n == 0 || n % 2 != 0 ) return;
int num_sweeps = log(n)/loge_of_2;
inPlaceFastHaarWaveletTransformWithNumSweeps(row_data, n, num_sweeps);
}

We can test this code as follows.

void displayArray(double* ary, int size) {
for(int i = 0; i < size; i++) { std::cout << ary[i] << " "; }
std::cout << std::endl << std::flush;
}

double row_data[32] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 255, 255, 255, 255, 255};
inPlaceFastHaarWaveletTransform(row_data, 32);
displayArray(row_data, 32);

This code produces the following output.

39.8438 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -39.8438 0 0 0 0 0 0 0 -79.6875 0 -63.75 -127.5 -95.625 0 0 0

cv::Mat Manipulation

The next technical problem we need to solve is how to extract rows from cv::Mat structures. The above Haar transform can then be applied to each row. Edges can be detected by looking at the absolute magnitudes of Haar wavelets at specific frequencies (e.g., 1/2^{n-1}, where n is the length of the image row which has to be an integral power of 2). We will use grayscale images. cv::Mat mat(n, n, CV_8U, cv::Scalar(0));

If you work in RGB images, check OpenCV documentation for appropriate conversions. The cv::Mat structures consist of uchar rows and columns.  Hence, a simple function to convert a row of uchars into an array of doubles.

double* ucharAryToDoubleAry(uchar* row_data, int size) {
double* rslt = new double[size];
for(int i = 0; i < size; i++) {
rslt[i] = (double)row_data[i];
}
return rslt;
}

Actual Detection of Edges
The function, given below, that detects edges uses the highest frequency Haar wavelets computed in the current image row. The function goes through the image specified in the first argument row by row. It takes the wavelet bandwidth defined by the next two arguments: lower_theta and upper_theta. The fourth parameter, rv, specifies the value on the scale [0, 255] that we want a pixel to contain if that pixel contains an edge. If the value of the wavelet is negative, then the current pixel contains an edge. If the value of the wavelet is positive, the previous pixel contains an edge.

void haarDetectHorEdges(cv::Mat &img, int lower_theta, int upper_theta, int rv)
{
const int nr = img.rows;
const int nc = img.cols * img.channels();
uchar* row_data;
double* converted_row_data;
double hd = 0.0;

for(int r = 0; r < nr; r++) {
row_data = img.ptr<uchar>(r);
converted_row_data = ucharAryToDoubleAry(row_data, nc);
inPlaceFastHaarWaveletTransform(converted_row_data, nc);
for(int c = 0; c < nc; c += 2) {
hd = converted_row_data[c];
if ( std::abs(hd) > lower_theta && std::abs(hd) < upper_theta ) {
if ( hd < 0 ) {
img.at<uchar>(r, c) = (uchar)rv;
}
else if ( hd > 0 && c > 0 ) {
img.at<uchar>(r, c-1) = (uchar)rv;
}
}
}
delete [] converted_row_data;
}
}

Tests
Now on to a few simple tools in OpenCV that we can use in our tests. One thing to keep in mind is that the length and width of each image must be equal to an integral power of 2. Otherwise, the Haar transform will not work.

cv::Mat createGrayscaleSquareMat(int n) {
cv::Mat mat(n, n, CV_8U, cv::Scalar(0));
return mat;
}

cv::Mat createGrayscaleSquareMatWithLine(int n, int y1, int x1, int y2, int x2) {
cv::Mat mat(n, n, CV_8U, cv::Scalar(0));
cv::Point p1(y1, x1);
cv::Point p2(y2, x2);
cv::line(mat, p1, p2, cv::Scalar(255), 1, 8);
return mat;
}

void drawCircle(cv::Mat &img, int center_y1, int center_x1, int rad, int thickness) {
circle(img,
cv::Point(center_y1, center_x1),
cv::Scalar(255),
thickness,
8);
}

void drawEllipse(cv::Mat &img, double angle)
{
int thickness = 2;
int lineType = 8;
const int nr = img.rows;
cv::ellipse(img,
cv::Point((int)nr/2, (int)nr/2),
cv::Size((int)nr/4, (int)nr/16),
angle,
0,
360,
cv::Scalar(255),
thickness,
lineType);
}

The first test, test01(), creates a grayscale 256 x 256 image and draws a line in it. The created image is shown in Figure 1 and is displayed with cv::imshow()

void test01() {
cv::Mat img = createGrayscaleSquareMatWithLine(256, 10, 10, 120, 120);
cv::imshow("orig", img);
haarDetectHorEdges(img, 125, 128, 255);
cv::imshow("edge", img);
cv::waitKey(0);
}

 Figure 1. A grayscale image with a line

When the function haarDetectHorEdges is applied, the resulting image is shown in Figure 2.

 Figure 2. Edges detected in Figure 1
The second test, test02(), draws three white circles in a 256 x 256 grayscale image, as shown in Figure 3, and detects edges, as shown in Figure 4.

void test02() {
cv::Mat img = createGrayscaleSquareMat(256);
drawCircle(img, 170, 170, 50, 5);
drawCircle(img, 120, 120, 50, 5);
drawCircle(img, 170, 120, 50, 5);
cv::imshow("orig", img);
haarDetectHorEdges(img, 125, 128, 255);
cv::imshow("edge", img);
cv::waitKey(0);
}

 Figure 3. Three white circles in a grayscale image

 Figure 4. Edges detected in Figure 3

The third test, test03(), draws four ellipses circles in a 256 x 256 grayscale image, as shown in Figure 5, and detects edges, as shown in Figure 6.

void test03() {
cv::Mat img = cv::Mat::zeros(256, 256, CV_8U);
drawEllipse(img, 90);
drawEllipse(img, 0);
drawEllipse(img, 45);
drawEllipse(img, -45);
cv::imshow("orig", img);
haarDetectHorEdges(img, 125, 128, 255);
cv::imshow("edge", img);
cv::waitKey(0);
}

 Figure 5. Four ellipses in a grayscale image
 Figure 6. Edges detected in Figure 5 One final note. We have deliberately simplified our images to show how this method works. Grayscale images with a greater variety of pixel values will require more sophisticated thresholds and wavelet frequencies.