알고리즘 - Edge 추출 방법 ( Phase Congruency )

영상에서 에지(edge)나 코너(corner) 같은 **특징(feature)**을 밝기(intensity) 값에 의존하지 않고 검출하는 방법입니다.
이는 인간 시각 시스템(HVS, Human Visual System)이 위상(phase) 정보를 통해 형태와 구조를 인식

픽셀의 절대 밝기 변화가 아니라 위상이 일치하는 정도를 측정합니다. 즉, 조명 변화나 콘트라스트 변화에 강건하며, 동일한 구조를 다양한 밝기 조건에서도 안정적으로 검출

 

2. 기본 원리
이미지를 주파수 영역(특히 다중 스케일 및 다중 방향)으로 변환했을 때, 각 주파수 성분은 **진폭(amplitude)**과 위상(phase) 정보

• 에지나 코너 부근에서는 여러 주파수 성분의 위상이 같은 위치에서 최대값

 

 

3. 장점
1. 조명 변화 불변성 — 밝기나 명암 대비가 달라도 특징 위치가 동일하게 검출됨.
2. 스케일/방향 분석 가능 — 멀티스케일, 멀티방향 필터 사용으로 정밀한 특징 분석 가능.
3. 정규화 필요 없음 — 이미지 명암값의 전처리(normalization)가 불필요.

 

4. 구현 방식

가장 널리 쓰이는 구현은 Log-Gabor 필터를 다중 스케일·다중 방향으로 적용한 뒤,
각 위치에서 위상 일치도를 계산하는 방법입니다.

구현 절차:

1. 이미지 → Log-Gabor 필터 뱅크 변환
2. 각 스케일·방향별 진폭과 위상 계산
3. 모든 스케일·방향에서 위상 평균 계산
4. 위상 일치도(phase congruency) 산출
5. 결과를 특징 맵(feature map)으로 변환

 

 

소스 코드

#include <opencv2/opencv.hpp>
#include <iostream>
#include <vector>
#include <cmath>

using namespace cv;
using namespace std;

// 간단한 Log-Gabor 필터 생성 (주파수 영역)
Mat createLogGaborKernel(Size size, double wavelength, double sigmaOnf, double theta) {
    Mat kernel(size, CV_32FC2, Scalar(0, 0));
    Point center(size.width / 2, size.height / 2);

    double thetaRad = theta * CV_PI / 180.0;
    double cosTheta = cos(thetaRad);
    double sinTheta = sin(thetaRad);

    for (int y = 0; y < size.height; y++) {
        for (int x = 0; x < size.width; x++) {
            double u = (x - center.x) / (double)size.width;
            double v = (y - center.y) / (double)size.height;

            // 회전 좌표계로 변환
            double uRot = u * cosTheta + v * sinTheta;
            double vRot = -u * sinTheta + v * cosTheta;

            double radius = sqrt(uRot * uRot + vRot * vRot);
            if (radius == 0) continue;

            // Log-Gabor 필터 방정식
            double logGabor = exp((-(log(radius / (1.0 / wavelength)) * log(radius / (1.0 / wavelength)))) /
                                  (2 * log(sigmaOnf) * log(sigmaOnf)));

            kernel.at<Vec2f>(y, x)[0] = logGabor; // 실수부
            kernel.at<Vec2f>(y, x)[1] = 0.0f;     // 허수부
        }
    }
    return kernel;
}

int main() {
    // 입력 이미지
    Mat imgGray = imread("test.png", IMREAD_GRAYSCALE);
    if (imgGray.empty()) {
        cerr << "이미지를 불러올 수 없습니다." << endl;
        return -1;
    }
    imgGray.convertTo(imgGray, CV_32F, 1.0 / 255.0);

    // FFT 준비
    Mat padded;
    int m = getOptimalDFTSize(imgGray.rows);
    int n = getOptimalDFTSize(imgGray.cols);
    copyMakeBorder(imgGray, padded, 0, m - imgGray.rows, 0, n - imgGray.cols, BORDER_CONSTANT, Scalar::all(0));

    Mat planes[] = {padded.clone(), Mat::zeros(padded.size(), CV_32F)};
    Mat complexImg;
    merge(planes, 2, complexImg);
    dft(complexImg, complexImg);

    // 다중 스케일 & 방향 설정
    vector<double> wavelengths = {3.0, 5.0, 8.0}; // 스케일
    vector<double> orientations = {0, 45, 90, 135}; // 방향

    Mat sumAmp = Mat::zeros(padded.size(), CV_32F);
    Mat sumPC  = Mat::zeros(padded.size(), CV_32F);

    for (double wl : wavelengths) {
        for (double angle : orientations) {
            // Log-Gabor 필터 생성
            Mat logGabor = createLogGaborKernel(padded.size(), wl, 0.55, angle);

            // 주파수 영역에서 곱셈
            Mat filtered;
            mulSpectrums(complexImg, logGabor, filtered, 0);

            // 역 FFT
            Mat invDFT;
            idft(filtered, invDFT, DFT_SCALE | DFT_REAL_OUTPUT);

            // 진폭과 위상 계산
            Mat magImg, phaseImg;
            magnitude(invDFT, Mat::zeros(invDFT.size(), CV_32F), magImg);
            phase(invDFT, Mat::zeros(invDFT.size(), CV_32F), phaseImg);

            // 위상 평균(단순 예시)
            sumAmp += magImg;
            sumPC += magImg.mul(cos(phaseImg - mean(phaseImg)[0]));
        }
    }

    // Phase Congruency 계산
    Mat PC = sumPC / (sumAmp + 1e-9);

    // 시각화
    normalize(PC, PC, 0, 1, NORM_MINMAX);
    imshow("Phase Congruency", PC);
    waitKey(0);

    return 0;
}