Vondrak滤波原理详解及Matlab实现

  • 一.Vondrak基本思想:
  • 二、Vondrak平滑法的原理:
  • 三、Vondrak滤波平滑公式:
  • 四、Vondrak滤波应用
  • 五、Matlab实现

一.Vondrak基本思想:

通过选择不同的而平滑因子控制数据平滑的程度,在观测数据的绝对拟合与绝对平滑之间选择一条折衷的曲线,保留有用信号,滤除噪声信号。
这种方法 对于等间距的观测数据、不等间距的观测数据均适用。
Vondrak滤波的本质是通过确定合理的平滑因子

在观测数据的绝对拟合和决定平滑之间寻求一条折衷的曲线。
它在不需知道观测资料的变化规律及其拟合函数的情况下,就能够对观测资料进行有效的平滑。
Vondrak滤波的关键是 选择合理的平滑因子,目前常用的方法有:频率响应法、观测误差法、平滑误差法以及交叉证认法。
前三种的应用都必须在一定的已知条件下采用,如观测序列的频率及先验标准差,而交叉证认法的主观随意性和计算量很大。

二、Vondrak平滑法的原理:

修匀数学中指出一个修匀序列满足光滑性要求,所谓光滑性要求,就是在满足

很小的情况下,所得出的修匀值位于一条光滑的曲线上。
S为平滑度,yi‘ 为修匀值,上式中z为差分阶数,n为待处理数据的个数。
但在实际处理中,一味的追求光滑性将会使得估计值偏离原始的观测值,使得处理的结果与真值的偏离程度增大。
在实际处理中,会预先采用一定的处理方法,去除粗大误差,使得观测值与真值的偏离程度很小。
考虑到整体加权及正负抵消的情况下,令

F为拟合度,yi为观测值,yi‘ 为修匀值,上式中z为差分阶数,n为观测值的个数。
同样的,当当拟合度时,说明修勾值与原观测值完全重合,此时并没有起到修匀效果。
Whittaker于1923年将2-1 、2-2 做了线性组合,发展为Whittaker修匀:

其中:pi为观测值的权序列,yi 为观测值序列,yi’ 为平滑值序列,h 是正常数。在满足上述准则下所获得的修匀值既满足位于一条光滑的曲线上,又与原观测序列具有一定的拟合度。
1967年捷克天文学家J.Vondrak将Whittaker修匀中原定一的平滑度改为以平滑之的3阶差分的平方和,发展为Whittaker-Vondrak平滑法,简称Vondrak滤波,
Vondrak平滑基本准则是:

1976年对此进行了改进,将拟合度和平滑度分别用他们的平均值来代替:

其中,pi为观测值的权序列,yi 为观测值序列,yi’ 为平滑序列,λ是正常数,

在观测数据的绝对拟合和绝对平滑之间起着平衡的作用,ε 越小,曲线的平滑程度越强。

三、Vondrak滤波平滑公式:


在满足平滑准则达到最小的情况下,利用拉格朗日多项式推到的Vondrak滤波的基本方程组为:

实际应用中,常用2-5 改进的Vondrak滤波。解方程组2-15也可获得满足准则2-5的平滑结果。

特点
Vondrak滤波可以在未知观测资料的变化规律且未知其拟合函数的情况下,就能够对观测资料进行有效的平滑,因此在许多领域都得到了广泛的应用。

四、Vondrak滤波应用

对于某一确定平滑因子周期信号被分成了三段:通过带(F=1)、过滤带(0《 F《 1)、压制带(F《 0)。Vondrak滤波器是很好的低通滤波器,具体应用案列详见:

吴芸芸. Vondrak滤波准则及应用研究[D]. 中南大学, 2012.

五、Matlab实现

function yy = vondrak(x,y,epsilon,sigma)
% Vondrak滤波程序 n = length(x);
B = zeros(n,1); p = ones(1,n);
p = p*sigma; a = zeros(1,n+3);
b = zeros(1,n+3);
c = zeros(1,n+3);
d = zeros(1,n+3);
for i=1:n-3 a(i+3) = 6*sqrt(x(i+2)-x(i+1))/((x(i)-x(i+1))*(x(i)-x(i+2))*(x(i)-x(i+3))); b(i+3) = 6*sqrt(x(i+2)-x(i+1))/((x(i+1)-x(i))*(x(i+1)-x(i+2))*(x(i+1)-x(i+3))); c(i+3) = 6*sqrt(x(i+2)-x(i+1))/((x(i+2)-x(i))*(x(i+2)-x(i+1))*(x(i+2)-x(i+3))); d(i+3) = 6*sqrt(x(i+2)-x(i+1))/((x(i+3)-x(i))*(x(i+3)-x(i+1))*(x(i+3)-x(i+2)));
end A = zeros(n,7);
for i = 1:n A(i,1) = a(i)*d(i); A(i,2) = a(i+1)*c(i+1)+b(i)*d(i); A(i,3) = a(i+2)*b(i+2)+b(i+1)*c(i+1)+c(i)*d(i); A(i,4) = epsilon*p(i)/(n-3)+a(i+3)^2+b(i+2)^2+c(i+1)^2+d(i)^2; A(i,5) = a(i+3)*b(i+3)+b(i+2)*c(i+2)+c(i+1)*d(i+1); A(i,6) = a(i+3)*c(i+3)+b(i+2)*d(i+2); A(i,7) = a(i+3)*d(i+3); B(i) = epsilon*p(i)/(n-3);
end
A(3,1:6) = A(3,2:7); A(3,7) = 0;
A(2,1:5) = A(2,3:7); A(2,6:7) = [0 0];
A(1,1:4) = A(1,4:7); A(1,5:7) = [0 0 0];
y = y.*B;
ls = 4;
for k = 1:n-1 max = 0; for i = k:ls t = abs(A(i,1)); if t>max max = t; is = i;   %列选主元 end end help1 = y(k); y(k) = y(is); y(is) = help1;  %常数向量交换行 help2 = A(k,:); A(k,:) = A(is,:); A(is,:)=help2;  %系数矩阵交换行 y(k) = y(k)/A(k,1);   A(k,:) = A(k,:)/A(k,1);  %系数归一化 for i = k+1:ls t = A(i,1); y(i) = y(i)-y(k)*t;  %常数向量消元 A(i,:) = A(i,:)-A(k,:)*t;  %系数矩阵消元 A(i,1:6) = A(i,2:7);  %系数矩阵左移一位 A(i,7) = 0; end if ls~=n ls = ls+1; end
end
q = A(n,1);
y(n) = y(n)/q;
ls = 2;
for i = n-1:-1:1 y(i) = y(i)-A(i,2:ls)*y(i+1:i+ls-1); if ls~=7 ls = ls+1; end
end
yy = y;

参考文献:
吴芸芸. Vondrak滤波准则及应用研究[D]. 中南大学, 2012.

Vondrak滤波原理详解及Matlab实现相关推荐

  1. 【RS码2】RS码的BM迭代译码原理详解及MATLAB实现(不使用MATLAB库函数-代码见CSDN同名资源)

    关注公号[逆向通信猿]更精彩!!! 理论基础 订阅<信道编码>专栏,首先查阅各子程序的详解 [有限域生成]本原多项式生成有限域的原理及MATLAB实现 [有限域除法]二元多项式除法电路原理 ...

  2. 【BCH码2】BCH码的快速BM迭代译码原理详解及MATLAB实现(不使用MATLAB库函数-代码见CSDN同名资源)

    关注公号[逆向通信猿]更精彩!!! 理论基础 订阅<信道编码>专栏,首先查阅各子程序的详解 [有限域生成]本原多项式生成有限域的原理及MATLAB实现 [有限域除法]二元多项式除法电路原理 ...

  3. 【BCH码3】BCH码的彼得森译码原理详解及MATLAB实现(不使用MATLAB库函数『需要完整代码请先私信』)

    结果预览 理论基础 订阅<信道编码>专栏,首先查阅各子程序的详解 [有限域生成]本原多项式生成有限域的原理及MATLAB实现 [有限域除法]二元多项式除法电路原理及MATLAB详解 [有限 ...

  4. opencv 高通滤波和低通滤波_一阶低通滤波原理详解

    在汽车标定中,使用最多的滤波算法即低通滤波,很多朋友可能知道怎么标定,但是不清楚具体的原理,因此本文将介绍一阶低通滤波的原理.算法建模仿真和优缺点: 一阶滤波算法的原理 一阶滤波,又叫一阶惯性滤波,或 ...

  5. 用Burg法估计AR模型的参数原理详解及matlab实现

    用Burg法估计AR模型的参数. 借助如图所示的格型预测误差滤波器,伯格法通过求出前向预测误差和后向预测误差的平均功率来选取最佳的反射系数k,使误差的平均功率取得最小值,进而通过反馈求出模型系数和噪声 ...

  6. 一阶低通滤波器方程_一阶低通滤波原理详解

    在汽车标定中,使用最多的滤波算法即低通滤波,很多朋友可能知道怎么标定,但是不清楚具体的原理,因此本文将介绍一阶低通滤波的原理.算法建模仿真和优缺点:一阶滤波算法的原理 一阶滤波,又叫一阶惯性滤波,或一 ...

  7. 短波信道模型--多径瑞利信道原理详解及matlab实现

    瑞利衰落是一种小尺度衰落效应,它总是叠加于大尺度衰落效应上如衰减.阴影等.发射端和接收端相对运动速度的大小决定了信道衰落的快慢.相对运对导致接收信号存在多普勒频移,即信道衰落的快慢与多普勒频移的大小有 ...

  8. 蚁群算法原理详解和matlab代码

    1原理: 蚂蚁在寻找食物源的时候,能在其走过的路径上释放一种叫信息素的激素,使一定范围内的其他蚂蚁能够察觉到.当一些路径上通过的蚂蚁越来越多时,信息素也就越来越多,蚂蚁们选择这条路径的概率也就越高,结 ...

  9. 快速引导滤波原理详解(fast guided filter)

    符号标记 q : 输出图像 p : 输入图像 I : 引导图 a k , b k : 线性变换模型 r : 滤波半径 s : 下采样因子 q: 输出图像 \\[2ex] p: 输入图像 \\[2ex] ...

最新文章

  1. postgresql 可调试
  2. 机器人瓦力 配乐_《WALL-E》机器人小王子
  3. C语言学习之有一个已排好序的数组,要求输入一个数后,按原来排序的规律将它插入数组中
  4. 哈夫曼编译器c语言程序,哪位大牛有哈夫曼编码的C语言源程序,麻烦帮帮忙啦!...
  5. spark ui 上schedulingDelay理解
  6. Linux内核分析 - 网络[十]:ARP杂谈
  7. 网易整合邮箱和博客 可通过邮件更新博客日志
  8. JDK11使用HSDB
  9. 《OpenGL编程指南》第8版 第9版 VS2015 VS2017配置方法
  10. 方舟原始恐惧mod生物代码_方舟:生存进化新手攻略
  11. java堆排序思想及代码实现
  12. 网络编程day1-本地信息的获取
  13. oracle求和函数的写法,Oracle over函数的学习笔记三 求和函数的使用
  14. 纬度和经度的最大长度是多少?
  15. 润乾报表 echarts统计图分类显示不全
  16. macbook更新系统服务器,为Mac新系统做准备 苹果更新OSX Server
  17. 益聚星荣:不打老婆的即时到账”?还呗贷款平台广告词惹争议
  18. CVPR2019目标检测方法进展
  19. 洗牌一副n张牌,每一张牌都用字母顺序标记。
  20. Android : TextView

热门文章

  1. 【大数据】大数据的特点
  2. 【Python】AttributeError: ‘DatetimeProperties’ object has no attribute ‘weekday_name’ 的解决方法
  3. 来,一起来实现一个符合Promise/A+的Promose(1.0.1版本)
  4. Microsoft Enterprise Library 5.0 系列教程(二) Cryptography Application Block (初级)
  5. 兔子--百度地图所需的jar+so下载地址
  6. 【日常填坑】之ajax请求laravel的api接口
  7. hibernate 使用别名查询
  8. 聊聊有哪些参与项目的好途径吧
  9. HONGJIN4 2013
  10. 《深入体验 飞鸽传书 开发内幕 核心基础》