用 OpenCV 去卷积？

Posted 2023-02-17

【中文标题】用 OpenCV 去卷积？

【英文标题】：Deconvolution with OpenCV?

【发布时间】：2013-11-28 12:56:07

【问题描述】：

有没有办法用 OpenCV 进行反卷积？

这里显示的改进给我留下了深刻的印象

并且想将此功能也添加到我的软件中。

编辑（赏金的其他信息。）

我还没有弄清楚如何实现反卷积。 这段代码可以帮助我锐化图像，但我认为反卷积可以做得更好。

void ImageProcessing::sharpen(QImage & img)

    IplImage* cvimg = createGreyFromQImage( img );
    if ( !cvimg ) return;

    IplImage* gsimg = cvCloneImage(cvimg );
    IplImage* dimg = cvCreateImage( cvGetSize(cvimg), IPL_DEPTH_8U, 1 );
    IplImage* outgreen = cvCreateImage( cvGetSize(cvimg), IPL_DEPTH_8U, 3 );
    IplImage* zeroChan = cvCreateImage( cvGetSize(cvimg), IPL_DEPTH_8U, 1 );
    cvZero(zeroChan);

    cv::Mat smat( gsimg, false );
    cv::Mat dmat( dimg, false );

    cv::GaussianBlur(smat, dmat, cv::Size(0, 0), 3);
    cv::addWeighted(smat, 1.5, dmat, -0.5 ,0, dmat);
    cvMerge( zeroChan, dimg, zeroChan, NULL, outgreen);

    img = IplImage2QImage( outgreen );
    cvReleaseImage( &gsimg );
    cvReleaseImage( &cvimg );
    cvReleaseImage( &dimg );
    cvReleaseImage( &outgreen );
    cvReleaseImage( &zeroChan );

希望得到有用的提示！

【问题讨论】：

要求代码的问题必须表明对所解决问题的最低理解。包括尝试的解决方案、为什么它们不起作用以及预期的结果。另请参阅：堆栈溢出问题清单

发现了一篇关于该主题的有趣文章：yuzhikov.com/articles/BlurredImagesRestoration1.htm

【参考方案1】：

当然，您可以使用 OpenCV 编写反卷积代码。但是还没有准备好使用的函数（还）。

要开始使用，您可以查看this Example，它显示了使用 OpenCV 在 Python 中实现 Wiener Deconvolution。

这里是 another Example 使用 C，但这是从 2012 年开始的，所以可能已经过时了。

【讨论】：

谢谢。它看起来比我预期的要复杂得多。第二个例子很简单，应该可以在 OpenCV2 中理解和重写。如果没有其他建议，我将等待一段时间再提出其他建议（可能是其他方法，例如最近邻解卷积）并选择您的建议。

最后我学到了足够多的python来运行第一个例子。汽车图像的结果是惊人的。我相信甚至可以从具有图像的一些上下文信息的运动分析（所有在同一方向的小线）中确定 psf。然而，在我的图像（显微镜中的细胞）上，结果并没有那么好。可能我没有正确配置 PFS，并且不像汽车示例那样明显。

了解您的 PSF 至关重要，但不知道您的 PSF 反卷积是一项非常艰巨的任务（您必须猜到）。如果您详细描述拍摄图像的设置，也许人们可以帮助您。

好的，我会添加信息。实际上，我的方向是 Z。这可能意味着，PSF 不能从一个图像中定义，而是从一堆捕获的图像中定义，同时在同一位置上下移动显微镜。我认为从当前的问题描述来看，这并不清楚。

【参考方案2】：

最近邻去卷积是一种技术，通常用于光学显微镜中 Z 平面中的一堆图像。这篇评论论文：让-巴蒂斯特·西巴里塔。 Deconvolution Microscopy。 Adv Biochem Engin/Biotechnol (2005) 95: 201–243 涵盖了很多使用的技术，包括您感兴趣的技术。这也是一个很好的介绍：http://blogs.fe.up.pt/BioinformaticsTools/microscopy/

这个 numpy+scipy python 示例展示了它是如何工作的：

from pylab import *
import numpy
import scipy.ndimage

width = 100
height = 100
depth = 10
imgs = zeros((height, width, depth))

# prepare test input, a stack of images which is zero except for a point which has been blurred by a 3D gaussian
#sigma = 3
#imgs[height/2,width/2,depth/2] = 1
#imgs = scipy.ndimage.filters.gaussian_filter(imgs, sigma)

# read real input from stack of images img_0000.png, img_0001.png, ... (total number = depth)
# these must have the same dimensions equal to width x height above
# if imread reads them as having more than one channel, they need to be converted to one channel
for k in range(depth):
    imgs[:,:,k] = scipy.ndimage.imread( "img_%04d.png" % (k) )

# prepare output array, top and bottom image in stack don't get filtered
out_imgs = zeros_like(imgs)
out_imgs[:,:,0] = imgs[:,:,0]
out_imgs[:,:,-1] = imgs[:,:,-1]

# apply nearest neighbor deconvolution
alpha = 0.4 # adjustabe parameter, strength of filter
sigma_estimate = 3 # estimate, just happens to be same as the actual

for k in range(1, depth-1):
    # subtract blurred neighboring planes in the stack from current plane
    # doesn't have to be gaussian, any other kind of blur may be used: this should approximate PSF
    out_imgs[:,:,k] = (1+alpha) * imgs[:,:,k]  \
        - (alpha/2) * scipy.ndimage.filters.gaussian_filter(imgs[:,:,k-1], sigma_estimate) \
        - (alpha/2) * scipy.ndimage.filters.gaussian_filter(imgs[:,:,k+1], sigma_estimate)

# show result, original on left, filtered on right
compare_img = copy(out_imgs[:,:,depth/2])
compare_img[:,:width/2] = imgs[:,:width/2,depth/2]
imshow(compare_img)
show()

【讨论】：

太酷了！谢谢你的链接！第二天我会阅读链接，测试代码，很可能会给你赏金。 P.S.：同时我考虑在我的应用程序中嵌入 Python 或将所有图像处理功能移动到用 python 实现的单独进程。有这么多强大的软件包......

对不起，我刚开始学习python。你能用真实的图像帮助测试代码吗？我试过 img = scipy.ndimage.imread("19756782_g.png") 其中 png 是灰色的 8 位图像。然后我把这些图像放在集合/（数组？）imgs 中。不幸的是，我在 k in range(1, imgs.shape[2]-1) 的循环中收到错误“'tuple' object has no attribute 'shape'”。

@ValentinHeinitz：这个样本对你有用吗？请查看上面的更新代码，如果您有更多问题，请告诉我。

@ValentinHeinitz：“然后我把这样的图像放在集合/（数组？）imgs” - imgs 它是一个 numpy 数组（ndarray），一个完整的图像堆栈在一个二进制结构中（所以它是 3D，宽度 x 高度 x 深度）。没有方便的方法来附加它 :) 但请参阅上文了解如何从单独的图像文件中读取它。

赏金的获胜者是...:-) 说真的，感谢您提供链接和代码示例。我的问题仍未如我所愿解决，但我现在明白了一些。在 Z 方向捕获图像时，我必须使用 PSF 并减少抖动。一旦我找到了如何在 OpenCV 中执行此操作，我将发布示例。

【参考方案3】：

您提供的示例图像实际上是 Lucy-Richardson 反卷积的一个很好的示例。 OpenCV 库中没有针对这种反卷积方法的内置函数。在 Matlab 中，您可以使用带有“deconvlucy.m”函数的反卷积。实际上，您可以通过键入“open”或“edit”来查看Matlab中某些函数的源代码。 下面，我尝试在 OpenCV 中简化 Matlab 代码。

// Lucy-Richardson Deconvolution Function
// input-1 img: NxM matrix image
// input-2 num_iterations: number of iterations
// input-3 sigma: sigma of point spread function (PSF)
// output result: deconvolution result

// Window size of PSF
int winSize = 10 * sigmaG + 1 ;

// Initializations
Mat Y = img.clone();
Mat J1 = img.clone();
Mat J2 = img.clone();
Mat wI = img.clone(); 
Mat imR = img.clone();  
Mat reBlurred = img.clone();    

Mat T1, T2, tmpMat1, tmpMat2;
T1 = Mat(img.rows,img.cols, CV_64F, 0.0);
T2 = Mat(img.rows,img.cols, CV_64F, 0.0);

// Lucy-Rich. Deconvolution CORE

double lambda = 0;
for(int j = 0; j < num_iterations; j++) 
       
    if (j>1) 
        // calculation of lambda
        multiply(T1, T2, tmpMat1);
        multiply(T2, T2, tmpMat2);
        lambda=sum(tmpMat1)[0] / (sum( tmpMat2)[0]+EPSILON);
        // calculation of lambda
    

    Y = J1 + lambda * (J1-J2);
    Y.setTo(0, Y < 0);

    // 1)
    GaussianBlur( Y, reBlurred, Size(winSize,winSize), sigmaG, sigmaG );//applying Gaussian filter 
    reBlurred.setTo(EPSILON , reBlurred <= 0); 

    // 2)
    divide(wI, reBlurred, imR);
    imR = imR + EPSILON;

    // 3)
    GaussianBlur( imR, imR, Size(winSize,winSize), sigmaG, sigmaG );//applying Gaussian filter 

    // 4)
    J2 = J1.clone();
    multiply(Y, imR, J1);

    T2 = T1.clone();
    T1 = J1 - Y;


// output
result = J1.clone();

以下是一些示例和结果。

Example results with Lucy-Richardson deconvolution

访问我的博客Here，您可以在其中访问整个代码。

【讨论】：

你的博客不错！！

【参考方案4】：

我不确定你是否理解什么是反卷积。反卷积背后的想法是从图像中移除检测器响应。这在天文学中很常见。

例如，如果您将 CCD 安装到望远镜上，那么您拍摄的任何图像都是您在天空中看到的东西和光学系统响应的卷积。望远镜（或相机镜头或其他任何东西）将具有一些点扩散功能（PSF）。也就是说，如果你看一个很远的点源，比如一颗星星，当你拍摄它的图像时，星星会模糊几个像素。这种模糊——点扩散——是你想要去除的。如果您非常了解光学系统的点扩散函数，那么您可以从图像中对 PSF 进行反卷积，从而获得更清晰的图像。

除非您碰巧知道光学器件的 PSF（测量起来很重要！），否则您应该寻找其他选项来锐化图像。我怀疑 OpenCV 是否内置了 Richardson-Lucy 算法。

【讨论】：

感谢您提及 Richardson-Lucy 算法。同时，我或多或少地按照您描述的方式理解了反卷积。不幸的是，Alex I 发布的示例不适用于我的图像，因为我的 PSF 不是高斯模糊。当我尝试使用一堆未锐化且模糊不清的图像来测试我的自动对焦算法时，我也遇到了 PSF。它没有用，SO 上的某个人将我指向“PSF”。

如果你想使用反卷积，你需要测量你的光学器件的 PSF。我想你可以猜到，但你会在黑暗中摸索。要对其进行测量，您需要在视场中的各个点对白光点光源进行成像。光学 PSF 通常近似为高斯分布，但您至少需要知道分布的宽度。真的，除非您愿意测量您的 PSF，否则我强烈建议您寻找另一种锐化算法。