c++ - 选择正确的内核大小以进行中值模糊处理以减少噪声

Question

我是图像处理的新手。我们需要从图像中获得具有亚像素精度的圆心。我使用了中值模糊来减少噪音。图像的一部分如下所示。下面给出了我获取圆边界的步骤

用medianBlur降低噪音

通过阈值API来应用OTSU阈值

使用findContours方法识别圆的边界。

当对 medianBlur使用不同的内核大小时，我得到不同的结果。我选择了 medianBlur来保持边缘。我尝试了内核大小3、5和7。现在，我对于 medianBlur使用正确的内核大小感到困惑。

如何确定正确的内核大小？

是否有任何科学方法来确定medianBlur的正确内核大小？

Answer 1

我将在此处为您提供有关如何找到这些磁盘质心的两个建议，您可以根据所需的精度级别选择一个。
首先，使用轮廓并不是最好的方法。轮廓线很大程度上取决于在阈值化时哪些像素恰好落在对象内，而噪声对这些像素的影响很大。
更好的方法是找到磁盘的质心(或更确切地说是一阶矩)。 Read Wikipedia了解有关图像分析中更多信息的信息。关于矩的一件好事是，我们可以将像素值用作权重，从而提高精度。
您可以从二进制形状的轮廓计算其矩，但在这种情况下，不能使用图像强度。 OpenCV有一个函数 cv::moments ，它可以计算整个图像的矩，但是我不知道有一个函数可以分别为每个对象执行此操作。因此，我将使用DIPlib进行这些计算(我是作者)。

关于过滤:
只要对象离图像边缘足够远，任何行为良好的线性平滑都不应影响对象的质心。靠近边缘会导致模糊在对象最靠近边缘的那一侧进行其他操作，这会产生偏差。
任何非线性平滑滤波器都可以更改质心。请避免使用中值过滤器。
因此，我建议您使用高斯滤波器，这是性能最好的线性平滑滤波器。

方法1:使用二进制形状的矩:
首先，我要达到阈值而不会出现任何形式的模糊。

import diplib as dip

a = dip.ImageRead('/Users/cris/Downloads/Ef8ey.png')
a = a(1)  # Use green channel only, simple way to convert to gray scale

_, t = dip.Threshold(a)
b = a<t

m = dip.Label(b)

msr = dip.MeasurementTool.Measure(m, None, ['Center'])
print(msr)

这个输出

  |                  Center | 
- | ----------------------- | 
  |       dim0 |       dim1 | 
  |       (px) |       (px) | 
- | ---------- | ---------- | 
1 |      18.68 |      9.234 | 
2 |      68.00 |      14.26 | 
3 |      19.49 |      48.22 | 
4 |      59.68 |      52.42 | 

现在我们可以对输入图像 a进行平滑处理，然后再次进行计算:

a = dip.Gauss(a,2)

_, t = dip.Threshold(a)
b = a<t
m = dip.Label(b)
msr = dip.MeasurementTool.Measure(m, None, ['Center'])
print(msr)

  |                  Center | 
- | ----------------------- | 
  |       dim0 |       dim1 | 
  |       (px) |       (px) | 
- | ---------- | ---------- | 
1 |      18.82 |      9.177 | 
2 |      67.74 |      14.27 | 
3 |      19.51 |      47.95 | 
4 |      59.89 |      52.39 | 

您会看到质心有一些小的变化。

方法2:使用灰度级矩:
在这里，我们使用误差函数将伪阈值应用于图像。这样做是将对象像素设置为1，背景像素设置为0，但是边缘周围的像素保留一些中间值。有人将此称为“模糊阈值”。这两个图像显示了正常(“硬”阈值)和错误功能片段(“模糊阈值”):

通过使用此模糊阈值，我们保留了有关边缘的精确(子像素)位置的更多信息，可在计算一阶矩时使用这些信息。

import diplib as dip

a = dip.ImageRead('/Users/cris/Downloads/Ef8ey.png')
a = a(1)  # Use green channel only, simple way to convert to gray scale

_, t = dip.Threshold(a)
c = dip.ContrastStretch(-dip.ErfClip(a, t, 30))

m = dip.Label(a<t)
m = dip.GrowRegions(m, None, -2, 2)

msr = dip.MeasurementTool.Measure(m, c, ['Gravity'])
print(msr)

这个输出

  |                 Gravity | 
- | ----------------------- | 
  |       dim0 |       dim1 | 
  |       (px) |       (px) | 
- | ---------- | ---------- | 
1 |      18.75 |      9.138 | 
2 |      67.89 |      14.22 | 
3 |      19.50 |      48.02 | 
4 |      59.79 |      52.38 | 

现在我们可以对输入图像 a进行平滑处理，然后再次进行计算:

a = dip.Gauss(a,2)

_, t = dip.Threshold(a)
c = dip.ContrastStretch(-dip.ErfClip(a, t, 30))
m = dip.Label(a<t)
m = dip.GrowRegions(m, None, -2, 2)
msr = dip.MeasurementTool.Measure(m, c, ['Gravity'])
print(msr)

  |                 Gravity | 
- | ----------------------- | 
  |       dim0 |       dim1 | 
  |       (px) |       (px) | 
- | ---------- | ---------- | 
1 |      18.76 |      9.094 | 
2 |      67.87 |      14.19 | 
3 |      19.50 |      48.00 | 
4 |      59.81 |      52.39 | 

您可以看到这次差异较小，因为测量更精确。

在二进制情况下，有和没有平滑的质心差异为:

array([[ 0.14768417, -0.05677508],
       [-0.256     ,  0.01668085],
       [ 0.02071882, -0.27547569],
       [ 0.2137167 , -0.03472741]])

在灰度情况下，差异为:

array([[ 0.01277204, -0.04444567],
       [-0.02842993, -0.0276569 ],
       [-0.00023144, -0.01711335],
       [ 0.01776011,  0.01123299]])

如果质心的度量单位是µm，而不是px，那是因为您的图像文件包含像素大小信息。测量功能将使用此功能为您提供真实的测量结果(质心坐标为w.r.t.左上像素)。如果您不希望这样做，可以重置图像的像素大小:

a.SetPixelSize(1)

C++中的两种方法
这是上述代码对C++的翻译，包括显示步骤，用于仔细检查阈值产生了正确的结果:

#include "diplib.h"
#include "dipviewer.h"
#include "diplib/simple_file_io.h"
#include "diplib/linear.h"       // for dip::Gauss()
#include "diplib/segmentation.h" // for dip::Threshold()
#include "diplib/regions.h"      // for dip::Label()
#include "diplib/measurement.h"
#include "diplib/mapping.h"      // for dip::ContrastStretch() and dip::ErfClip()

int main() {

   auto a = dip::ImageRead("/Users/cris/Downloads/Ef8ey.png");
   a = a[1];  // Use green channel only, simple way to convert to gray scale

   dip::Gauss(a, a, {2});

   dip::Image b;
   double t = dip::Threshold(a, b);
   b = a < t; // Or: dip::Invert(b,b);

   dip::viewer::Show(a);
   dip::viewer::Show(b); // Verify that the segmentation is correct
   dip::viewer::Spin();

   auto m = dip::Label(b);

   dip::MeasurementTool measurementTool;
   auto msr = measurementTool.Measure(m, {}, { "Center"});
   std::cout << msr << '\n';

   auto c = dip::ContrastStretch(-dip::ErfClip(a, t, 30));
   dip::GrowRegions(m, {}, m, -2, 2);

   msr = measurementTool.Measure(m, c, {"Gravity"});
   std::cout << msr << '\n';

   // Iterate through the measurement structure:
   auto it = msr["Gravity"].FirstObject();
   do {
      std::cout << "Centroid coordinates = " << it[0] << ", " << it[1] << '\n';
   } while(++it);
}