前言
本文集中前面主要介绍了离散数据的傅里叶变换,并且得到了较好的效果!那既然有了傅里叶变换这个工具,为什么还需要小波变换呢?因为:傅里叶变换只能告诉你原始信号中有哪些频率,但不能告诉你这些频率的信号出现在什么时间!也就说明:如果信号是"时变"的(频率随着时间是改变的),那么单纯用傅里叶变换所能反映的信息就十分有限了!因此,针对时变信号,我们使用小波变换。图1展示"时变信号"与"时不变信号"区别:

图1:时不变信号与时变信号
时不变与时变的区别,看下面的实现的代码就很轻易理解:

x = 0:0.001:1;

% 4个频率:
f1 = 50; f2 = 80; f3 = 110; f4 = 20;

% 时不变信号: 多频
y1 = sin(2pif1x) + sin(2pif2x) + sin(2pif3x) + sin(2pif4x);
% 时变信号: 多频
y2 = sin(2pif1x).(x<=0.3) + sin(2pif2x).(x>0.3 & x<=0.6)+…
sin(2pif3x).(x>0.6 & x<=0.8) + sin(2pif4x).(x>0.8);

figure(1);

subplot(2,1,1);
plot(y1);
grid on;
xlabel(‘采样点’); ylabel(‘振幅’);
axis([0 1000 -inf inf]);
title(‘时不变信号’);

subplot(2,1,2);
plot(y2);
grid on;
xlabel(‘采样点’); ylabel(‘振幅’);
axis([0 1000 -inf inf]);
title(‘时变信号’)

本文同样考虑的是离散数据的小波变换使用。通过手动matlab编程实现小波变换"塔式分解"与"重构"来深刻了解小波变换实现的内在含义。之后,借助matlab自带的一系列相关小波变换程序来实现"时频分析"和"小波去噪"。

说明:本文更加侧重详细介绍matlab自带各种小波功能函数的使用!除了小波分解与重构的程序我们手动实现外,其他的各种操作都建议用自带函数实现。

小波分解:
小波分解的流程总结为:先将信号对半分解成"低频近似"与"高频细节"2个部分;同样的操作每次将上一次的"低频近似"部分再分成低频近似和高频细节部分,逐次细分(最多分解到每个部分只有1个点)。每次分出的高频细节部分不做分解。因此:每次分出低频近似部分相当于对本次信号做"低通滤波",分出的高频细节部分相当于对本次信号做"高通滤波"。所以:每次小波分解就是用1个低通滤波器和1个高通滤波器对本次信号做1次低通滤波和1次高通滤波而已。

由上述说明可得:小波分解的关键在于2个(一组)滤波器。对于现实的离散数据而言,滤波器看上去很高大上其实就是很简单的数字而已,滤波听起来很难,其实就是做"点乘相加"而已。离散小波分解中最简单的一组滤波器为:

低通滤波器:[0.5, 0.5];高通滤波器:[0.5, -0.5]。

下面举例说明如何用上面这一组最简单滤波器对离散数据进行小波分解:

假设我们的离散数据为:[2,5,8, 9, 7, 4, -1, 1]

(1) 第一级分解:

低通滤波:
20.5 + 50.5 = 1 + 2.5 = 3.5
80.5 + 90.5 = 4 + 4.5 = 8.5
70.5 + 40.5 = 3.5 + 2 = 5.5
-10.5 + 10.5 = -0.5 + 0.5 = 0

得到低频近似部分:[3.5,8.5,5.5,0]

高通滤波:
20.5 + 5(-0.5) = 1 - 2.5 = -1.5
80.5 + 9(-0.5) = 4 - 4.5 = -0.5
70.5 + 4(-0.5) = 3.5 - 2 = 1.5
-10.5 + 1(-0.5) = -0.5 - 0.5 = -1

得到高频细节部分:[-1.5,-0.5,1.5,-1]

(2)第二级分解:

分解上一次分解得到的低频近似部分,高频细节不动。即待分解信号为:
[3.5,8.5,5.5,0]
低通滤波:
3.50.5 + 8.50.5 = 1.75 + 4.25 = 6
5.50.5 + 00.5 = 2.25 + 0 = 2.25
得到的低频近似部分为:[6,2.25]

高通滤波:
3.50.5 + 8.5(-0.5) = 1.75 - 4.25 = -2.5
5.50.5 + 0(-0.5) = 2.25 + 0 = 2.25
得到的高频细节部分为:[-2.5,2.25]

(3)第3级分解:
同理,待分解信号是上一次分解得到的低频近似部分:[6,2.25]
低通滤波:
60.5 + 2.250.5 = 3 + 1.125 = 4.125
得到的低频近似部分为:[4.125]

高通滤波:
60.5 + 2.25(-0.5) = 3 - 1.125 = 1.875
得到的高频细节部分为:[1.875]

到此为止,例子中给出的原始离散信号就达到了其最大分解级数(分解到元素只有1个)。整个过程很简单,只需注意每次点乘求和元素是无重合的。整个的多级分解过程如图2所示:
图2:离散信号小波多级分解示意图

注意:不同组的高通和低通滤波中都有这样的一个规律:两者的区别只是高通滤波器中第2个值是负数而已;数都是一样的。

小波多级分解清楚了,那怎么"重构/恢复"回去呢?塔式分解的逆向合成而已。根据滤波器的规律,我们可以设:

以2级分解到3级为例,我们知道:
6a + 2.25a = 3级低通近似
6a + 2.25(-a) = 3级高通细节

那么逆过程就是,我们知道了3级低通近似和高通细节的值,我们还知道滤波器的数值(a已知),然后反推2级低通近似和高通细节数值,即:
2级低通近似a + 2级高通细节a = 4.125
2级低通近似a + 2级高通细节(-a) = 1.875

所以:

这样我们就得到了2级低通近似和2解高通细节,然后逐级往上递推即可实现重构/恢复~ so easy.

整个分解过程我们清楚了,现在我们引入一些专业的名词:在离散数据中,一组低通高通滤波器,其实就是"小波基函数"!取不同的小波基函数其实就是滤波器里面的数值不同而已。最常用的"haar小波基"。下面我们就利用haar小波基,在matlab里手动实现小波分解与重构:

matlab手动实现小波分解程序:

clc ; clear;

% 每次修改这里的原始数据, 个数最好是2^n
% x = [9 7 3 5];
x = [2 5 8 9 7 4 -1 1];
% x = [2 5 8 9 7 4 -1 1 2 1 8 3 8 0 3 1];

order_max = log(length(x))/log(2);
fprintf(‘当前数据最多分解%d阶\n’,order_max);
order = double(input(‘自定义分解阶数( order<order_max ):’));

% matlab默认的haar小波基函数:
low = [1/sqrt(2) 1/sqrt(2)];
high = [1/sqrt(2) -1/sqrt(2)];

new = zeros(1,length(x)); % 最后结果 —— 均值 + 差值
ave = zeros(1,length(x)/2); % 均值(低频)记录
dec = zeros(1,length(x)/2); % 差值(高频)记录

% 小波循环分解部分: 其实就是低通和高通2种情况的卷积计算
m = 1;
xtmp = x;
for norder = 1:order
for n = 1:2:length(xtmp)
ave(m) = sum(xtmp(n:n+1).*low);
dec(m) = sum(xtmp(n:n+1).*high);
m = m + 1;
end
% 下面2句的赋值过程, 总体是从后往前赋值的~
% 进入过的高频就不会动了; 进入的低频再下次就会被自己分解出的高频和低频取代——总体从后往前√
new( length(xtmp)/2+1:length(xtmp) ) = dec( 1:2^(order_max-norder) ); % 记录每次更新的高频内容;
new( 1:length(xtmp)/2 ) = ave( 1:2^(order_max-norder) ); % 记录每次更新的低频内容;
xtmp = ave( 1:2^(order_max-norder) );
m = 1;
end

fprintf(‘手写%d级分解结果为:\n’,order);
new

fprintf(‘matlab自带%d级分解结果为:\n’,order);
wavedec(x,order,‘haar’)
matlab手动实现小波重构/恢复程序:

clc ; clear;

% 每次修改这里的原始数据, 个数最好是2^n
% x = [9 7 3 5];
x = [2 5 8 9 7 4 -1 1];
% x = [2 5 8 9 7 4 -1 1 2 1 8 3 8 0 3 1];

order_max = log(length(x))/log(2);
fprintf(‘当前数据最多分解%d阶\n’,order_max);
order = double(input(‘自定义分解阶数( order<order_max ):’));

% 直接用知自带的分解函数:
new = wavedec(x,order,‘haar’);

% 小波重构——任意haar基函数全能重构回去
% 计算的系数: 与基函数有关
tmp1 = 1/( low(1) + high(1) );

% 迭代重构开始:
xrec = zeros(1,length(new)); % 记录复原的数据
new_rec = new;

for norder = order_max+1-order : order_max % 这个规律让任意低阶都可以适用
m = 1; % 专门用来记录"前半段"位置——给奇数位用的
for n = 1:2^norder
half = 2^norder/2; % 当前重构区间的一半——也是用来给奇数位计算用的
if mod(n,2) ~= 0
% 奇数时操作:
xrec(n) = tmp1*( new_rec(m) + new_rec(m+half) );
m = m + 1;
else
% 偶数时操作:
xrec(n) = (new_rec(m-1) - low(1)*xrec(n-1))/low(2);
end
end
new_rec(1:n) = xrec(1:n); % 每次要更新的(后一次重构基于前一次结果): 从前往后
end

fprintf(’%d级分解后重构:\n’,order);
new_rec

fprintf(‘自带的%d级重构结构:\n’,order);
[C,S] = wavedec(x,order,‘haar’);
waverec(C,S,‘haar’)

链接:https://www.jianshu.com/p/8847e6eebe16

手动实现一维离散数据小波分解与重构相关推荐

  1. 哈儿小波分解和重构(降维和升维)实现算法

    [0]README 0.1)本文旨在讲解 哈儿小波变换(分解和重构)进行数据的降维和升维: [timestamp: 1703281610]时隔几个月再来review 哈儿小波变换算法的具体思路: 1) ...

  2. python小波分解与重构_小波分解和重构

    小波变换能够很好地表征一大类以低频信息为主要成分的信号, 小波包变换可以对高频部分提供更精细的分解 详见(http://www.cnblogs.com/welen/articles/5667217.h ...

  3. python小波分解与重构_python - 使用pyWavelets进行多级局部小波重构 - 堆栈内存溢出...

    我设法编写了自己的wrcoef函数版本,该版本似乎可以正常工作: import pywt import numpy as np def wrcoef(X, coef_type, coeffs, wav ...

  4. 多元经验模态分解_交通运输|基于小波分解和长短时记忆网络的地铁进站量短时预测...

    山东科学 ›› 2019, Vol. 32 ›› Issue (4): 56-63.doi: 10.3976/j.issn.1002-4026.2019.04.008 摘要: 针对城市地铁车站进站客流 ...

  5. 【小波变换】离散小波分解Discrete Wavelet Transform

    此篇博客记录自学离散小波分解的相关内容,以后若有更多理解在此篇更新. 一. 为什么需要离散小波分解    除离散变换外,还有连续小波分解,通过改变分析窗口大小,在时域上移动窗口和基信号相乘,最后在全时 ...

  6. 一维信号小波去噪算法C语言,[转载]一维小波分解与去噪重构

    对随机一维信号实现多尺度小波分解,长度M=256,层数N=3. 解:一维随机信号是用nelec函数作为信号源,产生0~256范围内的信号. 实现小波按层分解的函数是: [C,L]=wavedec(s, ...

  7. 小波分解处理非平稳时序数据(降噪)

    最近在搞毕业论文的预测部分,我收集到的销售数据波动性较大,尤其是存在尖峰数据,我导师给我提供了两种降噪的方法,一种是小波分解,一种是经验模式分解.这两种方法按照老师的话来说就是方法很旧,但是很有效.关 ...

  8. matlab小波变换数据少,一维离散数据小波变换实用案例

    前言 小波变换专业处理时变信号!其重要用途包含:突变点检测.时频分析.信号降噪等.本文将详细介绍小波变换的这3种主要用途,借助具体例子来说明并总结相关函数的使用. 间断点检测 现实信号中的间断点是较为 ...

  9. 基于小波分解与LSTM的城市轨道短时客流预测

    1.文章信息 文章题为<A novel prediction model for the inbound passenger flow of urban rail transit>,是一篇 ...

  10. C语言实现小波分解,提取近似与细节分量,包含详细例程

    C语言实现小波分解,提取近似与细节分量,包含详细例程 声明 本文的C语言实现小波分解非本人原创,均参考了网络上的文章(详见最后的参考资料),程序主要来自李承宇的文章和程序. 我只对程序进行了少量的修改 ...

最新文章

  1. 网易严选 x 网易有数:数据产品+数据中台双引擎模式实践
  2. epoll 版 高并发服务器
  3. save()、saveOrUpdate()、merge()的区别
  4. Nand Flash与Nor Flash
  5. 编程填空:学生信息处理程序_项目学生:业务层
  6. alexeyab darknet 编译_【目标检测实战】Darknet—yolov3模型训练(VOC数据集)
  7. 只用2000行代码实现google protocol buffer c++版的功能
  8. 第九篇:Spring Boot整合Spring Data JPA_入门试炼05
  9. CSS3-边框-外轮廓-文本-渐变-WEB字体
  10. 单片机modbus rtu通讯_Modbus-RTU通讯
  11. 在不推动提交的情况下触发Travis-CI重建?
  12. dreamweaver 正则表达式为属性值加上双引号_IT兄弟连 HTML5教程 HTML5表单 新增的表单属性3...
  13. java中rhino什么用_使用require.js和Java / Rhino解析模块
  14. 只要听说过电脑的人都能看懂的网上pdf全书获取项目
  15. Windows 7程序开发系列之一(任务栏篇)
  16. 選大學了﹖請看“上大學”網
  17. 基于SSM的校园二手交易平台
  18. java求航班飞行时间代码,基于JAVA的航班动态接口调用代码实例
  19. cdn回源php_简述回源原理和CDN缓存
  20. qq邮箱对方服务器退回,为什么我用QQ邮箱发邮件被退回来了?他说地 – 手机爱问...

热门文章

  1. Linux修改默认静态IP
  2. csdn头像修改失败的解决办法
  3. vant 个人中心头像修改
  4. iOS之iCloud云存档实现笔记
  5. 破解vba工程密码——VBA代码
  6. 安卓一键新机_知道华为手机变慢的罪魁祸首吗?用这四招两年旧机秒变新机
  7. JavaScript的js文件压缩和格式化工具
  8. Android抓包思想总结
  9. TRNSYS 内区之间通风原理试验
  10. 手游游戏资源提取 (破解、AssetStudio、VGMToolbox、disunity、Il2CppDumper、 .NET Reflector)...