sinc插值和NUFFT小结
篇一:sinc基函数分解
FFT要求时域数据是均匀采样,频率是均匀分布于(-π,π)的。
一、Sinc插值原理:
设有函数f(x),采样定理表明,在满足以下两个条件时,就可以从f(x)的等间隔离
式(1)可以理解为所有输入样本的加权叠加,根据此式就可以计算任意点处f(x)函数值。为精确计算某一点上的f(x)需要覆盖无限多个点。实际上这是无法做到的,而且使用大量的数据点会使得插值非常耗时,但精度提高却很小。计算表明核值随着与x的间隔增大而降低,这意味着可以在不过度损失精度的同时对卷积核进行截断。在现有的计算资源下个8点的加权sinc函数比较适合数据处理。
当使用截断后的sinc函数对存在陡峭边缘的函数进行插值时,会出现一种称为Gibbs效应的振铃现象。为减少这种影响,应对插值进行加窗锐化,通常加入Kaiser窗的处理效果比较理想。对与加窗后的插值核,需要进行归一化处理,使其增益单位化,否则采样点上的权值和不再等于1,并且不同插值点之间会出现较大误差。
1、功率谱估计的提出:主要是针对经典谱估计的分辨率低和方差性能不好的问题提出的。
二、非均匀离散傅里叶变换(NUDFT)
1、{tn}不均匀采样而{fn}均匀采样,即fmmf,则有F(wm)称为第一类非均匀DFT(nonuniform DFT,NUDFT)。
2、{fn}不均匀采样而{tn}均匀采样,即tnnt,则有F(wm)称为第二类NUDFT。
3、{fn}和{tn}都是均匀采样,即tnnt并且fmmf,则有
f(tn)e
n0
d
N1
j2mftn
,
f
n0
N1
j2nmtm
,(t)edn
F(wm)fd(tn)e
n0N1n0
N1
j2nmtf
,称为均匀DFT。特别是,如果f1/(Nt),则有
F(wm)fd(tn)ej2nm/N,这就是传统意义上的DFT。
为了方便,我们将两类NUFFT表示为完全离散的符号形式:
M1
Hn
m0
he
mN/2
n
j2vmn
Gm
nN/2
ge
j2vmn
其中n=-N/2,-N/2+1,…N/2-1,vm[1/2,1/2]表示非均匀采样位置。
三、NUFFT算法
1、第一类NUDFT可由NUFFT算法加速的流程
1)、由公式NUFFT算法采用最小均方误差的解x(c)(AHA)1AHb(c)计算插值系数xr(cm),这里r=-q/2,-q/2+1,…q/2,m=0,1,…M-1。需要O((q1)M)个乘数。 2)、计算Fourier系数l
mrm
r[q/2,q/2][ucm]rml
j2nl/uN
2
h
x(c)。需要O((q1)M)个乘法数。
。需要O(uNlog(uN))个乘法数。
luN/2
3)、使用FFT计算Tn
l
luN/2
e
4)、 Hsn.Tn。需要O(N)个乘法。sinc基函数分解。
2、第二类NUDFT可由NUFFT算法加速的流程
1)、由公式由公式NUFFT算法采用最小军方误差的解x(c)(AA)Ab(c)计算插值系数xr(cm),这里r=-q/2,-q/2+1,…q/2,m=0,1,…M-1。需要O((q1)M)个乘数。
1
2)、计算tnsn.gn。需要O(N)个乘法数。
luN/2
1
H1H
2
3)、使用FFT计算Tn
l
luN/2
e
j2nl/uN
。需要O(uNlog(uN))个乘法数。
4)、计算G
r
r[q/2,q/2]
x(c
m
)T[ucm]r,需要O((q1)M)个乘法数。
可以看出,NUFFT-1和NUFFT-2算法十分相似,并且具有相同的计算复杂度。他们的区别在于:NUFFT-1处理非均匀间隔的输入数据,所以为了使用FFT技术,需要首先将输入数据插值到均匀间隔上,如NUFFT-1流程第2步所示,NUFFT-2处理均匀间隔输入非均匀间隔输出变换,所以可以首先使用FFT得到在变换空间的均匀采样数据,然后将他们插值到期望的输出采样间隔上,该插值过程在NUFFT-2中的最后一步完成。
四、NUFFT 的优点:NUFFT
是一种时频非均匀采样的变换的快速算法,把数据直
接进行NUFFT变换,而不像传统方法需要先对数据进行插值才能用来成像,因此可以减少差值带来的误差以及所需花费的额外操作。
关于sinc函数第一副瓣的问题
篇二:sinc基函数分解
关于sinc函数第一副瓣的问题
sinc函数定义为,,函数图形如下:
图一
一些书中对sinc函数第一副班处的值有一些分歧,本文就此问题做出一点解答。 下面讨论(如图二),第一副瓣的取值问题。
第一副瓣最大值问题就是sinc函数极值问题,下面用导数方法来求第一副瓣极大值。
令
得,,
排除点,解此方程可得第一
副瓣处t= 1.2901;计算|sinc(1.2901)|= 0.1950.归一化后,用分贝表示第一副瓣相对于t=0处的值(1),。
五种常见小波基函数及其matlab实现
篇三:sinc基函数分解
与标准的傅里叶变换相比,小波分析中使用到的小波函数具有不唯一性,即小波函数 具有多样性。小波分析在工程应用中,一个十分重要的问题就是最优小波基的选择问题,因为用不同的小波基分析同一个问题会产生不同的结果。目前我们主要是通过用小波分析方法处理信号的结果与理论结果的误差来判定小波基的好坏,由此决定小波基。常用小波基有Haar小波、Daubechies(dbN)小波、Mexican Hat(mexh)小波、Morlet小波、Meyer小波等。
Haar小波
Haar函数是小波分析中最早用到的一个具有紧支撑的正交小波函数,也是最
简单的一个小波函数,它是支撑域在t的定义如下:
[0,1]范围内的单个矩形波。Haar
0t1
1t1 其他
函数
1(t)-10
Haar小波在时域上是不连续的,所以作为基本小波性能不是特别好。但它也
有自己的优点:
1. 计算简单。
2.
而且与自己的整数位移正交,因此,(t)不但与(2jt)[jz]正交,在a
2j的多分辨率系统中,Haar小波构成一组最简单的正交归一的
小波族。
(t)的傅里叶变换是:
()=j
4
sin2()ej/2
a
Daubechies(dbN)小波
Daubechies小波是世界著名的小波分析学者Inrid·Daubechies构造的小波函数,简写为dbN,N是小波的阶数。小波(t)和尺度函数
(t)中的支撑区为
2N1,(t)的消失矩为N。除N1(Harr小波)外,dbN不具有对称性
(即非线性相位)。除N
1(Harr小波)外,dbN没有明确的表达式,但转换
函数h的平方模是明确的:
令
p(y)C
k0
N1
N-1k
k
y
k
,其中Ck
N-1k
为二项式的系数,则有
m0()
其中:
2
(cos
2
2
)p(sinsinc基函数分解。
2
2
)
1
hem()2
2N1k0
k
-jk
Daubechies小波具有以下特点:
1. 在时域是有限支撑的,即(t)长度有限。 2. 在频域3.
()在=0处有N阶零点。
k
(t-k)dt(t)和它的整数位移正交归一,即(t)
为低通函数,长度有限,支撑域在t=0~2N-1的范围内。
。
4. 小波函数(t)可以由所谓“尺度函数”
(t)求出来。尺度函数(t)
db4的时域和频域波形:
Daubechies小波常用来分解和重构信号,作为滤波器使用:sinc基函数分解。
Mexican Hat(mexh)小波
Mexican Hat函数为Gauss函数的二阶导数:
(t)(1
t2)e
2
t22
()2e2
因为它的形状像墨西哥帽的截面,所以也称为墨西哥帽函数。 Mexihat小波的时域和频域波形:
Matlab图形作业
篇四:sinc基函数分解
MATLAB具有强大的图形处理功能。下面给出了3个m脚本文件,请在MATLAB环境下运行,观察其输出。
要求根据每个m文件输出的图形(共18个),用中文翻译并解释产生每个输出图形的函数具体是什么?其功能是什么?
文件1:
clear all
clf %清除图像
mfilename('fullpath') %返回调用的m文件的完整路径和名称,但不包括文件扩展名 echo on %显示程序
subplot(2,3,1) %使(2*3)幅子图中第1个子图成为当前图
t = 0:0.1:10; %定义变量t,取值为1~10,间距为0.1
z = impulse(1, [1 1 1], t); %分析以t为时间长度的单位脉冲响应,然后赋值给z
stairs(t(1:5:end),z(1:5:end)) %把1~5赋值给向量t,向量z的奇数位,绘制t,z阶梯图 hold on %保持当前坐标,继续绘图
plot(t,z,'r') %用红色线绘制t,z所成图形
plot([0 t(end)], [0 0], 'k:') %用黑色虚线绘制z=0的曲线
title('Impulse Response - (STAIRS)') %添加该图标题'Impulse Response - (STAIRS)'
subplot(2,3,2) %使(2*3)幅子图中第2个子图成为当前图
theta = 2*pi*(0:74)/75; %把0~2*pi作75等分,并把值赋给theta
x = cos(theta); %定义变量x,x=cos(theta)
y = sin(theta); %定义变量y,y=sin(theta)
z = abs(fft(ones(10,1), 75))';
%生成一个全为1的10*1的矩阵,将矩阵的末尾补零使其长度为75,
返回75点的不连续Fourier变换,然后将变换后的值取绝对值,赋
值给z
stem3(x, y, z) %在参数x,y的指定位置上绘出z的元素
title('Polar FFT - (STEM3)') %在当前坐标轴上方正中央放置字符串'Polar FFT - (STEM3)'作为标题
subplot(2,3,3) %使(2*3)幅子图中第3个子图成为当前图
[X,Y,Z] = peaks(-2:0.25:2); %分别创建X,Y,Z三个均为2*2阶的方阵,且间距为0.25
[U,V] = gradient(Z, 0.25); %求Z在X,Y方向上的差分,每个方向上点之间的间距为0.25 contour(X,Y,Z,10); %绘制Z的等高线,X,Y分别为为X、Y轴的限度。线条密度为10 hold on %保持当前坐标上继续绘图
quiver(X,Y,U,V); %在每个由X和Y中相应元素对指定的坐标位置处绘制向量,向量以箭头来表示
title('Surface Gradient - (CONTOUR & QUIVER)')
%添加该图标题为'Surface Gradient - (CONTOUR & QUIVER)'
theta = 0:0.1:4*pi; %将0~4*pi赋值给theta,间距为0.1
[x,y] = pol2cart(theta(1:5:end), theta(1:5:end));
%把存放在数组theta,theta相应元素中的极坐标数据转换为二维
笛卡尔坐标或者xy坐标
subplot(2,3,4) %使(2*3)幅子图中第4个子图成为当前图
polar(theta,theta) %创建一个幅角theta,半径theta的极坐标图
axis([-13 13 -12.5 14.5]) %设置坐标范围取为-13~13,-12.5~14.5
title('Spiral Plot - (POLAR)') %添加该图标题为'Spiral Plot - (POLAR)’
subplot(2,3,5) %使(2*3)幅子图中第5个子图成为当前图
compass(x,y) %显示指南针图形,图形具有n个箭头,n是x或者y的元素个数,每
一个箭头的起始点的位置都是原点
axis([-13 13 -12.5 14.5]) %设置坐标范围取为-13~13,-12.5~14.5
title('Direction Vectors - (COMPASS)') %添加该图标题为'Direction Vectors - (COMPASS)'
subplot(2,3,6) %使(2*3)幅子图中第6个子图成为当前图
feather(x(1:19),y(1:19)) %显示沿水平轴等间距点出发的羽毛状的向量图像,x和用都从1~19,
间距为1
axis([1 21 -5 10]) %坐标轴的范围为1~21,-5~10
title('Direction Vectors - (FEATHER)') %添加该图标题为'Direction Vectors - (FEATHER)'
set(gcf,'Position', [64 111 887 564]) %对所有由gcf确定的对象,将单元数组Position中指
定的属性设置为单元数组[64 111 887 564]里对应的值
echo off %不显示程序
shg %显示图形窗
文件2:
clear all
clf %清除图形
echo on %显示程序
data = [10 2 3 5; 5 8 10 3; 9 7 6 1; 3 5 7 2; 4 7 5 3];%定义一个5*4的矩阵data
subplot(2,3,1)%使(2*3)幅子图中第1个子图成为当前图
bar(data, 'stacked');%为’stacked’中的每一个元素在data设定的位置处绘制一个条带 title('Bar Graph - (BAR, ''stacked'')');%将该图的标题设为'Bar Graph - (BAR, ''stacked'’)’ subplot(2,3,2)%使(2*3)幅子图中第2个子图成为当前图
bar3h(data);%绘制三维条带图,data中的每个元素对应于一个条带
title('Horizontal Bar Graph - (BAR3H, ''grouped'')');
%将该图的标题设为'Horizontal Bar Graph - (BAR3H,
''grouped'')’
subplot(2,3,3)%使(2*3)幅子图中第3个子图成为当前图
hist(randn(1000,3));%产生一个1000*3的随机数组,将向量data中元素10等分,用柱状图
显示数据值得分布
title('Histogram - (HIST)');%将该图的标题设为'Histogram - (HIST)'
subplot(2,3,4)%使(2*3)幅子图中第4个子图成为当前图
area(data);%将data中的元素显示为多条曲线,并且将每一条曲线的下方填充
title('Area Plot - (AREA)');%将该图的标题设为'Area Plot - (AREA)'
subplot(2,3,5)%使(2*3)幅子图中第5个子图成为当前图
pie3(sum(data), [0 0 1 0]);%将data中数据按列求和,将求和数据绘制一个饼图,从饼图的中
心分离图中的第三个切片,其余三个切片位置不变。
title(['3-D Pie Chart';' (PIE3) ']);
%将该图的标题设为['3-D Pie Chart';' (PIE3) ']
subplot(2,3,6)%使(2*3)幅子图中第6个子图成为当前图
rose(5/3*randn(1000,1), 18); %生成一个1000*1的随机数组,该数组*5/3,生成一个18瓣的玫瑰图
title('Polar Histogram - (ROSE)');%将该图的标题设为'Polar Histogram - (ROSE)‘
set(gcf,'Position',[184 248 740 424])
%对所有由gcf确定的对象,将单元数组Position中指定的
属性设置为单元数组[184 248 740 424]里对应的值
echo off %不显示程序
shg %显示图形窗
文件3:
echo on %显示程序
subplot(2,3,1) %使(2*3)幅子图中第1个子图成为当前图
x = -3:0.3:3; y = x; %定义x,范围为-3~3,间距为0.3,定义y=x;
[X,Y]=meshgrid(x,y); %由向量x和y所指定的域变换为矩阵X和Y,得到的矩阵可用来计算两个
变量的函数和绘制三维图形
[theat,R] = cart2pol(X,Y);%将存储在数组X、Y相应元素中的三维笛卡尔坐标变换为圆柱坐标 Z = sinc(R);%定义Z=sinc(R)
contourf(peaks(30), 10) %定义一个30*30的矩阵,绘制该矩阵的10级等高线图
colorbar%显示表征颜色刻度的块
grid on %添加主要的网格线
title('Peaks Function - (CONTOURF & COLORBAR)') %将该图的标题设为'Peaks Function -
(CONTOURF & COLORBAR’
subplot(2,3,2)%使(2*3)幅子图中第2个子图成为当前图
plot3(X,Y,Z)%绘制X、Y、Z的三维曲线图
grid on%在当前坐标轴添加主要的网格线
axis([-3 3 -3 3 -1 1])%设置X、Y、Z坐标的范围分别为-3~3,-3~3,-1~1
title('Sinc Function - (PLOT3)')%将该图的标题设为'Sinc Function - (PLOT3)’
subplot(2,3,3)%使(2*3)幅子图中第3个子图成为当前图
waterfall(membrane(1));%绘制一个L型的薄片,然后使用x=1:size(membrane(1),1)和y=1:
size(membrane(1),1)来创建瀑布表面图。membrane(1)决定
着表面的颜色。
title('L-shaped Membrane - (WATERFALL)')%将该图的标题设为'L-shaped Membrane -
(WATERFALL)'
subplot(2,3,4) %使(2*3)幅子图中第4个子图成为当前图
contour3(peaks(30), 25);%创建一个30*30的矩阵,然后在三维视图下绘制该矩阵的25级等
高线图
title('Peaks Function - (CONTOUR3)')%将该图的标题设为'Peaks Function - (CONTOUR3)'
subplot(2,3,5) %使(2*3)幅子图中第5个子图成为当前图
mesh(X,Y,Z)%绘制三维网线图,Z定义绘制颜色
axis([-3 3 -3 3 -1 1])%定义X、Y、Z的坐标范围分别为-3~3,-3~3,-1~1
title('Sinc Function - (MESH)')%将该图的标题设为'Sinc Function - (MESH)'
subplot(2,3,6) %使(2*3)幅子图中第6个子图成为当前图
surf(membrane(1))%绘制一个L型的薄片,然后由生成的矩阵绘制三维阴影表面图 title('L-shaped Membrane - (SURF)')%将该图的标题设为'L-shaped Membrane - (SURF)'
set(gcf,'Position',[211 248 713 413])%对所有由gcf确定的对象,将单元数组Position中指定的
属性设置为单元数组[211 248 713 413]里对应的值
echo off%不显示程序
shg%显示图形窗
姓名:熊智秋
学号201205050229
http://m.myl5520.com/mingrenmingyan/125496.html