图像处理-线性滤波-1 基础(相关算子、卷积算子、边缘效应)
这里讨论利用输入图像中像素的小邻域来产生输出图像的方法,在信号处理中这种方法称为滤波(filtering)。其中,最常用的是线性滤波:输出像素是输入邻域像素的加权和。
1.相关算子(Correlation Operator)
步骤:
例:
A =[17
24 1 8 15 h=[8 1 6
23 5 7 14 16 3 5 7
4 6 13 20 22 4 9 2]
10 12 19 21 3
11 18 25 2 9] Matlab 函数:imfilter(A,h)
2.卷积算子(Convolution)
步骤:
Matlab 函数:Matlab 函数:imfilter(A,h,'conv')%imfilter默认是相关算子,因此当进行卷积计算时需要传入参数'conv'
3.边缘效应
当对图像边缘的进行滤波时,核的一部分会位于图像边缘外面。
常用的策略包括:
1)使用常数填充:imfilter默认用0填充,这会造成处理后的图像边缘是黑色的。
2)复制边缘像素:I3 = imfilter(I,h,'replicate');
4.常用滤波
fspecial函数可以生成几种定义好的滤波器的相关算子的核。
例:unsharp masking 滤波
12345I = imread(
'moon.tif'
);
h = fspecial(
'unsharp'
);
I2 = imfilter(I,h);
imshow(I), title(
'OriginalImage'
)
figure, imshow(I2), title(
'FilteredImage'
)
图像处理-线性滤波-2 图像微分(1、2阶导数和拉普拉斯算子)
更复杂些的滤波算子一般是先利用高斯滤波来平滑,然后计算其1阶和2阶微分。由于它们滤除高频和低频,因此称为带通滤波器(band-passfilters)。
在介绍具体的带通滤波器前,先介绍必备的图像微分知识。
1 一阶导数
对于离散情况(图像),其导数必须用差分方差来近似,有
,前向差分forwarddifferencing
1)前向差分的Matlab实现
123456789101112131415161718192021222324252627function dimg = mipforwarddiff(img,direction)
%MIPFORWARDDIFF
Finitedifferencecalculations %
%
DIMG= MIPFORWARDDIFF(IMG,DIRECTION) %
%
Calculates theforward-difference for
agiven direction
%
IMG :input image %
DIRECTION: 'dx'
or
'dy'
%
DIMG :resultant image %
%
Seealso MIPCENTRALDIFF MIPBACKWARDDIFFMIPSECONDDERIV %
MIPSECONDPARTIALDERIV
%
OmerDemirkaya, Musa Asyali, Prasana Shaoo, ...9/1/06 %
MedicalImage Processing Toolbox
imgPad = padarray(img,[1 1],
'symmetric'
,
'both'
);%将原图像的边界扩展
[row,col] = size(imgPad);
dimg = zeros(row,col);
switch
(direction)
case
'dx'
,
dimg(:,1:col-1)=imgPad(:,2:col)-imgPad(:,1:col-1);%x方向差分计算,
case
'dy'
,
dimg(1:row-1,:)=imgPad(2:row,:)-imgPad(1:row-1,:);
otherwise, disp(
'Directionis unknown'
);
end;
dimg = dimg(2:end-1,2:end-1);
2)中心差分的Matlab实现
12345678910111213141516171819202122232425262728function dimg = mipcentraldiff(img,direction)
%MIPCENTRALDIFF
Finitedifferencecalculations %
%
DIMG= MIPCENTRALDIFF(IMG,DIRECTION) %
%
Calculates thecentral-difference for
agiven direction
%
IMG :input image %
DIRECTION: 'dx'
or
'dy'
%
DIMG :resultant image %
%
Seealso MIPFORWARDDIFF MIPBACKWARDDIFFMIPSECONDDERIV %
MIPSECONDPARTIALDERIV
%
OmerDemirkaya, Musa Asyali, Prasana Shaoo, ...9/1/06 %
MedicalImage Processing Toolbox
img = padarray(img,[1 1],
'symmetric'
,
'both'
);
[row,col] = size(img);
dimg = zeros(row,col);
switch
(direction)
case
'dx'
,
dimg(:,2:col-1)= (img(:,3:col)-img(:,1:col-2))/2;
case
'dy'
,
dimg(2:row-1,:)= (img(3:row,:)-img(1:row-2,:))/2;
otherwise,
disp(
'Directionis unknown'
);
end
dimg = dimg(2:end-1,2:end-1);
1
实例:技术图像x方向导数
12I = imread(
'coins.png'
);figure; imshow(I);
Id = mipforwarddiff(I,
'dx'
);figure, imshow(Id);
原图像 x方向1阶导数
2 图像梯度(Image Gradient)
图像I的梯度定义为
Matlab函数
1)gradient:梯度计算
2)quiver:以箭头形状绘制梯度。注意放大下面最右侧图可看到箭头,由于这里计算横竖两个方向的梯度,因此箭头方向都是水平或垂直的。
实例:仍采用上面的原始图像
12345I =
double
(imread(
'coins.png'
));
[dx,dy]=gradient(I);
magnitudeI=sqrt(dx.^2+dy.^2);
figure;imagesc(magnitudeI);colormap(gray);%梯度幅值
hold on;quiver(dx,dy);%叠加梯度方向
3 二阶导数
3.1 普拉斯算子(laplacian operator)
3.1.2 概念
拉普拉斯算子是n维欧式空间的一个二阶微分算子。它定义为两个梯度向量算子的内积
其在二维空间上的公式为:
对于1维离散情况,其二阶导数变为二阶差分
2)因此,二阶差分为
对于2维离散情况(图像),拉普拉斯算子是2个维上二阶差分的和(见式3.3),其公式为:
上式对应的卷积核为
常用的拉普拉斯核有:
3.1.2 应用
拉普拉斯算子会突出像素值快速变化的区域,因此常用于边缘检测。
Matlab里有两个函数
1)del2
2)fspecial:图像处理中一般利用Matlab函数fspecial
h = fspecial('laplacian', alpha) returns a 3-by-3 filterapproximating the shape of the two-dimensional Laplacianoperator.
The parameter alpha controls the shape of the Laplacian and must bein the range 0.0 to 1.0. The default value for alpha is0.2.
3.1.3 资源
http://fourier.eng.hmc.edu/e161/lectures/gradient/node8.html
http://homepages.inf.ed.ac.uk/rbf/HIPR2/log.htm
sift算法
尺度不变特征转换(Scale-invariant featuretransform
Sift算法就是用不同尺度(标准差)的高斯函数对图像进行平滑,然后比较平滑后图像的差别,
差别大的像素就是特征明显的点。
sift可以同时处理亮度,平移,旋转,尺度的变化,利用特征点来提取特征描述符,最后在特征描述符之间寻找匹配
五个步骤
1构建尺度空间,检测极值点,获得尺度不变性
2特征点过滤并进行经确定位,剔除不稳定的特征点
3 在特征点处提取特征描述符,为特征点分配方向直
4声称特征描述子,利用特征描述符寻找匹配点
5计算变换参数
当2幅图像的sift特征向量生成以后,下一步就可以采用关键点特征向量的欧式距离来作为2幅图像中关键点的相似性判定量度
尺度空间:
尺度就是受delta这个参数控制的表示
而不同的L(x,y,delta)就构成了尺度空间,实际上具体计算的时候即使连续的高斯函数,都要被离散为矩阵来和数字图像进行卷积操作
L(x,y,delta)=G(x,y,e)*i(x,y)
尺度空间=原始图像(卷积)一个可变尺度的2维高斯函数G(x,y,e)
G(x,y,e) = [1/2*pi*e^2] * exp[ -(x^2 +y^2)/2e^2]
为了更有效的在尺度空间检测到稳定的关键点,提出了高斯差分尺度空间,利用不同尺度的高斯差分核与原始图像i(x,y)卷积生成
D(x,y,e)=(G(x,y,ke)-G(x,y,e))*i(x,y)
=L(x,y,ke)-L(x,y,e)
(为避免遍历每个像素点)
高斯卷积:
在组建一组尺度空间后,再组建下一组尺度空间,对上一组尺度空间的最后一幅图像进行二分之一采样,得到下一组尺度空间的第一幅图像,然后进行像建立第一组尺度空间那样的操作,得到第二组尺度空间,公式定义为
高斯差分
图像处理之卷积概念
我们来看一下一维卷积的概念.
连续空间的卷积定义是 f(x)与g(x)的卷积是 f(t-x)g(x)在t从负无穷到正无穷的积分值.t-x要在f(x)定义域内,所以看上去很大的积分实际上还是在一定范围的.
实际的过程就是f(x)先做一个Y轴的反转,然后再沿X轴平移t就是f(t-x),然后再把g(x)拿来,两者乘积的值再积分.想象一下如果g(x)或者f(x)是个单位的阶越函数.那么就是f(t-x)与g(x)相交部分的面积.这就是卷积了.
把积分符号换成求和就是离散空间的卷积定义了.
那么在图像中卷积卷积地是什么意思呢,就是图像f(x),模板g(x),然后将模版g(x)在模版中移动,每到一个位置,就把f(x)与g(x)的定义域相交的元素进行乘积并且求和,得出新的图像一点,就是被卷积后的图像.模版又称为卷积核.卷积核做一个矩阵的形状.
卷积定义上是线性系统分析经常用到的.线性系统就是一个系统的输入和输出的关系是线性关系.就是说整个系统可以分解成N多的无关独立变化,整个系统就是这些变化的累加.
如 x1->y1, x2->y2; 那么A*x1 + B*x2-> A*y1 + B*y2 这就是线性系统. 表示一个线性系统可以用积分的形式 如 Y =Sf(t,x)g(x)dt S表示积分符号,就是f(t,x)表示的是A B之类的线性系数.
看上去很像卷积呀,,对如果f(t,x) = F(t-x)不就是了吗.从f(t,x)变成F(t-x)实际上是说明f(t,x)是个线性移不变,就是说 变量的差不变化的时候,那么函数的值不变化.实际上说明一个事情就是说线性移不变系统的输出可以通过输入和表示系统线性特征的函数卷积得到.
古人曰:”说一堆大道理不如举一个好例子”,冲量这一物理现象很能说明”冲击函数”。在t时间内对一物体作用F的力,我们可以让作用时间t很小,作用力F很大,但让Ft的乘积不变,即冲量不变。于是在用t做横坐标、F做纵坐标的坐标系中,就如同一个面积不变的长方形,底边被挤的窄窄的,高度被挤的高高的,在数学中它可以被挤到无限高,但即使它无限瘦、无限高、但它仍然保持面积不变(它没有被挤没!),为了证实它的存在,可以对它进行积分,积分就是求面积嘛!于是”卷积”这个数学怪物就这样诞生了。说它是数学怪物是因为追求完美的数学家始终在头脑中转不过来弯,一个能瘦到无限小的家伙,竟能在积分中占有一席之地,必须将这个细高挑清除数学界。但物理学家、工程师们确非常喜欢它,因为它解决了很多当时数学家解决不了的实际问题。最终追求完美的数学家终于想通了,数学是来源于实际的,并最终服务于实际才是真。于是,他们为它量身定做了一套运作规律。于是,妈呀!你我都感觉眩晕的卷积分产生了。
有一个七品县令,喜欢用打板子来惩戒那些市井无赖,而且有个惯例:如果没犯大罪,只打一板,释放回家,以示爱民如子。
有一个无赖,想出人头地却没啥指望,心想:既然扬不了善名,出恶名也成啊。怎么出恶名?炒作呗!怎么炒作?找名人呀!他自然想到了他的行政长官——县令。
无赖于是光天化日之下,站在县衙门前撒了一泡尿,后果是可想而知地,自然被请进大堂挨了一板子,然后昂首挺胸回家,躺了一天,嘿!身上啥事也没有!第二天如法炮制,全然不顾行政长管的仁慈和衙门的体面,第三天、第四天……每天去县衙门领一个板子回来,还喜气洋洋地,坚持一个月之久!这无赖的名气已经和衙门口的臭气一样,传遍八方了!
县令大人噤着鼻子,呆呆地盯着案子上的惊堂木,拧着眉头思考一个问题:这三十个大板子怎么不好使捏?……想当初,本老爷金榜题名时,数学可是得了满分,今天好歹要解决这个问题:
——人(系统!)挨板子(脉冲!)以后,会有什么表现(输出!)?
——费话,疼呗!
——我问的是:会有什么表现?
——看疼到啥程度。像这无赖的体格,每天挨一个板子啥事都不会有,连哼一下都不可能,你也看到他那得意洋洋的嘴脸了(输出0);如果一次连揍他十个板子,他可能会皱皱眉头,咬咬牙,硬挺着不哼
(输出1);揍到二十个板子,他会疼得脸部扭曲,象猪似地哼哼(输出3);揍到三十个板子,他可能会象驴似地嚎叫,一把鼻涕一把泪地求你饶他一命(输出5);揍到四十个板子,他会大小便失禁,勉
强哼出声来(输出1);揍到五十个板子,他连哼一下都不可能(输出0)——死啦!
县令铺开坐标纸,以打板子的个数作为X轴,以哼哼的程度(输出)为Y轴,绘制了一条曲线:
——呜呼呀!这曲线象一座高山,弄不懂弄不懂。为啥那个无赖连挨了三十天大板却不喊绕命呀?
——呵呵,你打一次的时间间隔(Δτ=24小时)太长了,所以那个无赖承受的痛苦程度一天一利索,没有叠加,始终是一个常数;如果缩短打板子的时间间隔(建议Δτ=0.5秒),那他的痛苦程度可就迅速叠加了;等到这无赖挨三十个大板(t=30)时,痛苦程度达到了他能喊叫的极限,会收到最好的惩戒效果,再多打就显示不出您的仁慈了。
——还是不太明白,时间间隔小,为什么痛苦程度会叠加呢?
——这与人(线性时不变系统)对板子(脉冲、输入、激励)的响应有关。什么是响应?人挨一个板子后,疼痛的感觉会在一天(假设的,因人而异)内慢慢消失(衰减),而不可能突然消失。这样一来,只要打板子的时间间隔很小,每一个板子引起的疼痛都来不及完全衰减,都会对最终的痛苦程度有不同的贡献:
t个大板子造成的痛苦程度=Σ(第τ个大板子引起的痛苦*衰减系数)
[衰减系数是(t-τ)的函数,仔细品味]
数学表达为:y(t)=∫T(τ)H(t-τ)
——拿人的痛苦来说卷积的事,太残忍了。除了人以外,其他事物也符合这条规律吗?
——呵呵,县令大人毕竟仁慈。其实除人之外,很多事情也遵循此道。好好想一想,铁丝为什么弯曲一次不折,快速弯曲多次却会轻易折掉呢?
——恩,一时还弄不清,容本官慢慢想来——但有一点是明确地——来人啊,将撒尿的那个无赖抓来,狠打40大板!
卷积(convolution, 另一个通用名称是德文的Faltung)的名称由来,是在于当初定义它时,定义成integ(f1(v)*f2(t-v))dv,积分区间在0到t之间。举个简单的例子,大家可以看到,为什么叫”卷积”了。比方说在(0,100)间积分,用简单的辛普生积分公式,积分区间分成100等分,那么看到的是f1(0)和f2(100)相乘,f1(1)和f2(99)相乘,f1(2)和f2(98)相乘,……… 等等等等,就象是在坐标轴上回卷一样。所以人们就叫它”回卷积分”,或者”卷积”了。
为了理解”卷积”的物理意义,不妨将那个问题”相当于它的时域的信号与系统的单位脉冲响应的卷积”略作变化。这个变化纯粹是为了方便表达和理解,不影响任何其它方面。将这个问题表述成这样一个问题:一个信号通过一个系统,系统的响应是频率响应或波谱响应,且看如何理解卷积的物理意义。
假设信号函数为f,响应函数为g。f不仅是时间的函数(信号时有时无),还是频率的函数(就算在某一固定时刻,还有的地方大有的地方小);g也是时间的函数(有时候有反应,有时候没反应),同时也是频率的函数(不同的波长其响应程度不一样)。那我们要看某一时刻t 的响应信号,该怎么办呢?
这就需要卷积了。
要看某一时刻 t 的响应信号,自然是看下面两点:
1。你信号来的时候正赶上人家”系统”的响应时间段吗?
2。就算赶上系统响应时间段,响应有多少?
响 应不响应主要是看 f 和 g两个函数有没有交叠;响应强度的大小不仅取决于所给的信号的强弱,还取决于在某频率处对单位强度响应率。响应强度是信号强弱和对单位强度信号响应率的乘积。”交叠”体现在f(t1)和g(t-t1)上,g之所以是”(t-t1)”就是看两个函数错开多少。
由于 f 和 g两个函数都有一定的带宽分布(假若不用开头提到的”表述变化”就是都有一定的时间带宽分布),这个信号响应是在一定”范围”内广泛响应的。算总的响应信号,当然要把所有可能的响应加起来,实际上就是对所有可能t1积分了。积分范围虽然一般在负无穷到正无穷之间;但在没有信号或者没有响应的地方,积也是白积,结果是0,所以往往积分范围可以缩减。
这就是卷积及其物理意义啊。并成一句话来说,就是看一个时有时无(当然作为特例也可以永恒存在)的信号,跟一个响应函数在某一时刻有多大交叠。
*********拉普拉斯*********
拉普拉斯(1729-1827) 是法国数学家,天文学家,物理学家。他提出拉普拉斯变换(Laplace Transform)的目的是想要解决他当时研究的牛顿引力场和太阳系的问题中涉及的积分微分方程。
拉普拉斯变换其实是一个数学上的简便算法;想要了解其”物理”意义 — 如果有的话 — 请看我举这样一个例子:
问题:请计算十万乘以一千万。
对于没学过指数的人,就只会直接相乘;对于学过指数的人,知道不过是把乘数和被乘数表达成指数形式后,两个指数相加就行了;如果要问究竟是多少,把指数转回来就是。
“拉 普拉斯变换” 就相当于上述例子中把数转换成”指数” 的过程;进行了拉普拉斯变换之后,复杂的微分方程(对应于上例中”复杂”的乘法)就变成了简单的代数方程,就象上例中”复杂”的乘法变成了简单的加减法。再把简单的代数方程的解反变换回去(就象把指数重新转换会一般的数一样),就解决了原来那个复杂的微分方程。
所以要说拉普拉斯变换真有” 物理意义”的话,其物理意义就相当于人们把一般的有理数用指数形式表达一样。
另外说两句题外话:
1 。拉普拉斯变换之所以现在在电路中广泛应有,根本原因是电路中也广泛涉及了微分方程。
2。拉普拉斯变换与Z变换当然有紧密联系;其本质区别在于拉氏变换处理的是时间上连续的问题,Z变换处理的是时间上分立的问题。
Download from here
卷积的物理意义?
卷积(convolution, 另一个通用名称是德文的Faltung)的名称由来,是在于当初定义它时,定义成integ(f1(v)*f2(t-v))dv,积分区间在0到t之间。举个简单的例子,大家可以看到,为什么叫“卷积”了。比方说在(0,100)间积分,用简单的辛普生积分公式,积分区间分成100等分,那么看到的是f1(0)和f2(100)相乘,f1(1)和f2(99)相乘,f1(2)和f2(98)相乘,.........等等等等,就象是在坐标轴上回卷一样。所以人们就叫它“回卷积分”,或者“卷积”了。
卷积是个啥?我忽然很想从本质上理解它。于是我从抽屉里翻出自己珍藏了许多年,每每下决心阅读却永远都读不完的《应用傅立叶变换》。
3.1 一维卷积的定义
函数f(x)与函数h(x)的卷积,由函参量的无穷积分
定义虽然找到了,但我还是一头雾水。卷积是个无穷积分吗?那它是干啥用的?再往后翻:几何说明、运算举例、基本性质,一堆的公式,就是没有说它是干啥用的。我于是坐在那呆想,忽然第二个困扰我的问题冒了出来:傅立叶变换是个啥?接着就是第三个、第四个、……、第N个问题。
傅立叶变换是个啥?听说能将时域上的东东变到频域上分析?哎?是变到频域上还是空间域上来着?到底啥是时域,频域,空间域?
上网查傅立叶变换的物理意义,没发现明确答案,倒发现了许多和我一样晕着问问题的人。结果又多出了许多名词,能量?功率谱?图像灰度域?……没办法又去翻那本教材。
1.1 一维傅立叶变换的定义与傅立叶积分定理
设f(x)是实变量x的函数,该函数可实可复,称积分
为函数f(x)的傅立叶变换。
吐血,啥是无穷积分来着?积分是啥来着?还能记起三角函数和差化积、积化和差公式吗?我忽然有种想把高中课本寻来重温的冲动。
信号的时域的卷积等于频域的乘积。
利用这个性质以及特殊的δ函数可以通过抽样构造简单的调制电路
匹配滤波器最简单的形式就是原信号反转移位相乘积分得到的近似=相关
相关性越好得到的信号越强
还有解调中一些东西本质就是相关
卷积公式 解释 卷积公式是用来求随机变量和的密度函数(pdf)的计算公式。 定义式: z(t)=x(t)*y(t)=∫x(m)y(t-m)dm. 已知x,y的pdf,x(t),y(t).现在要求z=x+y的pdf.我们作变量替显,令 z=x+y,m=x. 雅可比行列式=1.那么,z,m联合密度就是f(z,m)=x(m)y(z-m)*1.这样,就可以很容易求Z的在(z,m)中边缘分布 即fZ(z)=∫x(m)y(z-m)dm.....由于这个公式和x(t),y(t)存在一一对应的关系。为了方便,所以记 ∫x(m)y(z-m)dm=x(t)*y(t) 长度为m的向量序列u和长度为n的向量序列v,卷积w的向量序列长度为(m+n-1), u(n)与v(n)的卷积w(n)定义为:w(n)=u(n)@v(n)=sum(v(m)*u(n-m)),m from 负无穷到正无穷; 当m=n时w(1) =u(1)*v(1) w(2) = u(1)*v(2)+u(2)*v(1) w(3) =u(1)*v(3)+u(2)*v(2)+u(3)*v(1) … w(n) = u(1)*v(n)+u(2)*v(n-1)+ …+u(n)*v(1) … w(2*n-1) = u(n)*v(n) 当m≠n时,应以0补齐阶次低的向量的高位后进行计算 这是数学中常用的一个公式,在概率论中,是个重点也是一个难点。
卷积公式是用来求随机变量和的密度函数(pdf)的计算公式。
定义式:
z(t)=x(t)*y(t)= ∫x(m)y(t-m)dm.
已知x,y的pdf,x(t),y(t).现在要求z=x+y的pdf. 我们作变量替显,令
z=x+y,m=x. 雅可比行列式=1.那么,t,m联合密度就是f(z,m)=x(m)y(z-m)*1.这样,就可以很容易求Z的在(z,m)中边缘分布
即fZ(z)=∫x(m)y(z-m)dm..... 由于这个公式和x(t),y(t)存在一一对应的关系。为了方便,所以记∫x(m)y(z-m)dm=x(t)*y(t)
高斯变换就是用高斯函数对图像进行卷积。高斯算子可以直接从离散高斯函数得到:
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
{
g[i*N+j]=exp(-((i-(N-1)/2)^2+(j-(N-1)/2)^2))/(2*delta^2));
sum += g[i*N+j];
}
}
再除以 sum 得到归一化算子
N是滤波器的大小,delta自选
信号与线性系统,讨论的就是信号经过一个线性系统以后发生的变化(就是输入输出和所经过的所谓系统,这三者之间的数学关系)。所谓线性系统的含义,就是,这个所谓的系统,带来的输出信号与输入信号的数学关系式之间是线性的运算关系。
因此,实际上,都是要根据我们需要待处理的信号形式,来设计所谓的系统传递函数,那么这个系统的传递函数和输入信号,在数学上的形式就是所谓的卷积关系。
卷积关系最重要的一种情况,就是在信号与线性系统或数字信号处理中的卷积定理。利用该定理,可以将时间域或空间域中的卷积运算等价为频率域的相乘运算,从而利用FFT等快速算法,实现有效的计算,节省运算代价。