同态滤波是属于图像增强的一个小算法,其原理和代码实现在众多博客中均有提及,再此,只对学习中一些自认为有用的知识点进行总结。
实现和学习过程中的一些总结:
- 同态滤波类似于灰度变换,都是对灰度级进行动态范围压缩和对比度增强来改善图像的外观。
- 图像
f
(
x
,
y
)
f(x,y)
f(x,y)可以表示为其照射分量
i
(
x
,
y
)
i(x,y)
i(x,y)和反射分量
r
(
x
,
y
)
r(x,y)
r(x,y)的乘积,即:
f
(
x
,
y
)
=
i
(
x
,
y
)
r
(
x
,
y
)
f(x,y)=i(x,y)r(x,y)
f(x,y)=i(x,y)r(x,y)
- 照射分量对应图像中的低频部分,在光源中来说对应漫反射光,反射分量对应图像中高频部分,在光源中对于镜面反射光。
- 在算法开始,先要对图像进行对数处理,若
f
(
x
,
y
)
f(x,y)
f(x,y)中有任何零值,则必须将图像加1,以避免
I
n
(
0
)
In(0)
In(0)的出现,然后从最终结果中减去1。
一些未出现的工具函数可以参考我之前的博客:图像基础频率域滤波
代码实现如下:
cv::Mat calc_Homomorphic(cv::Size size, double gamma_L, double gamma_H, double c, double D0)
{
cv::Mat result(size, CV_64FC1);
int cx = size.width / 2;
int cy = size.height / 2;
for (int i = 0; i < result.rows; ++i)
{
double* p = result.ptr<double>(i);
for (int j = 0; j < result.cols; ++j)
{
double d_2 = std::pow(i - cy, 2) + std::pow(j - cx, 2);
p[j] = (gamma_H - gamma_L) * (1 - std::exp(-c * d_2 / (D0 * D0))) + gamma_L;
}
}
cv::Mat planes[] = { result.clone(), result.clone() };
cv::Mat Homomorphic;
merge(planes, 2, Homomorphic);
return Homomorphic;
}
int main()
{
std::string path = "PET_image.tif";
cv::Mat src = cv::imread(path, cv::IMREAD_GRAYSCALE);
if (!src.data) {
std::cout << "Could not open or find the image" << std::endl;
return -1;
}
cv::Mat srclogMat;
src.convertTo(srclogMat, CV_64FC1);
srclogMat = srclogMat + cv::Scalar::all(1);
cv::log(srclogMat, srclogMat);
cv::Mat paddedMat;
int m = cv::getOptimalDFTSize(src.rows);
int n = cv::getOptimalDFTSize(src.cols);
cv::copyMakeBorder(srclogMat, paddedMat, 0, m - src.rows, 0, n - src.cols, cv::BORDER_CONSTANT, cv::Scalar::all(0));
cv::Mat paddedMatdouble;
paddedMat.convertTo(paddedMatdouble, CV_64FC1);
Dftshift(paddedMatdouble, paddedMatdouble);
cv::Mat complexF;
cv::dft(paddedMatdouble, complexF, cv::DFT_COMPLEX_OUTPUT);
double gamma_L = 0.25, gamma_H = 2.0;
double c = 1.0, D0 = 160;
cv::Mat Homomorphic = calc_Homomorphic(paddedMat.size(), gamma_L, gamma_H, c, D0);
cv::Mat complexFH;
cv::Mat iDft;
cv::multiply(complexF, Homomorphic, complexFH);
cv::idft(complexFH, iDft, cv::DFT_REAL_OUTPUT | cv::DFT_SCALE);
Dftshift(iDft, iDft);
cv::exp(iDft, iDft);
iDft = iDft - cv::Scalar::all(1);
cv::Mat result = iDft(cv::Rect(0, 0, src.cols, src.rows)).clone();
MinusToZero(result, result);
normalize(result, result, 0, 1, cv::NORM_MINMAX);
result.convertTo(result, CV_8UC1, 255);
cv::imwrite("ascend_ET.png", result);
cv::waitKey(0);
return 0;
}
测试结果:
|