机器学习的数学-sobel算子

前言：

文章分为四部分来介绍：

• 梯度
• 方向导数
• Sobel边缘检测
• C语言实现Sobel算子

Part1-梯度(Gradient)

        梯度是对多元函数的推广，是多元函数对各自变量偏导数形成的向量，其数学意义类似于一元函数的导数。函数的梯度是对自身所有的自变量的偏导数形成的向量。

        公式注解：称为梯度算子，梯度算子作用于多元函数可得到一个向量。

        与一元函数的导数类似，梯度决定了多元函数的单调性与极值。梯度与单调性的关系存在如下两种情况：

        情况一：与的夹角不超过，则在方向上函数单调递增，即有：，当二者夹角为0度时，沿着梯度方向时函数单调递增。

        情况二： 与的夹角超过， 则在方向上函数单调递减，即有：，当二者夹角为时， 沿着梯度反方向时函数单调递减。

        可导函数在极值点处的梯度为0，同样地，梯度为0的点称为多元函数的驻点。

Part2-方向导数(Directional Derivative)

        偏导数沿着坐标轴方向求导，而方向导数可以沿着任意方向求导。方向导数是自变量沿着某一方向变化时的导数。对于多元函数，其在点处沿着方向的方向导数定义为：

        即沿着方向趋于时，函数值的改变量与此方向的步长的比值的极限。如果函数在点处的偏导数存在，则其方向导数可以根据梯度计算得出：

公式注解："  "表示向量内积运算。

        如果指向各个自变量的坐标轴方向，此时方向导数等于偏导数。我们进一步考虑二元函数的方向导数，对于函数，如果，沿着直线：

趋向于，方向导数为：

Part3-sobel边沿检测

        边缘检测的目的是找出图像中所有边缘点，也就是变化剧烈的地方，导数反映了函数在某一点的变化速度，因此导数绝对值较大的点可能是边缘点。

        把一副灰度图像看作目标函数，函数值对应图像坐标处的像素值，计算图像在(水平)，(垂直)方向的偏导数的绝对值，记为，如果(threshold表示设定的阈值)，可将此图像坐标位置标记为边缘点。

        由于图像存储的信号为离散数字信号，因此只能取整数，图像的原点位于图像的左上角，偏导数用差分近似代替，sobel边沿检测方向的偏导数可用如下中心差分公式近似计算：

        为了抵抗噪声，用三行的中心差分值进行综合，中间一行权重为2，第一行和第三行的权重为1，可得到对应的卷积核为：

        同理，在y方向的偏导数的计算公式为：

对应卷积核为：

Part4-C语言实现Sobel算子

void sobel(unsigned char imageIn[img_H][img_W], unsigned char imageOut[img_H][img_W], unsigned char Threshold)
{
  //卷积核大小
  short KERNEL_SIZE = 3;
  short xStart = KERNEL_SIZE / 2;
  short xEnd = img_W - KERNEL_SIZE / 2;
  short yStart = KERNEL_SIZE / 2;
  short yEnd = img_H - KERNEL_SIZE / 2;
  short i, j, k;
  short temp[4];
  for (i = yStart; i < yEnd; i++)
  {
    for (j = xStart; j < xEnd; j++)
    {
      //计算不同方向梯度幅值
      temp[0] = -(short) imageIn[i - 1][j - 1] + (short) imageIn[i - 1][j + 1]     //{{-1, 0, 1},
        - (short) imageIn[i][j - 1] + (short) imageIn[i][j + 1]        // {-1, 0, 1},
          - (short) imageIn[i + 1][j - 1] + (short) imageIn[i + 1][j + 1];    // {-1, 0, 1}};
      
      temp[1] = -(short) imageIn[i - 1][j - 1] + (short) imageIn[i + 1][j - 1]     //{{-1, -1, -1},
        - (short) imageIn[i - 1][j] + (short) imageIn[i + 1][j]       // { 0,  0,  0},
          - (short) imageIn[i - 1][j + 1] + (short) imageIn[i + 1][j + 1];    // { 1,  1,  1}};
      
      temp[2] = -(short) imageIn[i - 1][j] + (short) imageIn[i][j - 1]       //  0, -1, -1
        - (short) imageIn[i][j + 1] + (short) imageIn[i + 1][j]       //  1,  0, -1
          - (short) imageIn[i - 1][j + 1] + (short) imageIn[i + 1][j - 1];    //  1,  1,  0
      
      temp[3] = -(short) imageIn[i - 1][j] + (short) imageIn[i][j + 1]       // -1, -1,  0
        - (short) imageIn[i][j - 1] + (short) imageIn[i + 1][j]       // -1,  0,  1
          - (short) imageIn[i - 1][j - 1] + (short) imageIn[i + 1][j + 1];    //  0,  1,  1
      
      temp[0] = abs(temp[0]);
      temp[1] = abs(temp[1]);
      temp[2] = abs(temp[2]);
      temp[3] = abs(temp[3]);
      
      //找出梯度幅值最大值
      for (k = 1; k < 4; k++)
      {
        if (temp[0] < temp[k])
        {
          temp[0] = temp[k];
        }
      }
      
      if (temp[0] > Threshold)
      {
        imageOut[i][j] = 255;
      }
      else
      {
        imageOut[i][j] = 0;
      }
    }
  }
}