MATLAB从零开始实现离散小波变化DWT

一、基础目标

在MATLAB中从零开始实现离散小波变换（DWT）是一项能让你深入理解小波分析核心思想的实践。下面将梳理其基本原理、具体的实现步骤、关键代码以及一些重要的注意事项。

二、离散小波变换基础

离散小波变换（DWT）的核心思想是通过一组低通和高通滤波器，将信号分解成不同频率的子带。这种多分辨率分析方法让我们能够同时观察信号的时域和频域特性。

核心概念与Mallat算法

DWT的数学本质是通过滤波器组对信号进行多尺度分解。Mallat算法是实现这一过程的高效方法，它包含两个核心步骤：

1. 分解（正向DWT）：信号通过低通滤波器（获得近似系数cA）和高通滤波器（获得细节系数cD），然后进行下采样（通常为2倍）。
2. 重构（逆向DWT）：对近似系数和细节系数进行上采样，然后通过重构滤波器组，最终恢复原始信号。

下表对比了DWT分解中的两种主要系数：

系数类型	表示信息	滤波器类型	物理意义
近似系数 (cA)	信号的低频成分	低通滤波器	信号的整体趋势、轮廓
细节系数 (cD)​	信号的高频成分	高通滤波器	信号的细节、突变点、噪声

小波滤波器系数获取

在实现DWT前，我们需要先获得小波滤波器的系数。MATLAB的wfilters函数可以帮我们获取标准小波的滤波器系数。

function [Lo_D, Hi_D, Lo_R, Hi_R] = get_wavelet_filters(wavelet_name)
% 获取小波滤波器系数
% 输入：wavelet_name - 小波名称，如'db4'
% 输出：Lo_D - 分解低通滤波器, Hi_D - 分解高通滤波器
% Lo_R - 重构低通滤波器, Hi_R - 重构高通滤波器

% 检查输入的小波名称是否有效
valid_wavelets = {
   
   'haar', 'db1', 'db2', 'db3', 'db4', 'db5', ...
                 'sym2', 'sym3', 'sym4', 'coif1', 'coif2'};
if ~any(strcmpi(wavelet_name, valid_wavelets))
    error('不支持的小波名称。请使用如db4、sym2等标准小波。');
end

% 获取滤波器系数
[Lo_D, Hi_D, Lo_R, Hi_R] = wfilters(wavelet_name);

fprintf('小波%s的滤波器系数：\n', wavelet_name);
fprintf('低通分解滤波器 (Lo_D): [%s]\n', num2str(Lo_D, '%.4f '));
fprintf('高通分解滤波器 (Hi_D): [%s]\n', num2str(Hi_D, '%.4f '));
fprintf('低通重构滤波器 (Lo_R): [%s]\n', num2str(Lo_R, '%.4f '));
fprintf('高通重构滤波器 (Hi_R): [%s]\n', num2str(Hi_R, '%.4f '));
end

三、从零实现DWT分解与重构

单层DWT分解实现

下面是单层DWT分解的具体实现，重点关注卷积和下采样过程：

function [cA, cD] = my_dwt(x, Lo_D, Hi_D, mode)
% 单层一维离散小波变换分解
% 输入：x - 输入信号, Lo_D - 分解低通滤波器, Hi_D - 分解高通滤波器
% mode - 边界延拓模式 ('sym'对称延拓, 'zpd'零延拓, 'ppd'周期延拓)
% 输出：cA - 近似系数, cD - 细节系数

if nargin < 4
    mode = 'sym'; % 默认对称延拓
end

N = length(x);
L = length(Lo_D);

% 边界延拓处理
switch mode
    case 'sym' % 对称延拓 (MATLAB默认方式)
        x_ext = [fliplr(x(1:L-1)), x, fliplr(x(end-L+2:end))];
    case 'zpd' % 零延拓
        x_ext = [zeros(1, L-1), x, zeros(1, L-1)];
    case 'ppd' % 周期延拓
        x_ext = [x(end-L+2:end), x, x(1:L-1)];
    otherwise
        error('不支持的边界延拓模式');
end

% 卷积运算（不使用MATLAB的conv函数，直接实现）
cA_conv = zeros(1, N+L-1);
cD_conv = zeros(1, N+L-1);

for i = 1:N+L-1
    for j = 1:L
        if (i-j+1) >= 1 && (i