[ 영상 처리 ] Part3-2. OpenCV Diagonal Edge Detection & Image Pyramid(C++)

이번시간에는 Gaussian Filter를 약간 변형 해 보고 Salt and Pepper Noise에 적용하여 그 결과를 분석한다.

또한, Sobel Filter를 변형하여 대각 Edge를 추출하며, Gaussian&Laplacian Pyramid를 구현한다.

---

Before This Episode

https://jinwoo-jung.com/64

[ 영상 처리 ] Part3-1. OpenCV Edge Detection(C++)

---

Gaussian Filtering

9x9 Gaussian Filter를 생성한 후 Gaussian Filtering을 진행 해 보자.

MakeGaussianFilter

Mat MakeGaussianFilter(int32_t s32_Size) {
	Size size(s32_Size, s32_Size);
	Mat kn = Mat::zeros(size, CV_32FC1);
	double sigma = (s32_Size - 1) / 6.0;
	float* kn_data = (float*)kn.data;

	for (int c = 0; c < kn.cols; c++) {
		for (int r = 0; r < kn.rows; r++) {
			kn_data[r * kn.cols + c] = (float)gaussian2D((float)(c - kn.cols / 2), (float)(r - kn.rows / 2), sigma);
		}
	}

	return kn;
}

float gaussian2D(float c, float r, double sigma) {
	return exp(-(pow(c, 2) + pow(r, 2)) / (2 * pow(sigma, 2))) / (2 * CV_PI * pow(sigma, 2));
}

 

MakeGaussianFilter()는 Kernel Size를 인자로 받아서 2D Gaussian 수식에 근거하여 Kernel Weight를 계산하고, 생성한 Kernel을 반환하는 함수이다.

 

Kernel Size만큼의 빈 Mat 객체 kn을 선언한다. $\sigma$는 Kernel의 Width의 반이 $3*\sigma$가 되면 가장 좋기 때문에 위와 같이 구현하였다. 2중 for문을 돌면서 Kernel 중심으로부터 각 픽셀의 떨어진 거리를 인자로 하여 gaussian2D()로 부터 Gaussian Weight를 계산한다. gaussian2D()는 이전에도 사용한 2D Gaussian 수식에 근거한다.

 

9x9 Gaussian Kerenl

 

생성한 9x9 Gaussian Kernel은 위와 같다. 

 

myKernelConv()

void myKernelConv(const Mat& src_img, Mat& dst_img, const Mat& kn) {
	dst_img = Mat::zeros(src_img.size(), CV_8UC1);

	int wd = src_img.cols; int hg = src_img.rows;
	int kwd = kn.cols; int khg = kn.rows;
	int rad_w = kwd / 2; int rad_h = khg / 2;

	float* kn_data = (float*)kn.data;
	uchar* src_data = (uchar*)src_img.data;
	uchar* dst_data = (uchar*)dst_img.data;

	float wei, tmp, sum;

	for (int c = rad_w; c < wd - rad_w; c++) {
		for (int r = rad_h; r < hg - rad_h; r++) {
			tmp = 0.f;
			sum = 0.f;
			for (int kc = -rad_w; kc <= rad_w; kc++) {
				for (int kr = -rad_h; kr <= rad_h; kr++) {
					wei = (float)kn_data[(kr + rad_h) * kwd + (kc + rad_w)];
					tmp += wei * (float)src_data[(r + kr) * wd + (c + kc)];
					sum += wei;
				}
			}
			if (sum != 0.f) tmp = abs(tmp) / sum;
			else tmp = abs(tmp);

			if (tmp > 255.f) tmp = 255.f; 

			dst_data[r * wd + c] = (uchar)tmp;
		}

 

이전에 사용했던 코드와 동일하다. 간략하게 설명하자면 가장자리 문제를 고려하여 유효한 픽셀에 대하여 Gaussian Filtering을 진행한다. 인접 픽셀과 Weight의 곱을 다 더한 후 Weight의 합이 0이 아닐 경우 나눠 줌으로써 Normalization을 진행한다. 또한, 해당 결과가 255보다 클 경우 255로 설정하는데, 이는 Result Image의 Intensity Resolution이 8byte 이기 때문이다. 이후 계산된 Intensity를 Processing 하는 픽셀의 Intensity로 설정한다.

 

Original Image , myGaussianFilter , Opencv GaussianBlur

 

Gaussian Filtering을 진행한 결과 2번째 이미지와 같이 Smoothing 효과가 잘 일어남을 확인할 수 있다. 이는

cv::GaussianBlur()와 비교해 봤을 때 유사한 결과를 보이는 것을 확인할 수 있다.

 

Gaussian Filtering을 통해 Smoothing 효과가 일어남을 Histogram 분포를 통해 분석 해 보자.

GetHistogram()

Mat GetHistogram(Mat src) {
	Mat histogram;
	const int* channel_numbers = { 0 };
	float channel_range[] = { 0.0, 255.0 };
	const float* channel_ranges = channel_range;
	int number_bins = 256;

	// GrayScale  Ch 1, 히스토그램을 구할 체널 번호(Gray -> 0)
	calcHist(&src, 1, channel_numbers, Mat(), histogram, 1, &number_bins, &channel_ranges);

	int hist_w = src.cols;
	int hist_h = src.rows;
	int bin_w = cvRound((double)hist_w / number_bins);

	Mat histImage(hist_h, hist_w, CV_8UC1, Scalar(0, 0, 0));
	normalize(histogram, histogram, 0, histImage.rows, NORM_MINMAX, -1, Mat());

	for (int i = 1; i < number_bins; i++) {
		line(histImage, Point(bin_w * (i - 1), hist_h - cvRound(histogram.at<float>(i - 1))),
			Point(bin_w * (i), hist_h - cvRound(histogram.at<float>(i))),
			Scalar(255, 0, 0), 2, 8, 0);
	}

	return histImage;
}

 

코드에 대한 설명은 이전 포스팅으로 대체한다.

https://jinwoo-jung.com/48

[ 영상 처리 ] Part1-1. OpenCV Image Processing(C++)

 

Original Image Histogram                                                                     Gaussian Filtered Image Histogram

 

Gaussian Filtering 결과 Smoothing 효과가 발생한다. Original Image의 Histogram의 경우 왼쪽과 같이 Intensity가 낮은 곳이 높은 분포를 보임을 확인할 수 있다. Gaussian Filtering 결과 Intensity 분포가 여전히 낮은 곳에 높지만 상대적으로 0~255 전체 범위에 대하여 골고루 분포됨을 확인할 수 있다. 이는 Smoothing을 통해 인접 픽셀의 평균 값으로 픽셀 Intensity가 변하였기 때문이다.

 

Gaussian Filtering은 Noise 제거에 최적의 Filter라고 배웠다. 하지만 이는 Noise가 Gaussian 분포를 따를 때에 한정적이다. Salt and Pepper Noise 처럼 Outlier(튀는 값)이 존재하면 오히려 역효과를 낸다. 이러한 경우에는 Median Filter가 더 효과적임을 이전 포스팅에서 확인 하였다. Salt and Pepper Noise를 임의로 생성하여 Filtering 결과를 분석 해 보자.

 

MakeSaltandPepperNoise()

void MakeSaltandPepperNoise(const Mat& st_OriginalImage, Mat& st_ResultImage, int32_t s32_Num)
{
	st_ResultImage = st_OriginalImage.clone();

	int32_t s32_I;

	for (s32_I = 0; s32_I < s32_Num; s32_I++)
	{
		int x = rand() % st_OriginalImage.cols;
		int y = rand() % st_OriginalImage.rows;

		if (s32_I < s32_Num/2)
		{
			// Salt
			st_ResultImage.at<uchar>(y, x) = 255;
		}
		else
		{
			// Pepper
			st_ResultImage.at<uchar>(y, x) = 0;
		}
	}
}

 

MakeSaltandPepperNoise는 인자로 받은 개수만큼 random하게 선정된 Pixel 위치에 흑,백의 Noise를 생성하는 함수이다. Mat::clone()을 통해 원본 이미지를 복사한다. 이후 s32_Num 만큼 반복하여 무작위로 생성되는 픽셀에 대하여 Intensity를 255(Salt), 0(Pepper)로 생성한다. 각각의 개수는 인자로 받은 개수의 반으로 설정하였다.

 

위에서 생성한 9x9 Kernel을 이용하여 Gaussian Filtering을 진행하였고, 이전 포스팅에서 구현한 Median Filtering을 Kernel Size를 9로 하여 두 결과를 비교 해 보자.

Original Image , Gaussian Filtering , Median Filtering    -    10000개의 Noise

 

Salt and Pepper  Noise를 10000개 생성하였을 때 Gaussian Filtering과 Median Filtering의 결과는 다음과 같다. Gaussian Filtering은 Kernel 내의 인접 픽셀의 Intensity를 모두 반영한다. 따라서 Outlier를 제거하지 못하고 결국 Weight에 따른 평균을 내기 때문에 완전히 제거되지 않고 얼룩져 보이는 것을 확인할 수 있다. 하지만 Median Filter의 경우 Kernel 내의 인접 픽셀의 Intensity의 중간값을 해당 픽셀의 Intensity로 설정하기 때문에 Outlier는 중간값으로 올 확률이 적어 효과적으로 Salt and Pepper Noise를 제거한 것을 확인할 수 있다.

Original Image , Gaussian Filtering , Median Filtering    -    100개의 Noise

 

반면, Salt and Pepper Noise가 100개로 Resolution에 비해 상대적으로 적으면 Gaussian Kernel 내에 존재할 확률 역시 적어지기 때문에 평균을 구한다 해도 효과적으로 Smoothing을 통한 Noise 제거가 가능함을 확인할 수 있다.

 

Diagonal Edge Detection

Sobel Filter를 학습할 때 각 방향에 대한 차분을 계산하기 위하여 다음과 같이 Filter를 설계하였다.

 

 

 

그렇다면 대각선의 Edge를 추출하기 위해선 어떻게 Filter를 설계해야 할까?

 

 

135도의 대각선 Edge를 검출하고 한다면 위와 같은 상황일 것이다. Correlation Filter는 Filter Weight와 픽셀의 Intensity 분포가 유사한 부분에서 높은 Response를 보인다. 따라서, (0,0) , (2,2)의 부호가 반대이고 Weight의 절댓값이 같으면 Edge를 찾아낼 수 있다. 동시에 (1,0) , (1,2) & (0,1) , (2,1) 역시 동일한 관계를 가진다면 Edge일 때 Pixel Intensity를 고려했을 때 대각선 Edge를 가지는 부분에서 높은 Response를 보일 것이다. 

 

float f32_KernelData135[]  = {  2.f,  1.f,  0.f,
							   1.f,  0.f, -1.f,
							   0.f, -1.f, -2.f };
float f32_KernelData45[] = {  0.f,  1.f,  2.f,
							  -1.f,  0.f,  1.f,
							  -2.f, -1.f,  0.f };

Mat st_Kernel45(Size(3, 3), CV_32FC1, f32_KernelData45);
Mat st_Kernel135(Size(3, 3), CV_32FC1, f32_KernelData135);

 

따라서 45, 135도 대각선 방향에 대한 Edge를 검출하기 위한 Filter는 위와 같이 만들 수 있다. 이때, 45, 135도의 대각선 방향을 이미지 좌표계를 고려하여 생각해야 함을 주의하자.

 

앞서 구현한 Gaussian Filtering 과정을 이용하여 5x5 Gaussian Kernel를 생성하여 Fitlering을 진행한 후 위에서 생성한 Filter를 기반으로 Convolution 연산을 한 결과 아래와 같다. 

GaussianKernel = MakeGaussianFilter(5);
myKernelConv(st_OriginalImage, st_GaussianImage45, GaussianKernel);
myKernelConv(st_OriginalImage, st_GaussianImage135, GaussianKernel);

myKernelConv(st_GaussianImage45, st_ResultImage45, st_Kernel45);
myKernelConv(st_GaussianImage135, st_ResultImage135, st_Kernel135);

 

 

Original Image , 45 Diagonal Edge  Detection , 135 Diagonal Edge Detection

 

긴 흰색 봉 부분을 보면 확실히 구별 가능한데, 해당 부분은 135도 방향의 Edge가 존재하는 곳이다. 따라서 오른쪽 Edge 검출 결과를 보면 해당 부분에서 Response가 매우 큼을 확인할 수 있다. 반면, 사각 볼트 쪽을 보면 45도 방향의 Edge가 많이 존재하여 두번째 Edge 검출 결과에서 해당 부분이 높은 Response를 보임을 확인할 수 있다.

 

Gaussian Pyramid

Gaussian Pyramid는 Bluring과 DownSampling을 통해 이미지를 1/2씩 축소시키면서 Pyramid를 생성한 것이다. Bluring 이 후 DownSampling을 진행하기에 원본 이미지의 형태를 축소시킨 후에도 유지시킬 수 있다. 동일한 이미지를 다양한 크기로 조절하여 요구되는 작업을 수행할 때 사용된다. 이를 직접 구현 해 보자.

 

myGaussianPyramid()

vector<Mat> myGaussianPyramid(Mat src_img) {
	vector<Mat> Vec;

	Vec.push_back(src_img);
	for (int i = 0; i < 4; i++) {
		src_img = myGaussianFilter(src_img);
		src_img = mySampling(src_img);
		Vec.push_back(src_img);
	}

	return Vec;
}

 

전체적인 동작 과정은 다음과 같다. Input으로 Orginal Image를 받아서 Vector Vec에 저장한다. 이후 for문을 돌면서 Gaussain Filtering을 통해 Bluring을 진행한 후 DownSampling 한다. 이후 축소된 이미지를 다시 Vector Vec에 저장한다. 이때, src_img를 계속적으로 변형하기에 진행할수록 1/2가 축소(전체 이미지 사이즈는 1/4 축소)된 이미지가 Vec에 저장된다.

 

myGaussianFilter()

Mat myGaussianFilter(Mat src_img) {
	int width = src_img.cols ;
	int height = src_img.rows;
	int kernel[3][3] = { 1, 2, 1,
						2, 4, 2,
						1, 2, 1 };

	Mat dst_img(src_img.size(), CV_8UC3);

	uchar* srcData = src_img.data;
	uchar* dstData = dst_img.data;

	for (int y = 0; y < height; y++) {
		for (int x = 0; x < width; x++) {
			// 각 체널의 픽셀에 대하여 Gaussian Filtering 진행
			for (int c = 0; c < 3; c++) {
				dst_img.at<cv::Vec3b>(y, x)[c] = myKernerConv3x3(src_img, kernel, x, y, width, height, c);
			}
		}
	}

	return dst_img;
}

 

이전까지는 GrayScale Image에 대해서만 진행하였지만, Gaussain Pyramid는 Color Image를 대상으로 진행된다. 따라서 Gaussian Filtering 결과를 저장할 Mat 객체를 CV_8UC3로 설정하였다. 이후 동일하게 유요한 픽셀에 대하여 접근하는데, 3채널 각각에 접근하기 위해서 for문이 하나 더 필요하다. 

 

for (int c = 0; c < 3; c++) {
    dst_img.at<cv::Vec3b>(y, x)[c] = myKernerConv3x3(src_img, kernel, x, y, width, height, c);
}

 

Mat::at()으로 해당 픽셀의 Intensity에 접근한다. 이때, Color Image이기 때문에 Datatype은 Vec3b이며, [0]은 B, [1]은 G, [2] R Channel을 의미한다. 각각의 Channel의 Intensity를 myKernerConv3x3()의 반환값으로 설정한다.

 

myKernerConv3x3

int myKernerConv3x3(const cv::Mat& src_img, int kernel[][3], int x, int y, int width, int height, int c) {
	int sum = 0;
	int sumKernel = 0;

	for (int j = -1; j <= 1; j++) {
		for (int i = -1; i <= 1; i++) {
			if ((y + j) >= 0 && (y + j) < height && (x + i) >= 0 && (x + i) < width) {
				sum += src_img.at<cv::Vec3b>(y+j, x+i)[c] * kernel[i + 1][j + 1];
				sumKernel += kernel[i + 1][j + 1];
			}
		}
	}

	if (sumKernel != 0) { return sum / sumKernel; }
	else return sum;
}

 

이전에 구현한 myKernerConv3x3()과 유사하지만, 현재 Channel을 나타내는 Parameter c가 추가되었다. 

$(x,y)$는 현재 Filtering을 진행할 픽셀의 좌표이다. 따라서 Kernel내의 인접 픽셀에 접근하기 위해 위와 같이 2중  for문을 설계하였다. 이때, 가장자리 문제를 고려하기 위해 유효한 픽셀의 위치를 판단하고, 유효한 픽셀에 대하여 원본이미지의 해당 체널의 Intensity를 Kernel Weight와 곱하여 sum에 저장한다. 해당 과정을 Kernel 내부의 모든 픽셀에 대하여 진행한 후 Kernel Weight의 합이 0이 아니면 나눠 줌으로써 Normalization을 진행하고 아니면 그대로 sum을 반환한다.

 

결국 Gray Image에서 진행한 Gaussian Filtering을 Color Image에서는 각 체널에 접근하여 동일한 작업을 진행한 것으로 이해하면 된다. Gaussian Filtering을 진행한 Image에 대하여 DownSampling을 진행한다.

mySampling()

Mat mySampling(Mat src_img) {

	// src_img에 대하여 가로, 세로 각각 1/2 Resize => 전체 크기는 1/4
	int width = src_img.cols / 2;
	int height = src_img.rows / 2;

	Mat dst_img(height, width, CV_8UC3);

	for (int y = 0; y < height; y++) {
		for (int x = 0; x < width; x++) {

			// 각 체널에 대하여 DownSampling
			for (int c = 0; c < 3; c++) {
				dst_img.at<cv::Vec3b>(y, x)[c] = src_img.at<cv::Vec3b>(y * 2, x * 2)[c];
			}
		}
	}

	return dst_img;
}

 

가로, 세로 1/2 Scale로 Resize하기 때문에 결국 원본 이미지에서 2x2의 Sub Image가 DonSampling된 이미지의 1x1 픽셀이 되며, 해당 픽셀의 Intensity는 Sub Image의 왼쪽 상단 픽셀의 Intensity가 된다고 이해하면 된다.

 

dst_img.at<cv::Vec3b>(y, x)[c] = src_img.at<cv::Vec3b>(y * 2, x * 2)[c];

 

왜 Index가 위와 같이 설정되어야 되는지는 위 그림을 통해 쉽게 이해 가능하다.

 

따라서 실행 결과 원본 이미지에 대하여 1/2 축소된 이미지가 피라미드의 각 층을 이루게 된다.

 

Laplacian Pyramid

Laplacian Pyramid는 Gaussian Pyramid의 인접한 층끼리의 차이(Difference)를 구한 것이다. 즉, Gaussian Pyramid를 통해 생성된 다중 해상도의 이미지의 차분을 계산하는 것이다.  Laplacian Pyramid에서 밝은 부분이 Edge의 위치임을 확인할 수 있다. 따라서 Laplacian은 결국 다양한 해상도 영상에서의 Edge 정보를 가지고 있는 것이다.

 

 

이를 직접 구현 해 보자.

myLaplacianPyramid()

vector<Mat> myLaplacianPyramid(Mat src_img) {
	vector<Mat> Vec;

	for (int i = 0; i < 4; i++) {
		if (i != 3) {
			Mat high_img = src_img;						// 이전 층의 Image
			src_img = myGaussianFilter(src_img);
			src_img = mySampling(src_img);
			Mat low_img = src_img;						// 현재 층의 Image -> Resolution이 감소(Resize)
			resize(low_img, low_img, high_img.size());	// 이전 층과 동일하게 Resize
			Vec.push_back(high_img - low_img + 128);	// 차분값을 128을 더하여 보기 좋게 만듦
		}
		else {
			Vec.push_back(src_img);
		}
	}
	return Vec;
}

 

기본적으로 Gaussian Pyramid를 구현해야 한다. 따라서 Gaussian Filtering과 DownSampling 과정은 동일하다. 인접한 Pyramid 층의 차분을 구해야 하는데, DownSampling되면서 두 이미지의 Size가 달라지기 때문에 cv::resize()를 통해 높은 해상도의 이미지와 동일하게 크기를 재설정하였다. 이후 두 이미지의 뺀 뒤 128을 더하여 Edge부분이 더 잘 보이도록 하였다.

 

실행결과 다양한 해상도에서의 Edge 정보를 Laplacian Pyramid를 통해 취득할 수 있음을 알 수 있다. 그렇다면 원본 이미지를 복원할 수 있을까?

 

// 4. Laplacian Pyramid
Mat st_OriginalImage, ResultImage;

st_OriginalImage = imread("gear.jpg", 1);
vector<Mat> VecLap = myLaplacianPyramid(st_OriginalImage);
reverse(VecLap.begin(), VecLap.end());

// 복원
for (int i = 0; i < VecLap.size(); i++) {
    if (i == 0) {
        ResultImage = VecLap[i];
    }
else {
    resize(ResultImage, ResultImage, VecLap[i].size());
    ResultImage = ResultImage + VecLap[i] - 128;
}

string windowName = "LaplacianPyramid" + to_string(i); 
imshow(windowName, ResultImage);

 

myLaplacianPyramid()에서 인접한 Gaussian Pyramid Image를 빼고 128을 더했으므로 해당 결과에 대하여 인접한 Laplacian Pyramid Image를 더한 후 128을 빼주면 된다.

 

 

따라서 생성한 Gaussian Pyramid를 복원할 수 있다.