2019-09-30 21:35:38 qq_41413256 阅读数 54

目标检测(object detection)

目标检测的目的

目标检测是找到图像中感兴趣的目标(物体),确定他们的位置和大小。我的理解就是将目标图中的物体用算法(传统算法或则是基于深度学习的算法)实现,用一个个的框框框出来,然后识别框框中的物体是什么,也是为下一的目标识别做准备。
机器视觉中的图像处理有四大任务:

  1. 分类(classification):假定给定一幅图像或者一段视频,要判断里面有什么类别的目标。
  2. 定位(location):定位出目标的位置,一般用边框将物体圈出来。
  3. 检测(detection):定位出目标的位置并检测出目标是什么。

传统的目标检测方法

  1. Viola-Jones(VJ):采用积分图特征+AdaBoost方法进行人脸检测等。
  2. HOG+SVM用于行人检测,通过对行人目标候选区域提取HOG特征,结合SVM分类器进行判定。
  3. DPM:基于HOG特征检测的一个变种,该方法是非深度学习方法中检测效果最好,性能最优的。

基于深度学习的目标检测

  1. One-stageYOLO和SSD系列):直接回归初目标框位置,不用产生候选框。
  2. Two-stageFaster RCNN系列):利用RPN网络对候选框进行推荐。
    上面写了一大堆,其实我也不是很清楚要做什么,这只是一个概述,方便我今后学习用的。就酱
2019-08-21 17:06:49 baidu_34971492 阅读数 259

FPGA图像处理之边缘检测算法的实现

1.背景知识

边缘检测是图像处理和计算机视觉中的基本问题,边缘检测的目的是标识数字图像中亮度变化明显的点。图像属性中的显著变化通常反映了属性的重要事件和变化。 这些包括(i)深度上的不连续、(ii)表面方向不连续、(iii)物质属性变化和(iv)场景照明变化。 边缘检测是图像处理和计算机视觉中,尤其是特征提取中的一个研究领域。

2.边缘检测算子

一阶:Roberts Cross算子,Prewitt算子,Sobel算子, Kirsch算子,罗盘算子;二阶: Marr-Hildreth,在梯度方向的二阶导数过零点,Canny算子,Laplacian算子。今天我们要讲的是基于Sobel算子的边缘检测的FPGA算法的实现。

3.Sobel算子实现

Sobel算法是像素图像边缘检测中最重要的算子之一,在机器学习、数字媒体、计算机视觉等信息科技领域起着举足轻重的作用。在技术上,它是一个离散的一阶差分算子,用来计算图像亮度函数的一阶梯度之近似值。在图像的任何一点使用此算子,将会产生该点对应的梯度矢量或是其法矢量
Soble边缘检测算法比较简,实际应用中效率比canny边缘检测效率要高,但是边缘不如Canny检测的准确,但是很多实际应用的场合,sobel边缘却是首选,尤其是对效率要求较高,而对细纹理不太关心的时候。
Soble边缘检测通常带有方向性,可以只检测竖直边缘或垂直边缘或都检测。
-1 0 +1
-2 0 +2
-1 0 +1
Sobel算子 x方向
+1 +2 +1
0 0 0
-1 -2 -1

Y方向

(i-1,j-1) ( i,j-1) (i+1,j-1)
(i-1,j) (i,j) (i+1,j)
(i-1,j+1) (i,j+1) (i+1,j+1)

原始图像P
实现步骤:

 1.Gx = P ★Sobelx   -- 原始图像与Sobel算子X方向卷积;
 2.  Gy = P★Sobely   -- 原始图像与Sobel算子Y方向卷积;
  1. 在这里插入图片描述
  2. 阈值比较形成边缘查找后的二值图像。

4.C语言实现

/* Sobel template
a00 a01 a02
a10 a11 a12
a20 a21 a22
*/
unsigned char a00, a01, a02;
unsigned char a10, a11, a12;
unsigned char a20, a21, a22;
void MySobel(IplImage* gray, IplImage* gradient)
{
CvScalar color ;
for (int i=1; i<gray->height-1; ++i)
{
for (int j=1; j<gray->width-1; ++j)
{
a00 = cvGet2D(gray, i-1, j-1).val[0];
a01 = cvGet2D(gray, i-1, j).val[0];
a02 = cvGet2D(gray, i-1, j+1).val[0];
a10 = cvGet2D(gray, i, j-1).val[0];
a11 = cvGet2D(gray, i, j).val[0];
a12 = cvGet2D(gray, i, j+1).val[0];
a20 = cvGet2D(gray, i+1, j-1).val[0];
a21 = cvGet2D(gray, i+1, j).val[0];
a22 = cvGet2D(gray, i+1, j+1).val[0];
// x方向上的近似导数  卷积运算
double ux = a20 * (1) + a10 * (2) + a00 * (1)
+ (a02 * (-1) + a12 * (-2) + a22 * (-1));
// y方向上的近似导数  卷积运算
double uy = a02 * (1) + a01 * (2) + a00 * (1)
+ a20 * (-1) + a21 * (-2) + a22 * (-1);
color.val[0] = sqrt(ux*ux + uy*uy);
cvSet2D(gradient, i, j, color);
}
}
}
//注释:该程序需要在安装Opencv软件下运行。

5.Matlab边缘检测的实现

ps=imread('lena.jpg'); %读取图像
subplot(1,3,1)
imshow(ps);
title('原图像');
ps=rgb2gray(ps);
[m,n]=size(ps); %用Sobel微分算子进行边缘检测
pa = edge(ps,'sobel');
subplot(1,3,2);
imshow(pa);
title('Sobel边缘检测得到的图像');

结果:
在这里插入图片描述

效果图

6.FPGA实现

我将在FPGA程序中注释,表示实现过程。我们使用的图像为480x272。

/* 
Filename    : Sobel.v
Compiler    : Quartus II 13.0
Description : implement Sobel Edge Detector 
Release     : 
*/

module sobel (
  input            iCLK,
  input            iRST_N,
  input      [7:0] iTHRESHOLD,
  input            iDVAL,
  input      [9:0] iDATA,
  output reg       oDVAL,
  output reg [9:0] oDATA
);
//----------------------------------------------------
// 将Sobel算子换算成有符号数(signed)
//----------------------------------------------------
// mask x
parameter X1 = 8'hff, X2 = 8'h00, X3 = 8'h01;
parameter X4 = 8'hfe, X5 = 8'h00, X6 = 8'h02;
parameter X7 = 8'hff, X8 = 8'h00, X9 = 8'h01;

// mask y
parameter Y1 = 8'h01, Y2 = 8'h02, Y3 = 8'h01;
parameter Y4 = 8'h00, Y5 = 8'h00, Y6 = 8'h00;
parameter Y7 = 8'hff, Y8 = 8'hfe, Y9 = 8'hff;

wire  [7:0] Line0;
wire  [7:0] Line1;
wire  [7:0] Line2;

wire  [17:0]  Mac_x0;
wire  [17:0]  Mac_x1;
wire  [17:0]  Mac_x2;

wire  [17:0]  Mac_y0;
wire  [17:0]  Mac_y1;
wire  [17:0]  Mac_y2;

wire  [19:0]  Pa_x;
wire  [19:0]  Pa_y;

wire  [15:0]  Abs_mag;
//---------------------------------------------
// 实现3x3矩阵原始图像 P
//---------------------------------------------
LineBuffer LineBuffer_inst (
  .clken(iDVAL),
  .clock(iCLK),
  .shiftin(iDATA[9:2]),
  .taps0x(Line0),
  .taps1x(Line1),
  .taps2x(Line2)
);
//--------------------------------------------
// Gx = P ★Sobelx
// x方向卷积运算实现
//---------------------------------------------
MAC_3 x0 (
  .aclr3(!iRST_N),
  .clock0(iCLK),
  .dataa_0(Line0),
  .datab_0(X9),
  .datab_1(X8),
  .datab_2(X7),
  .result(Mac_x0)
);

MAC_3 x1 (
  .aclr3(!iRST_N),
  .clock0(iCLK),
  .dataa_0(Line1),
  .datab_0(X6),
  .datab_1(X5),
  .datab_2(X4),
  .result(Mac_x1)
);

MAC_3 x2 (
  .aclr3(!iRST_N),
  .clock0(iCLK),
  .dataa_0(Line2),
  .datab_0(X3),
  .datab_1(X2),
  .datab_2(X1),
  .result(Mac_x2)
);
PA_3 pa0 (
  .clock(iCLK),
  .data0x(Mac_x0),
  .data1x(Mac_x1),
  .data2x(Mac_x2),
  .result(Pa_x)
);

//---------------------------------------------------
// Gy = P★Sobely
// y方向卷积运算的实现
//---------------------------------------------------
// Y
MAC_3 y0 (
  .aclr3(!iRST_N),
  .clock0(iCLK),
  .dataa_0(Line0),
  .datab_0(Y9),
  .datab_1(Y8),
  .datab_2(Y7),
  .result(Mac_y0)
);

MAC_3 y1 (
  .aclr3(!iRST_N),
  .clock0(iCLK),
  .dataa_0(Line1),
  .datab_0(Y6),
  .datab_1(Y5),
  .datab_2(Y4),
  .result(Mac_y1)
);

MAC_3 y2 (
  .aclr3(!iRST_N),
  .clock0(iCLK),
  .dataa_0(Line2),
  .datab_0(Y3),
  .datab_1(Y2),
  .datab_2(Y1),
  .result(Mac_y2)
);
PA_3 pa1 (
  .clock(iCLK),
  .data0x(Mac_y0),
  .data1x(Mac_y1),
  .data2x(Mac_y2),
  .result(Pa_y)
);
//-----------------------------------------------
// 得到G
//-----------------------------------------------
SQRT sqrt0 (
  .clk(iCLK),
  .radical(Pa_x * Pa_x + Pa_y * Pa_y),
  .q(Abs_mag)
);
//-------------------------------------------------
// 阈值比较
//-------------------------------------------------
always@(posedge iCLK, negedge iRST_N) begin
  if (!iRST_N)
    oDVAL <= 0;
  else begin
    oDVAL <= iDVAL;
    
    if (iDVAL)
      oDATA <= (Abs_mag > iTHRESHOLD) ? 0 : 1023;
    else
      oDATA <= 0;
  end
end

endmodule

IP设置
在这里插入图片描述
LineBuffer IP设置
在这里插入图片描述

LineBuffer IP的设置

在这里插入图片描述

MAC_3 IP的设置
在这里插入图片描述
PA_3 IP的设置
在这里插入图片描述
SQRT IP的设置
FPGA基于Sobel算子图像边缘检测的实现结果:
在这里插入图片描述
lena原图
在这里插入图片描述
阈值3
在这里插入图片描述
阈值5
在这里插入图片描述
阈值7

欢迎关注微信公众号:FPGA开源工作室
获取更多学习资料。
FPGA开源工作室

2018-10-03 15:46:39 u013921430 阅读数 4187

【fishing-pan:https://blog.csdn.net/u013921430 转载请注明出处】

灰度图边缘检测

   在学习图像处理时,首先接触到的就是灰度图像的边缘检测,这是图像处理最基础的也是最重要的一环,熟悉图像边缘检测有助于我们学习其他的数字图像处理方法。由于图像的边缘区域会存在明显的像素值阶跃,因此边缘检测主要是通过获得图像灰度梯度,进而通过梯度大小和变化来判断图像边缘的。
  
在这里插入图片描述
   因此,我们可以通过一阶差分;
Δfx(x,y)=f(x+1,y)f(x,y)Δfy(x,y)=f(x,y+1)f(x,y) \Delta f_{x}(x,y)=f(x+1,y)-f(x,y)\\ \Delta f_{y}(x,y)=f(x,y+1)-f(x,y)
   或者二阶差分对边缘区域进行判断;
Δfxx(x,y)=f(x+1,y)+f(x1,y)2f(x,y)Δfyy(x,y)=f(x,y+1)+f(x,y1)2f(x,y) \Delta f_{xx}(x,y)=f(x+1,y)+f(x-1,y)-2f(x,y)\\ \Delta f_{yy}(x,y)=f(x,y+1)+f(x,y-1)-2f(x,y)
   其中一阶差分可以判断边缘是否存在,二阶差分还可以根据正负号判断像素点在图像边缘亮的一侧还是暗的一侧。
   其他的边缘检测方法还包括一些梯度算子,例如Prewitt算子、Sobel算子,Canny算子,LOG边缘检测算子等,在此不做说明。

彩色图边缘检测

   RGB 图像使用三个通道存储像素信息,我们可以将这三个通道的信息看作是一个矢量,而矢量是不存在梯度的概念的,我们无法直接将上诉方法或算子直接用于RGB 图像,而且RGB图像单个通道的梯度信息又无法反映整体的梯度信息。
   在《数字图像处理》(冈萨雷斯)中提到了一种针对彩色图像的边缘检测方法,这种方法由 Di Zenzo 等人在1986年提出,下面就一起看看这种方法如何得出。

Di Zenzo’s gradient operator

   在图像多通道图像f(x,y)f(x,y) 中的某一点P(x,y)P(x,y) 处,假设其梯度方向为θ\theta
Δf=f(x+εcosθ,y+εsinθ)f(x,y) \Delta f=\left \| f(x+\varepsilon cos\theta,y+\varepsilon sin\theta) -f(x,y)\right \|
   为了便于计算,将计算绝对值换为计算平方,令
Δf2=f(x+εcosθ,y+εsinθ)f(x,y)2 \Delta f^{2}=\left \| f(x+\varepsilon cos\theta,y+\varepsilon sin\theta) -f(x,y)\right \|^{2}
   对f(x+εcosθ,y+εsinθ)f(x+\varepsilon cos\theta,y+\varepsilon sin\theta)进行二元泰勒展开;
f(x+εcosθ,y+εsinθ)=f(x,y)+i=1m(εcosθfi(x,y)x+εsinθfi(x,y)y)+onf(x,y)+i=1m(εcosθfi(x,y)x+εsinθfi(x,y)y) \begin{aligned} f(x+\varepsilon cos\theta,y+\varepsilon sin\theta)&amp;=f(x,y)+\sum_{i=1}^{m}(\varepsilon cos\theta\cdot\frac{\partial f_{i}(x,y)}{\partial x} +\varepsilon sin\theta\cdot\frac{\partial f_{i}(x,y)}{\partial y} )+o^{n}\\ &amp;\approx f(x,y)+\sum_{i=1}^{m}(\varepsilon cos\theta\cdot\frac{\partial f_{i}(x,y)}{\partial x} +\varepsilon\cdot sin\theta\cdot\frac{\partial f_{i}(x,y)}{\partial y} ) \end{aligned}
   其中mm表示图像通道数目,为了方便表述使用fix\frac{\partial f_{i}}{\partial x}代替fi(x,y)x\frac{\partial f_{i}(x,y)}{\partial x},而在求导时各个通道之间是相互独立的,则有;
Δf2i=1m(εcosθfix+εsinθfiy)2 \Delta f^{2}\approx\sum_{i=1}^{m}(\varepsilon cos\theta\cdot\frac{\partial f_{i}}{\partial x} +\varepsilon sin\theta\cdot\frac{\partial f_{i}}{\partial y} )^{2}
  重新定义一个函数G(θ)G(\theta),令
G(θ)=i=1m(εcosθfix+εsinθfiy)2=ε2(cosθ2i=1mfix2+sinθ2i=1mfiy2+2sinθcosθi=1mfixfiy) \begin{aligned} G(\theta)&amp;=\sum_{i=1}^{m}(\varepsilon cos\theta\cdot\frac{\partial f_{i}}{\partial x} +\varepsilon sin\theta\cdot\frac{\partial f_{i}}{\partial y} )^{2}\\ &amp;=\varepsilon ^{2}(cos\theta^{2}\sum_{i=1}^{m}\left \|\frac{\partial f_{i}}{\partial x}\right \|^{2}+sin\theta^{2}\sum_{i=1}^{m}\left \|\frac{\partial f_{i}}{\partial y}\right \|^{2}+2sin\theta cos\theta \sum_{i=1}^{m}\frac{\partial f_{i}}{\partial x}\frac{\partial f_{i}}{\partial y}) \end{aligned}
   进一步舍去式子中的ε\varepsilon 项,令
G(θ)=cosθ2i=1mfix2+sinθ2i=1mfiy2+2sinθcosθi=1mfixfiy G(\theta)=cos\theta^{2}\sum_{i=1}^{m}\left \|\frac{\partial f_{i}}{\partial x}\right \|^{2}+sin\theta^{2}\sum_{i=1}^{m}\left \|\frac{\partial f_{i}}{\partial y}\right \|^{2}+2sin\theta cos\theta \sum_{i=1}^{m}\frac{\partial f_{i}}{\partial x}\frac{\partial f_{i}}{\partial y}
   为了进一步方便表述;令
E=i=1mfix2;F=i=1mfiy2;H=i=1mfixfiy E=\sum_{i=1}^{m}\left \| \frac{\partial f_{i}}{\partial x} \right \|^{2}; F=\sum_{i=1}^{m}\left \|\frac{\partial f_{i}}{\partial y}\right \|^{2}; H=\sum_{i=1}^{m}\frac{\partial f_{i}}{\partial x}\frac{\partial f_{i}}{\partial y}
G(θ)=cosθ2E+sinθ2F+2sinθcosθH G(\theta)=cos\theta^{2}E+sin\theta^{2}F+2sin\theta cos\theta H
   现在θ\theta成为了式子中唯一的变量,再回到边缘的定义上,边缘的方向是图像像素梯度最大的方向。也就是说梯度的方向θmax\theta_{max} 会使G(θ)G(\theta) 取最大值,则;
θmax=G(θ)argmax \theta_{max}=\underset{argmax}{G(\theta )}
   对G(θ)G(\theta) 进行求导;
G(θ)=Esin2θ+Fsin2θ+2Hcos2θ G(\theta )^{&#x27;}=-Esin2\theta +F sin2\theta+2 H cos2\theta
   令G(θ)=0G(\theta )^{&#x27;}=0,得;
tan 2θmax=2HEFθmax=12arctan(2HEF+kπ) tan ~2\theta_{max} =\frac{2H}{E-F}\\ \theta_{max}=\frac{1}{2}arctan(\frac{2H}{E-F}+k\pi)
   很明显G(θ)G(\theta ) 是一个以π\pi 为周期的周期函数,如果只考虑区间[0,π)\left [ 0 ,\pi\right ),且θmax\theta_{max} 落到该区间内,则会有另一个让G(θ)G(\theta )取极值的解也落在该区域内,这个值是θmax+π2\theta_{max}+ \frac{\pi}{2}或者θmaxπ2\theta_{max}-\frac{\pi}{2}。但是不论如何这两个解有一个让G(θ)G(\theta )取极大值,另一个让其取极小值,两个角度相差 90°。
  
   说到这里大家应该都明白了,两个角度对应的方向是相互垂直的,一个是垂直于边缘的方向,也就是边缘的法向,此时梯度取最大值。另一个是平行于边缘的方向,是边缘的切向,此时梯度取极小值。
  

RGB图像的边缘检测

   在RGB图像中,m=3m=3,再令;
u=Rxr+Gxg+Bxb \overset{\rightarrow }{u}=\frac{\partial R}{\partial x}\overset{\rightarrow }{r}+\frac{\partial G}{\partial x}\overset{\rightarrow }{g}+\frac{\partial B}{\partial x}\overset{\rightarrow }{b}
v=Ryr+Gyg+Byb \overset{\rightarrow }{v}=\frac{\partial R}{\partial y}\overset{\rightarrow }{r}+\frac{\partial G}{\partial y}\overset{\rightarrow }{g}+\frac{\partial B}{\partial y}\overset{\rightarrow }{b}
   其中r\overset{\rightarrow }{r}g\overset{\rightarrow }{g}b\overset{\rightarrow }{b}分别代表不同颜色分量的单位向量,则
gxx=E=uTu=Rx2+Gx2+Bx2 g_{xx}=E=\overset{\rightarrow }{u}^{\tiny{T}}\overset{\rightarrow }{u}=\left \| \frac{\partial R}{\partial x} \right \|^{2}+\left \| \frac{\partial G}{\partial x} \right \|^{2}+\left \| \frac{\partial B}{\partial x} \right \|^{2}
gyy=F=vTv=Ry2+Gy2+By2 g_{yy}=F=\overset{\rightarrow }{v}^{\tiny{T}}\overset{\rightarrow }{v}=\left \| \frac{\partial R}{\partial y} \right \|^{2}+\left \| \frac{\partial G}{\partial y} \right \|^{2}+\left \| \frac{\partial B}{\partial y} \right \|^{2}
gxy=H=uTv=RxRy+GxGy+BxBy g_{xy}=H=\overset{\rightarrow }{u}^{\tiny{T}}\overset{\rightarrow }{v}=\frac{\partial R}{\partial x}\frac{\partial R}{\partial y}+\frac{\partial G}{\partial x}\frac{\partial G}{\partial y}+\frac{\partial B}{\partial x}\frac{\partial B}{\partial y}
   在利用Di Zenzo 提出的方法求得θmax\theta_{max} 后,将以上符号带入到G(θ)G(\theta),可以计算出像素点梯度大小为
G(θ)={12[(gxx+gyy)+(gxxgyy)cos 2θmax+2gxysin 2θmax]}12 G(\theta)=\left \{ \frac{1}{2}\left [ (g_{xx}+g_{yy}) +(g_{xx}-g_{yy})cos~2\theta_{max} +2g_{xy}sin~2\theta_{max} \right ]\right \}^{\frac{1}{2}}
   进而可以根据梯度大小进行边缘检测。

参考资料

  1. S. Di Zenzo, A note on the gradient of a multi-image, Computer Vision, Graphics, and Image Processing 33 (1) (1986) 116–125.
  2. 《数字图像处理》 (冈萨雷斯)
2016-03-12 12:41:39 HorizontalView 阅读数 1017

图像处理笔记 —— 边缘检测



一. 边缘检测的基本概念

  1. 图像的边缘相当于二元函数的梯度幅度比较大的位置点,而且梯度还能反映图像局部变化最快的方向,在某指定的方向上,把梯度模超过指定阈值的位置数据记录下来,就是沿该方向的图像边缘。

  2. 正确地描绘边缘应当包括两个方面:边缘点的位置以及局部的边缘方向。边缘一般是指灰度变化率最大的位置。图像边缘的成因:

    • 图像灰度在表面法向变化的不连续性造成的边缘。

    • 图像对象 在空间上的深度不一致形成的边缘。

    • 在光顺表面上由于颜色的不一致性形成的边缘

    • 物体的光影造成的边缘

    常用术语

    • 阶跃边:在剖面上看,图像灰度在奇异点,两侧由明显的差异

    • 高台边:图像灰度除去奇异点位置之外保持平坦,只在奇异点附近突然跃起或者降低。

    • 边缘点:图像灰度的显著变化点,又称奇异点

    • 边缘段:边缘的坐标及其方向

    • 边缘检测:从图像中检测边缘点和边缘段,并且描述边缘方向的过程。在边缘检测中,存在两种类型的错误:漏检了边缘点和非边缘点误认为边缘点。

    • 轮廓:图像边缘的列表

    • 边缘连接:将边缘列表元素有序化的过程

    • 边缘跟踪:在图像或者图像序列中,确定轮廓的过程。


二. 基于微分算子的边缘检测原理

  1. 微分算子的两个重要性质,由于这两个性质,在使用微分算子之前,应该采用预处理抗噪声技术。

    • 定域性: 即局部性,每一点的导数只与函数在该点邻近的信息有关

    • 敏感性:对局部的函数值变化敏感,抗噪性差

  2. 微分算子的定域性:
    设函数 f(x,y),(x,y)Ω 有足够高阶的连续可导性,则,Dm,n:f(x,y)Dm,nf(x,y) 是有限阶偏导数,
    Dm,nf(x,y)=m+nf(x,y)xmyn
    那么要计算一点(a,b)处的导数Dm,n:f(a,b) ,只需要f(x,y)(a,b) 的一个小邻域B={(x,y)|(xa)2+(yb)2ε} 中的信息{f(x,y)|(x,y)B} 即可。
  3. 边缘检测技术的基本步骤:

    • 将相应的微分算子简化为离散差分格式,进而简化成模板(记为 T)。

    • 利用模板对图像 f(x,y),获得模板作用后的结果 Tf(m,n)

    • 提出阈值 h,在采用一阶微分算子情形,记录下高于某个阈值 h 的位置坐标 sh={(m,n)||Tf(m,n)|h} ,采用二阶微分算子时,一般时对某个阈值 ε>0 确立 sh={(m,n)||Tf(m,n)|ε}

    • 对集合 sh 进行整理,同时调整阈值 h

  4. 阈值 h 很重要,如果过小,则结果呈现为到处都是边缘点,如果过大,则使边缘点特别稀疏。

三. 算子

  1. 方向导数与梯度之间具有以下关系εf(a,b)=<ε,f(a,b)> 并且 0εf(a,b)εf(a,b) 假定 f0,则该不等式右端取等号的充分必要条件:εf 同方向,即 ε=kf; 该不等式左端取等号的充分必要条件:εf 垂直。这表明,图像f(x,y)在任意一点 (a,b) 处,变化最快的方向就是梯度方向,同时,变化最慢的方向是与梯度相互垂直的方向。
  2. Roberts 差分公式:
    • |T1f(i,j)|=|f(i,j)f(i+1,j+1)|2+|f(i+1,j)f(i,j+1)|2
    • |T2f(i,j)|=|f(i,j)f(i+1,j+1)|+|f(i+1,j)f(i,j+1)|
  3. 若只对 Rxf进行判决,则选择的是与 x 方向垂直的边缘。若只对 Ryf 进行判决,则选择的是与 y 方向垂直的边缘。
  4. Sobel 算子差分公式为: |f(i,j)|=|Sx.f|+|Sy.f| 矩阵之间使用的是点乘)

    • 其中 f

      f(i1,j1)f(i1,1)f(i1,j+1)f(i,j1)f(i,j)f(i,j+1)f(i+1,j1)f(i+1,j)f(i+1,j+1)

    • Sx

      121000121

    • Sx

      101202101

  5. Prewitt算子的模板为

    • Sx

      111000111

    • Sx

      101101101

  6. Kirsch算子总共有8种模板,这里我们仅列出一种

    • Sx

      533503533

    • Sx

      333503553

  7. Robinson算子总共有8种模板,这里我们仅列出一种

    • Sx

      101202121

    • Sx

      012101210


四. 二阶微分算子

  1. 二阶微分算子可以找到一个一致的阈值进行选择来进行那个边缘检测定位。
  2. 在二维边缘检测中,经常采用 Laplace 微分算子
    f(x,y)=2f2x+2f2y
  3. 现实生活中采用前后向差分策略
    • 2fx2f(i+1,j)2f(i,j)+f(i1,j)
    • 2fy2f(i,j+1)2f(i,j)+f(i,j1)

相关源码及图像

%%matlab 
clc;clear all;close all;

%读取图片
img = imread('4.1.06.tiff');
figure('NumberTitle','off','Name','边缘检测');
subplot(331);imshow(img);title('原图');

%转换成灰度图
img_g = rgb2gray(img);
subplot(332);imshow(img_g);title('灰度图');

%基于 Roberts 算子的一阶微分算子检测
roberts_img = edge(img_g,'roberts',0.1);
subplot(333);imshow(roberts_img);title('Roberts 算子');

%基于 Sobel 算子的一阶微分算子检测
roberts_img = edge(img_g,'sobel',0.1);
subplot(334);imshow(roberts_img);title('Sobel算子');

%基于 Prewitt 算子的一阶微分算子检测
roberts_img = edge(img_g,'prewitt',0.1);
subplot(335);imshow(roberts_img);title('Prewitt 算子');

%二阶微分算子检测
roberts_img = edge(img_g,'log',0.01);
subplot(336);imshow(roberts_img);title('二阶微分算子');

%基于 zerocross 算子的一阶微分算子检测
roberts_img = edge(img_g,'zerocross',0.01);
subplot(337);imshow(roberts_img);title('zerocross 算子');

这里写图片描述


参考文献

[1] 王桥. 数字图像处理[M]. 北京:科学出版社, 2009. 113-131

2018-04-18 20:02:11 Gentleman_Qin 阅读数 7793

高光谱图像处理之目标检测技术

一、高光谱图像处理之目标检测


1、高光谱图像目标检测的发展趋势和研究现状:

   20世纪80年代末,美国的一些研究机构开始利用高光谱图像数据进行目标检测方面的研究。自上世纪九十年代,国外出现了进行高光谱图像目标检测算法理论研究的研究组。由Reed和Yu提出了基于广义似然比检验的恒虚警RX 检测器(RXD)。Chang课题组提出了基于正交子空间投影的OSP检测方法, Harsanyi提出了基于约束能量最小化的CEM算法。未来高光谱目标检测的发展将会越来越重视实用性,算法的性能将会进一步提高,同时更也加适合使用FPGA硬件对其进行加速从而具有更高的实际应用价值。

2、高光谱目标检测技术的应用范围:

   高光谱目标检测具有较强的实用性,可应用于公共安全、环境检测、城市规划、食品卫生、地质岩矿的识别等众多方面。


                         图1.1 高光谱图像数据结构

3、高光谱图像数据的特点:

   高光谱图像数据“图谱合一”,具有丰富的光谱维信息。高光谱图像数据共有三个维度,如图2.1所示,其中,图像空间维信息x、y用于表示物体的实际空间分布,而光谱波段L用于表示每个像素的光谱属性。

4、高光谱图像目标检测原理:

   高光谱图像的各波段在成像范围内都是连续成像,因此高光谱的光谱曲线一般是平滑的、连续的曲线。高光谱图像的波段L中涵盖了物质的光谱信息,而每种物质的光谱信息都不一样,我们可以利用图像像素的光谱波段L所包含的特定的光谱信息来判断该像素所代表的特定的物质种类。

5、高光谱图像目标检测流程:


                    图1.2 高光谱图像目标检测流程图

   如图所示,对于拍摄得到的原始高光谱图像数据,需要先对数据进行预处理,包括数据格式化、无用数据剔除以及亮度到反射率的转化等。同时,对于遥感仪拍摄的高光谱图像还需要进行辐射校正,在目标检测前,应对数据进行调整,包括数据归一化等。最后根据已知的先验信息选择相应的目标检测算法进行检测。

6、现场可编程门阵列(FPGA)在高光谱图像处理中的应用:

   一个FPGA可以大致定义为一系列互连的逻辑块。这些器件的主要优点之一是,为了实现不同的组合或时序逻辑功能,可以根据需要多次重新配置逻辑块及其互连。这一特性为FPGA提供了软件和硬件系统的优势。因为FPGA具有比专用集成电路(ASIC)更多的灵活性和更短的开发时间,更接近GPU提供的性能,但同时功耗要比GPU低得多。FPGA供应商通过改进FPGA体系结构(包括优化的硬件模块)并利用最新的硅技术不断提升FPGA的功耗和能效。可重构性、低功耗的特点以及FPGA对空间电离辐射耐受性的提高,这些因素已经使得FPGA成为目前板载高光谱遥感的最佳选择之一。


二、算法分类、比较和选择:


1、方法分类:

   高光谱目标检测方法按照先验信息是否已知分为监督方法和非监督方法。前者用于目标光谱已知的情况下,利用目标光谱与图像像元光谱进行匹配,从而完成目标检测,比如CEM算法、OSP算法;后者多用于异常目标检测,一般不需要目标和背景的先验信息,根据高光谱图像数据获取目标检测所需要的数据,然后根据数据的大小来判断是否为异常目标,比如RXD算法。

2、CEM、OSP、RXD算法的区别:

(1)CEM(Constrained Energy Minimization)算法:

   CEM算法主要思想是设计一个FIR线性滤波器,使得在满足约束条件式的情况下滤波器的输出能量最小。该算法不需要图像的背景信息,只需要知道要检测的先验光谱信息(目标向量)即可,具体方法是通过高光谱图像数据和先验已知的待检测目标确定一个滤波向量,让图像经过该滤波向量即可得到检测结果,其中滤波向量的作用是滤除图像中的非目标像素,让感兴趣的目标能够通过,同时抑制由其他信号带来的滤波器输出能量。

(2)OSP(Orthogonal Subspace Projection)算法:

   OSP算法与CEM算法相比,最大区别在于不仅需要目标的先验知识,还需要图像中背景的先验知识,但在实际中中这些先验信息很难全部得到。在高光谱检测中我们一般用其来检测异常。该算法需要前提条件:图像信息、目标像元、非目标像元(异常目标)信息。

(3)RXD(Reed-XiaoliDetector)算法:

   RXD算法是异常目标检测领域中最基础的算法,不需要目标光谱的先验知识,而是基于背景服从多元正态分布的假设,通过检测与背景分布中心相比属于异常像元,并在这些感兴趣区域进一步查找可能存在的目标。该算法主要针对的是小目标检测问题。

(4)确定所采用的目标检测算法:

   由于我们的应用场景多为有特定目标的目标检测,CEM正是针对未知场景中可能存在的特定目标的检测,只需要知道目标的光谱信息即可,而RXD算法适应于对特定场景的异常(未知目标)检测,而OSP算法除了需要已知目标光谱还需要背景信息。综上,我决定采用CEM算法进行高光谱目标检测的实现。


三、CEM算法分析:


1、算法步骤

(1)对高光谱图像进行预处理,得到二维化和归一化后的数据r(L*N);

(2)根据图像数据r,求得图像的自相关矩阵:

               

(3)确定目标光谱向量dd大小为L*1);

(4)根据公式:


                                   

   设计FIR线性滤波器:


(5)将归一化后的数据经过FIR滤波器,得到输出信号y

           

2、问题分解:

   CEM算法的实现可分为三部分:自相关矩阵、矩阵求逆、线性FIR滤波器。


图 1.3 CEM算法分解流程图

   如图所示,首先根据高光谱图像r求得自相关矩阵,再利用矩阵求逆模块求得自相关矩阵的逆矩阵,结合从光谱库获取到的目标向量的先验信息求得FIR滤波器的滤波向量,最后将高光谱图像r通过FIR滤波器即可得到最终的检测结果。


四、CEM算法实现:


1、算法流程:

   在MATLAB和C语言中实现CEM算法的具体流程如图所示,因为语言特性是串行执行命令,所以在编写程序时与硬件设计比较更加直接明了。

图1.4 CEM算法流程图

2、数据预处理:

   对前期得到的高光谱图像在MATLAB平台上进行预处理。这一过程主要对原始的200*200*189大小的高光谱图像进行操作:

(1)二维化:通过调用MATLAB里面的reshape()函数实现。

(2)归一化:采用“min-max”方法。 具体步骤是先找到数据的最大值
(maxA)和最小值(minA), 通过 y = (x-minA) /(maxA-minA) 计算的 y 即为归一化后的数据。高光谱数据量大,也造成了数据的存储比较困难,数据在程序中的存储等处理

3、CEM算法的MATLAB实现:

   首先,根据 2 对高光谱图像数据进行预处理,而后求其自相关矩阵。高光谱图像数据的大小为200*200*189,其中200*200为像素数,189为波段数。

              

   首先用reshape()函数和transpose()函数将200*200*189的三维高光谱图像数据转为189*40000的二维矩阵,再用自己写的Normalize()函数对矩阵数据用“min-max”方法进行归一化,然后根据公式式利用MATLAB中的矩阵求转置函数transpose()和矩阵相乘操作得到一个189*189的矩阵,最后对189*189的矩阵除以像素数N即可得到自相关矩阵。

最后,根据先验已知的目标向量和自相关矩阵R的逆矩阵(在MATLAB中矩阵R求逆即为(1/R))求得FIR滤波器,将高光谱图像模拟数据通过FIR滤波器,即可得到最终的检测结果。(本次测试中选取的目标向量 d 是模拟图像中的 C 物质)

4、CEM算法的C语言实现:

   在Visual Studio平台上完成CEM算法C程序的编写。相比于MATLAB实现,CEM算法的C语言实现主要难点在于矩阵运算。因为在C语言中,矩阵的转置、相乘、求逆等操作均没有现有的函数,需要编写相应的函数。在实现矩阵基本运算过程中,通过动态分配内存运用二维指针传递参数,完成矩阵加减法运算、矩阵转置运算、矩阵相乘运算以及矩阵求逆运算,这样可以节省存储空间,使用完后释放空间即可。

CEM算法的C语言实现主要包括如下步骤:

(1)在MATLAB中,将归一化后的高光谱图像数据转为189*40000的二维形式,保存为CEM.mat;

(2)编写矩阵初始化、矩阵转置、矩阵相乘、矩阵求逆、内存释放函数;

(3)编写main函数,读CEM.mat,调用上述函数,求得FIR滤波器,进而获得输出信号y,将其写入CEM.txt;

(4)在MATLAB中显示CEM.txt中数据所代表的图像。

5、难点解决:

(1)矩阵转置

   矩阵转置函数的输入是大小为m*n的矩阵A,输出是大小n*m为的矩阵B。采用嵌套的for循环分别遍历矩阵的行和列,将输入的二维矩阵A按列读出,重新按行写入新矩阵B中,即可实现矩阵的转置。

(2)矩阵相乘

   矩阵相乘函数的输入是大小m*n为的矩阵A和n*k的矩阵B,输出是大小为m*k的矩阵C。采用嵌套的for循环分别遍历矩阵的行和列,首先将矩阵A按行读出,矩阵B按列读出,然后将读出的矩阵A的第i行和矩阵B的第j列对应位相乘求和,将计算的结果写入矩阵C第i行第j列元素中。

(3)矩阵求逆——QR分解求逆

   QR分解求逆的原理是:对于可逆矩阵A,首先利用QR分解将A矩阵分解为Q矩阵和R矩阵。即A=QR,其中Q是正交矩阵,R是上三角矩阵。然后将公式左右同时求逆,可以求得A的逆矩阵,其中Q矩阵的逆矩阵和转置矩阵相同,R求逆有特定的公式。QR分解求逆的运算较为简单,且数据稳定度较高,可以得到误差较小的求逆结果。

6、算法实现结果和对比分析:

1)MATLAB和C语言实现CEM算法的检测结果如图所示:

       

      (a)               (b)            (c)

图1.6 原图CEM算法检测结果:(a)原图(每一行是一种物质)(b)MATLAB检测C物质(c)C语言检测C物质

 

通过对MATLAB和C语言实现结果进行对比,发现二者均成功实现了CEM算法并完成了对目标的检测,观察检测结果基本无差异。

######### 转载请注明出处:https://blog.csdn.net/Gentleman_Qin/article/details/79992314 #########


没有更多推荐了,返回首页