视频压缩 运动估计算法

1、课程设计任务书 题目: 通信工程应用技术 初始条件： MATLAB 软件，电脑 要求完成的主要任务: 设计视频压缩系统中的运动估计算法：全搜索法（FS: Full Search）和三步法(TSS: Three Step Search)，比较二种方法的搜索点和每帧的峰值信噪比（PSNR： peak signal to noise ratio）要求：编制算法代码；对视频进行运动估计；计算PSNR时间安排：序号设 计 内 容所用时间1根据设计任务确定实验方案2天2根据实验条件进行电路的测试，并对结果进行分析7天3撰写课程设计报告1天合 计2周指导教师签名： 2013 年 月 日系主任（或责任教师）签

2、名： 2013 年 月 日 目录摘要IAbstractII1 设计任务12 实验原理及基本思想22.1 实验原理22.2基本思想23 运动估计算法33.1全搜索算法33.1.1全搜索算法及程序流程图33.1.2全搜索的特点43.2三步法43.2.1三步法的方法43.2.2三步法的特点41.搜索范围为-7,7；44 仿真结果64.1全搜索算法仿真结果64.2三步法算法仿真结果64.3全搜索算法和三步法指标对比124.3.1全搜索算法指标124.3.2三步法指标124.4仿真结果分析125 心得体会136 参考文献14试验程序15附录：28 摘要在视频编码和处理系统中，运动估计和运动补偿技术对降低

3、视频序列时间冗余度、提高编码效率起着非常关键的作用。运动估计的准确程度将直接决定视频编码器的编码效率。它极大地消除了视频序列的帧间相关性。运动估计算法的复杂性将直接决定视频压缩编码系统的复杂性，如何提高运动估计的效率，使运动估计算法的搜索过程更快速、更高效一直是人们研究的热点。掌握运动估计的块匹配算法，以及快速运动估计算法。关键词：运动估计、运动补偿技术、位移（运动）矢量 AbstractIn video coding and processing system, motion estimation and motion compensation technology to reduce th

4、e video sequence time redundancy, improve coding efficiency plays a key role. Motion estimation accuracy will directly decide video encoder coding efficiency. It greatly eliminates the video sequence correlation between the frame. The complexity of motion estimation algorithm will directly determine

5、 the complexity of video compression coding system, how to improve the efficiency of the motion estimation, make the motion estimation algorithm of the search process more quickly and more efficient has been a research hot spot. Master the block matching motion estimation algorithm, and fast motion

6、estimation algorithm.Keywords: motion estimation and、motion compensation technology and displacement（motion） vectorII1 设计任务随着视频压缩技术的发展,迫切需要解决视频数据的高效压缩问题、网络带宽变化适应性问题以及提高容错能力问题。视频数据压缩主要通过三种技术手段来实现:利用DCT变换和矢量量化来消除视频帧内的空间冗余;利用熵编码来消除编码冗余;利用运动估计来消除帧间的时间冗余。运动估计(Motion Estimation)是视频编码中广泛使用的一种技术，由于视频原始

7、图像间存在着大量的信息冗余，所以需要采用运动估计和运动补偿(Motion Compensation)技术来消除视频信号的冗余以提高编码效率。运动估计指压缩编码过程中计算当前帧和参考帧间的运动矢量。运动补偿是解码过程中根据运动矢量和当前帧恢复出下一帧。在帧间预测编码中，由于活动图像邻近帧中的景物存在着一定的相关性。因此，可将活动图像分成若干块或宏块，并设法搜索出每个块或宏块在邻近帧图像中的位置，并得出两者之间的空间位置的相对偏移量，得到的相对偏移量就是通常所指的运动矢量，得到运动矢量的过程被称为运动估计。运动矢量和经过运动匹配后得到的预测误差共同发送到解码端，在解码端按照运动矢量指明的位置，从已

8、经解码的邻近参考帧图像中找到相应的块或宏块，和预测误差相加后就得到了块或宏块在当前帧中的位置。通过运动估计可以去除帧间冗余度，使得视频传输的比特数大为减少，因此，运动估计是视频压缩处理系统中的一个重要组成部分。本节先从运动估计的一般方法入手，重点讨论了运动估计的三个关键问题：将运动场参数化、最优化匹配函数定义以及如何寻找到最优化匹配。然而运动估计的计算复杂度在整个视频数据压缩编码系统中最大。因此研究高效的运动估计算法对提高视频数据压缩编码的效率有着非常重要的意义。2 实验原理及基本思想 2.1 实验原理 在帧间预测编码中，由于活动图像邻近帧中的景物存在着一定的相关性。因此，可将活动图像分成若干

9、块或宏块，并设法搜索出每个块或宏块在邻近帧图像中的位置，并得出两者之间的空间位置的相对偏移量，得到的相对偏移量就是通常所指的运动矢量，得到运动矢量的过程被称为运动估计。 运动矢量和经过运动匹配后得到的预测误差共同发送到解码端，在解码端按照运动矢量指明的位置，从已经解码的邻近参考帧图像中找到相应的块或宏块，和预测误差相加后就得到了块或宏块在当前帧中的位置。 运动估计的准确程度往往用补偿图像与原图像比较的PSNR来衡量表示。2.2基本思想运动估计的基本思想是尽可能准确地获得序列图像帧间的运动位移，即运动矢量。运动估计算法的目标是效率和准确性。由于在成象的场景中一般有多个物体作不同的运动，如果直接按

10、照不同类型的运动将图像分割成复杂的区域是比较困难的。最直接和不受约束的方法是在每个像素都指定运动矢量，这就是所谓基于像素表示法。这种表示法是对任何类型图像都是适用的，但是它需要估计大量的未知量，并且它的解时常在物理上是不正确，除非在估计过程中施加适当的物理约束。这在具体实现时是不可能的，通常采用基于块的物体运动表示法。因为运动估计越准确，预测补偿的图像质量越高，补偿的残差就越小，补偿编码所需位数越少，需要传输的比特率就越小。3 运动估计算法 运动估计指压缩编码过程中计算当前帧和参考帧间的运动矢量。运动估计的基本思想是将图像序列的每一帧分成许多互不重叠的宏块，并认为宏块内所有象素的位移量都相同，

11、然后对每个宏块到参考帧某一给定特定搜索范围内根据一定的匹配准则找出与当前块最相似的块，即匹配块，匹配块与当前块的相对位移即为运动矢量。视频压缩的时候，只需保存运动矢量和残差数据就可以完全恢复出当前块。3.1全搜索算法3.1.1全搜索算法及程序流程图将图像的k+1帧分为不重叠，紧靠的N*N大小的块，然后依次对每一块进行处理。在处理某一块时，以该块的中心点为中心点，在源帧k中的窗口内的每个像素进行一次匹配，计算所有的MAD，然后找出最小的MAD的点，就是与之相匹配的点。开始获取图像分成N*N大小的块计算MSE进行块匹配计算MAD比较结束图3.1-1程序流程图3.1.2全搜索的特点全搜索具有以下一些

12、特点：1.对搜索区域的所有位置进行穷尽搜索；2.精度最高计算复杂，难以实时处理；3.必须研究相应的运动估计快速算法。3.2三步法三步法是一种较好的搜索算法，快速而且高效，它是在对数法的基础上对其进行了改进，提出在每一步搜索后搜索步长均减半的算法。它基本保持了FS的 性能，但其计算量只有FS的10左右。TSS在会议电视和可视电话中应用较多， 它通过三步搜索，逐步较小搜索步长。若最大搜索范围为士7，搜索精度取一个像 素，则步长为4、2、1，共需三步即可满足要求，从而得名。3.2.1三步法的方法三步法是一种经典的快速方法，如图2所示，按八方向（19）搜索，从中选取候选矢量得指向3的矢量，0为前一帧搜

13、索值（或初值）。紧接搜索八方向（AH），以此降低运算量，提高运算速度。根据这样的思路，针对不同情况可以做进一步的改进。 图3.2-1十字法取样3.2.2三步法的特点1.搜索范围为-7,7；2.搜索模板半径依次减半；3.3 峰值信噪比 峰值信噪比（经常缩写为PSNR）是一个表示信号最大可能功率和影响它的表示精度的破坏性噪声功率的比值的工程术语。由于许多信号都有非常宽的动态范围，峰值信噪比常用对数分贝单位来表示。PSNR是“Peak Signal to Noise Ratio”的缩写。peak的中文意思是顶点。而ratio的意思是比率或比列的。整个意思就是到达噪音比率的顶点信号，psnr是一般是用

14、于最大值信号和背景噪音之间的一个工程项目。通常在经过影像压缩之后，输出的影像通常都会有某种程度与原始影像不一样。为了衡量经过处理后的影像品质，我们通常会参考PSNR 值来认定某个处理程序够不够令人满意。PSNR计算公式如下：Peak 就是指8 bits 表示法的最大值255。MSE 指 Mean Square Error（均方误差，各值相差的n次方和的平均值的n次平方根(这几个字应该没有)），I(角标n)指原始影像第n个pixel 值，P(角标n)指经处理后的影像第n个pixel 值。PSNR 的单位为dB。所以PSNR值越大，就代表失真越少。PSNR 是最普遍，最广泛使用的评鉴画质的客观量测

15、法，不过许多实验结果都显示，PSNR 的分数无法和人眼看到的视觉品质完全一致，有可能 PSNR 较高者看起来反而比 PSNR 较低者差。这是因为人眼的视觉对于误差的敏感度并不是绝对的，其感知结果会受到许多因素的影响而产生变化（例如：人眼对空间频率较低的对比差异敏感度较高，人眼对亮度对比差异的敏感度较色度高，人眼对一个区域的感知结果会受到其周围邻近区域的影响）。4 仿真结果4.1全搜索算法仿真结果 图4.1-1 the first frame 图4.1-2 the second frame 图4.1-3 帧间差值图4.1-4 DFD图4.1-5 恢复后的第二帧图像 图4.1-6 第一帧图像的运动

16、矢量图图4.1-7 第二帧图像的运动矢量图4.2三步法算法仿真结果 图4.2-1 the first frame 图4.2-2 the second frame 图4.2-3 帧间差值图4.2-4 DFD图4.2-5 恢复后的第二帧图像图4.2-6 第一帧图像的运动矢量图图4.2-7 第二帧图像的运动矢量图4.3全搜索算法和三步法指标对比4.3.1全搜索算法指标全搜索算法：耗时长为Elapsed time is 1.053320 seconds.4.3.2三步法指标三步法：耗时长为Elapsed time is 0.881349 seconds.4.4仿真结果分析 根据峰值信噪比可知全搜索算法

17、的效果比三步法的效果要好；但全搜索所消耗的时间要比三步法所消耗的时间长，其运行较慢，因为全搜索算法的搜索点数要比三步法的搜索点数多很多。 5 心得体会通过本次对全搜索块匹配算法课程设计，让我更深刻的理解了多媒体信息的获取、搜索和保存的基本方法。基于块的运动估计，是先将视频图像分成一个个规则的图像块，然后对每个图像块估计运动矢量。基于块的运动估计和运动补偿已经广泛应用于各种视频压缩编码标准。因此，本次课程设计是非常有作用的，让我们清醒的认识到自己还有好大的差距，它让我了解了多媒体处理的基本原理。课程设计中会遇到很多很多的困难，但在同学的帮助和自己思考下最终还是没有解决问题。所以，本次课程设计我个

18、人的收获不大，但是让我知道自己还在图像编程方面有很多欠缺的地方，不仅对均方误差(MSE)、绝对值误差(MAD)更深的理解，并且采用全搜索块匹配方法得到所有宏块的运动矢量，把c语言的知识又重新的过了一遍，同时还理解了些多媒体信息处理的基本原理。这次课程设计为以后毕业设计敲响了一次警钟。6 参考文献1一种适用于H_263建议的块匹配运动估计算法2 H_263中运动估计算法的TMS320C80实现3快速块匹配运动估计算法的探索与分析 浙江大学 第13期srtp研究报告4几种经典快速块匹配运动估计算法的比较研究 肖敏连5 黎洪松. 数字视频处理M. 北京邮电大学出版社，2006试验程序FS：clear

19、 all;I1=imread('claire1.bmp'); %read the first frameI2=imread('claire2.bmp'); %read the second framedm=7;%¸给图像扩边，每个边都扩大dm大小% |-|-|-|% | B | C | D |% |-|-|-|% | | | |% | | | |% | E | A | F | -> 图像% | | | |% | | | |% | | | |% |-|-|-|% | G | H | I |% |-|-|-|% I1=double(I1);I2=do

20、uble(I2);rownum colnum = size(I1); II=zeros(rownum+2*dm,colnum+2*dm);II(dm+1:dm+rownum,dm+1:dm+colnum)=I1;for i=1:dm II(i,dm+1:dm+colnum)=II(dm+1,dm+1:dm+colnum); II(rownum+dm+i,dm+1:dm+colnum)=II(dm+rownum,dm+1:dm+colnum);endfor j=1:dm II(1:rownum+2*dm,j)=II(1:rownum+2*dm,dm+1); II(1:rownum+2*dm,co

21、lnum+dm+j)=II(1:rownum+2*dm,dm+colnum);end %下面进行全搜索算法 blocksize=16; rowblocks =rownum/blocksize; colblocks =colnum/blocksize; A=99999999999999999999; %为了找到最小的均方误差，A用于设定一个很大的初值 Eij=0; xrecord=ones(16,16); %xrecord，yrecord用于存放匹配快的块号，即运动矢量 yrecord=ones(16,16); diff=zeros(256,256); %这幅图的大小为256*256 tic f

22、or x=0:(rowblocks-1) %x表示行中第几个子块 row=x*blocksize; for y=0:(colblocks-1) %y表示列中第几个子块 col=y*blocksize;% tempx=x*blocksize+1:(x+1)*blocksize;% tempy=y*blocksize+1:(y+1)*blocksize; for p=-dm:dm for q=-dm:dm %（p,q）表示x,y对应子块在前一帧所的搜索位置 Eij=0; Eij=sum(sum(I2(row+1:row+blocksize,col+1:col+blocksize)-II(row+d

23、m+p+1:row+dm+p+blocksize,col+dm+q+1:col+dm+q+blocksize).2)/(blocksize2); if Eij<A A=Eij; xrecord(x+1,y+1)=p; yrecord(x+1,y+1)=q; end end end A=999999999999999999; for mx=1:blocksize for ny=1:blocksize diff(row+mx,col+ny)=I2(row+mx,col+ny)-II(row+mx+dm+xrecord(x+1,y+1),col+ny+dm+yrecord(x+1,y+1);

24、end end end endtoc figure,imshow(I1,); title('the first frame'); figure,imshow(I2,); title('the second frame'); IIII=I2-I1; figure,imshow(IIII,); title(帧间差值'); figure,imshow(diff,); title('DFD'); %title('利用全搜索算法匹配后的帧间差'); for x=0:(rowblocks-1) % row=x*blocksize; f

25、or y=0:(colblocks-1) % col=y*blocksize; III(row+1:row+blocksize,col+1:col+blocksize)=II(row+dm+xrecord(x+1,y+1)+1:row+dm+xrecord(x+1,y+1)+blocksize,col+dm+yrecord(x+1,y+1)+1:col+dm+yrecord(x+1,y+1)+blocksize)+diff(row+1:row+blocksize,col+1:col+blocksize); end end %III=I1+abs(diff); figure,imshow(III

26、,); title('»Ö¸´ºóµÄµÚ¶þ֡ͼÏñ'); ERR=diff;% figure,imshow(ERR,);% title('DFD'); numberarray=0:1:255; for m=1:255 numberarray(m+1)=0; end; zeronumber=0; for n=1:rownum for m=1:colnum dif=abs(ERR(

27、m,n); if(dif=0) temp=zeronumber; zeronumber=temp+1; else numberarray(dif)=numberarray(dif)+1; end; end; end; figure;plot(0,zeronumber,'k*');hold on;plot(numberarray,'r*'),title('DFD distribution');hold off; ERR1=zeros(16,16); for i=0:15 for j=0:15 ERR1(i+1,j+1)=round(sum(sum(

28、ERR(i*blocksize+1:i*blocksize+blocksize,j*blocksize+1:j*blocksize+blocksize)/(blocksize*blocksize); end end numberarray=0:1:255; for m=1:255 numberarray(m+1)=0; end; zeronumber=0; for n=1:16 for m=1:16 dif=abs(ERR1(m,n); if(dif=0) temp=zeronumber; zeronumber=temp+1; else numberarray(dif)=numberarray

29、(dif)+1; end; end; end; figure;plot(0,zeronumber,'k*');hold on;plot(numberarray,'r*'),title('DFD(block average) distribution');hold off; %figure;mesh(diff);figure;contour(diff,15); figure; for i=1:16 for j=1:16 quiver(i,j,xrecord(i,j)/16,yrecord(i,j)/16); hold on; end end gri

30、d on;figure;quiver(1:16,1:16,yrecord,xrecord); grid on; S_3SS：clear all; I1=imread('claire1.bmp'); %read the first frameI2=imread('claire2.bmp'); %read the second frame dm=7;I1=double(I1);I2=double(I2);rownum colnum = size(I1);II=zeros(rownum+2*dm,colnum+2*dm);II(dm+1:dm+rownum,dm+1:

31、dm+colnum)=I1;for i=1:dm II(i,dm+1:dm+colnum)=II(dm+1,dm+1:dm+colnum); II(rownum+dm+i,dm+1:dm+colnum)=II(dm+rownum,dm+1:dm+colnum);endfor j=1:dm II(1:rownum+2*dm,j)=II(1:rownum+2*dm,dm+1); II(1:rownum+2*dm,colnum+dm+j)=II(1:rownum+2*dm,dm+colnum);end ¨ blocksize=16; rowblocks =rownum/blocksize;

32、 colblocks =colnum/blocksize; A=99999999999999999999; % Eij=0; xrecord=ones(16,16); % yrecord=ones(16,16); diff=zeros(256,256); % tic for x=0:(rowblocks-1) % row=x*blocksize; for y=0:(colblocks-1) % col=y*blocksize;% tempx=x*blocksize+1:(x+1)*blocksize;% tempy=y*blocksize+1:(y+1)*blocksize; for p1=-

33、4:4:4 %µÚÒ»²½ for q1=-4:4:4 % Eij=0; Eij=sum(sum(I2(row+1:row+blocksize,col+1:col+blocksize)-II(row+dm+p1+1:row+dm+p1+blocksize,col+dm+q1+1:col+dm+q1+blocksize).2)/(blocksize2); if Eij<A A=Eij; xrecord(x+1,y+1)=p1; yrecord(x+1,y+1)=q1; end end end p1=xrecord(x+1,y+1)

34、; q1=yrecord(x+1,y+1); for p2=p1-2:2:p1+2 %µÚ¶þ²½ for q2=q1-2:2:q1+2 if p2=p1 | q2=q1 Eij=0; Eij=sum(sum(I2(row+1:row+blocksize,col+1:col+blocksize)-II(row+dm+p2+1:row+dm+p2+blocksize,col+dm+q2+1:col+dm+q2+blocksize).2)/(blocksize2); if Eij<A A=Eij; xrecord(x+1,y+1)=

35、p2; yrecord(x+1,y+1)=q2; end end end end p2=xrecord(x+1,y+1); q2=yrecord(x+1,y+1); for p3=p2-1:1:p2+1 %µÚÈý²½ for q3=q2-1:1:q2+1 if p3=p2 | q3=q2 Eij=0; Eij=sum(sum(I2(row+1:row+blocksize,col+1:col+blocksize)-II(row+dm+p3+1:row+dm+p3+blocksize,col+dm+q3+1:col+dm+q3+bloc

36、ksize).2)/(blocksize2); if Eij<A A=Eij; xrecord(x+1,y+1)=p3; yrecord(x+1,y+1)=q3; end end end end A=999999999999999999; for mx=1:blocksize for ny=1:blocksize diff(row+mx,col+ny)=I2(row+mx,col+ny)-II(row+mx+dm+xrecord(x+1,y+1),col+ny+dm+yrecord(x+1,y+1); end end end endtoc figure,imshow(I1,); titl

37、e('the first frame'); figure,imshow(I2,); title('the second frame'); IIII=I2-I1; figure,imshow(IIII,); title('Ö¡¼ä²îÖµ'); figure,imshow(diff,); title('DFD'); for x=0:(rowblocks-1) row=x*blocksize; for y=0:(colblocks-1) col=y*blocksize; III(row+1:row+blocksize,col+1:col+blocksize)=II(row+dm+xrecord(x+1,y+1)+1:row+dm+xrecord(x+1,y+1)+blocksize,col+dm+yrecord(x+1,y+1)+1:col+dm+yrecord(