来自 DFT/IDFT 计算的负值

Negative values from DFT/IDFT calculation

我正在频域中创建一个截止频率为 0.015 的低通滤波器和相应的高通滤波器,如下所示:

auto getGaussianFilter(int32_t size) {
    const auto sigmaf = 0.015f*size; // cutoff freq.

    cv::Mat kernel = cv::Mat::zeros(size, size, CV_32FC1);

    for (auto fy = -size/2; fy < size/2; ++fy)
        for (auto fx = -size/2; fx < size/2; ++fx) 
            kernel.at<float>(fy + size/2, fx + size/2) = std::exp(-(fy*fy+fx*fx)/(2*sigmaf*sigmaf));
    return kernel;
} 

...

lowpassFilter  = getGaussianFilter(256);
highpassFilter = 1 - lowpassFilter; 

,并将这些应用于图像,如以下代码段。

std::vector<cv::Mat> planes { src, cv::Mat::zeros(src.size(), src.type()) }; 
std::vector<cv::Mat> fplanes { filter, filter };
cv::Mat complex, complexFilter;
cv::merge(planes, complex);
cv::merge(fplanes, complexFilter);

cv::dft(complex, complex, cv::DFT_COMPLEX_OUTPUT);
shift(complex); // swapping quadrants
cv::mulSpectrums(complex, complexFilter, complex, cv::DFT_ROWS);
shift(complex); 
cv::idft(complex, complex, cv::DFT_SCALE | cv::DFT_REAL_OUTPUT); 

输入图像为 256x256 且为对数尺度,低通滤波器应用于高通滤波输入图像的平方。

然后我想对最终 idft 的结果求平方根, 但它包含负值;因此,出现了几个 NaN 值。

applyFilter(src, dst, highpassFilter);
cv::pow(dst, 2, dst);
applyFilter(dst, dst, lowpassFilter);
cv::sqrt(dst, dst);  // NaN !

为什么会有负值,如何处理负值才能求平方根?


编辑:添加了shift

的代码
void shift(cv::Mat& src) {
    src = src(cv::Rect(0, 0, src.cols & -2, src.rows & -2)); 
    const auto cy = src.rows/2, cx = src.cols/2;
    cv::Mat q0(src, cv::Rect(0, 0, cx, cy));
    cv::Mat q1(src, cv::Rect(cx, 0, cx, cy));
    cv::Mat q2(src, cv::Rect(0, cy, cx, cy));
    cv::Mat q3(src, cv::Rect(cx, cy, cx, cy));
    cv::Mat tmp;
    q0.copyTo(tmp);
    q3.copyTo(q0);
    tmp.copyTo(q3);

    q1.copyTo(tmp);
    q2.copyTo(q1);
    tmp.copyTo(q2);
 } 

直到第三次或第四次阅读你的代码我才明白这一点...

std::vector<cv::Mat> fplanes { filter, filter };
cv::merge(fplanes, complexFilter);

您在这里创建了一个复值滤波器,它不是共轭对称的。应用此滤波器会创建一个非共轭对称的频谱,因此您会得到一个非实数逆变换。

你想要做的是拥有一个纯粹真实的过滤器:

std::vector<cv::Mat> fplanes { filter, cv::Mat::zeros(filter(), filter()) };
cv::merge(fplanes, complexFilter);

根据幅度和相位(而不是实部和虚部)来考虑频域滤波器。幅度会改变图像中每个频率分量幅度的变化方式,相移会改变每个频率分量。一般来说(除了非常特殊的情况)你不想移动图像中的频率分量,因为它会弄乱你的图像。因此,将频域滤波器的相位保持在 0。零相位意味着每个值都是非负且实数。

也就是说,您移动图像频谱、应用滤镜和移回的顺序:

cv::dft(complex, complex, cv::DFT_COMPLEX_OUTPUT);
shift(complex); // swapping quadrants
cv::mulSpectrums(complex, complexFilter, complex, cv::DFT_ROWS);
shift(complex); 
cv::idft(complex, complex, cv::DFT_SCALE | cv::DFT_REAL_OUTPUT); 

与简单地移动过滤器相同:

cv::dft(complex, complex, cv::DFT_COMPLEX_OUTPUT);
shift(complexFilter); // swapping quadrants
cv::mulSpectrums(complex, complexFilter, complex, 0);
cv::idft(complex, complex, cv::DFT_SCALE | cv::DFT_REAL_OUTPUT); 

另请注意,我将 cv::DFT_ROWS 替换为 0。根据 the documentation,该标志表示每个图像行都是独立的一维变换。我的猜测是,这只对 CCS 打包矩阵很重要,但最好还是将其排除在外。