DSP 项目1 - 基于DCT的音频压缩代码模板

▼ 📄 ex0_all.m（主函数，无需修改）

clc;
clear;
close all;

% ex0文件作为总文件，可在ex全部完成后运行测试，可对参数部分进行修改体验区别
global path path_out fs;  % 定义全局变量方便调用 (增加fs)

% 定义音频文件路径
path = './audio/greatest_work.flac';    % 设定path为原始音频路径，格式flac
path_out = './audio/greatest_work_compressed.flac';   % 设定path_out为重建音频路径，格式flac

% 窗函数大小设置
window = 512;  % 设定窗函数大小为512

% 调用 ex1_read_audio_cut 函数读取音频并进行分段
[x_cut, N, nb] = ex1_read_audio_cut(window, path);        

% 调用 ex2_dct_x_cut 对分段音频作DCT变换
[X_CUT, row, col] = ex2_dct_x_cut(x_cut);                       

% 设置阈值百分比
perc = 0.99;

% 调用 ex3_energy_cut 对DCT结果进行阈值处理
[X_CUT_99, ind, i] = ex3_energy_cut(X_CUT, perc);             

% 调用 ex4_quantization 对阈值化结果作量化和解量化
[q, tmp, quant] = ex4_quantization(X_CUT_99);                 

% 调用 ex5_idct_quant 对解量化结果作IDCT变换
I_X_quant = ex5_idct_quant(quant);                  

% 调用 ex6_recombine_I_X 对IDCT结果作合并重建
[dis, reconstruct] = ex6_recombine_I_X(I_X_quant);         

% 调用 ex7_analyze 对比原始信号和重建信号，默认显示图形
ratio = ex7_analyze(reconstruct, true);             

disp('音频文件已保存');

▼ 📄 ex1_read_audio_cut.m（读取音频并分段）

function [x_cut, N, nb] = ex1_read_audio_cut(window, path)
%EX1_READ_AUDIO_CUT 读取音频文件并对其进行分段

global fs; % 声明全局变量fs
[x, fs_local] = audioread(path);  % 使用 audioread 函数读取音频文件，得到音频数据 x 和采样频率 fs
fs = fs_local; % 将读取到的采样率存入全局变量

% N 表示音频数据的总长度，即采样点的数量。
N = ?;  % 使用 length 函数

% nb 表示可以分成多少段音频，即音频数据长度 N 除以 window 的值。
nb = ?;  % 使用 floor 函数

% x_cut 是一个二维矩阵，每列表示一段音频，每段的长度为 window 个采样点。
x_cut = ?;  % 使用 zeros 函数初始化

% 循环遍历每一段音频数据，并将其填充到 x_cut 的每一列中。
for n = 1:nb
    % 在此处将每一段音频数据填充到 x_cut 的第 n 列。
    x_cut(:, n) = ?;  % 使用 n*window 及其变形定义起点和终点
end
end

▼ 📄 ex2_dct_x_cut.m（DCT变换）

function [X_CUT, row, col] = ex2_dct_x_cut(x_cut)
%EX2_DCT_X_CUT 将分割后的音频信号分别作DCT计算

[row, col] = ?;  % 获取x_cut的行数和列数
X_CUT = ?;       % zeros初始化dct计算结果，应和时域信号大小一致

for c = 1:col                   % 循环，对x_cut的每一列作dct计算
    X_CUT(:,c) = ?;  % 使用dct函数
end
end

▼ 📄 ex3_energy_cut.m（能量阈值处理）

function [X_CUT_99, ind, i] = ex3_energy_cut(X_CUT, perc)
%EX3_ENERGY_CUT 将DCT计算结果进行阈值化处理

[row, col] = size(X_CUT);         % 取得DCT系数矩阵的行数和列数
X_CUT_99 = X_CUT;                 % 初始化阈值化后的结果矩阵
ind = cell(col, 1);               % 初始化索引单元数组
i = zeros(col, 1);                % 初始化保留系数数量向量

for c = 1:col                     % 对每一列进行操作
    % 对每一列系数的绝对值进行降序排列，并获取对应的索引
    [sorted_coeffs, sorted_idx] = ?;
    
    % 计算每列的总能量
    total_energy = sum(sorted_coeffs.^2);
    cumulative_energy = ?;     % 初始化累积能量
    count = 1;                    % 初始化阈值位置

    % 当累积能量不足总能量的 perc 部分时，增加 count 的值
    while ?
        cumulative_energy = cumulative_energy + sorted_coeffs(count)^2;
        count = ?;
    end
    
    % 保留前 count-1 个最重要的系数
    i(c) = count - 1;
    ind{c} = sorted_idx(1:i(c));
    
    % 获取小于这些能量的位置并置零
    zero_indices = sorted_idx(?); % 从第 count 个到最后的索引
    X_CUT_99(?) = 0;  % 将这些位置的值置为零
end
end

▼ 📄 ex4_quantization.m（量化与解量化）

function [q, tmp, quant] = ex4_quantization(X_CUT_99)
%EX4_QUANTIZATION 对DCT系数进行量化和解量化

% 设置量化位数
b = 6;  % 6位量化，可表示2^6=64个级别

% 找出最大值作为量化区间
L = ?;  % 使用max和abs函数

% 计算量化步长
q = ?;  % 量化步长公式

% 量化：将连续值映射到离散级别
tmp = ?;  % 使用round函数四舍五入

% 解量化：将离散值映射回连续值
quant = ?;  % 乘以量化步长

end

▼ 📄 ex5_idct_quant.m（IDCT逆变换）

function I_X_quant = ex5_idct_quant(quant)
%EX5_IDCT_QUANT 对量化后的DCT系数进行逆变换

[row, col] = size(quant);
I_X_quant = zeros(row, col);

for c = 1:col
    I_X_quant(:, c) = ?;  % 使用idct函数对quant的第c列进行逆变换
end
end

▼ 📄 ex6_recombine_I_X.m（合并重建）

function [dis, reconstruct] = ex6_recombine_I_X(I_X_quant)
%EX6_RECOMBINE_I_X 将IDCT后的片段合并重建完整音频

global path_out fs; % 声明全局变量fs

[row, col] = size(I_X_quant);

% 计算重建信号的总长度
dis = ?;  % 行数乘以列数

% 初始化重建信号
reconstruct = ?;  % 创建一维零向量

% 将每个片段拼接起来
for n = 1:col
    start_idx = ?;  % 起始索引
    end_idx = ?;    % 结束索引
    reconstruct(start_idx:end_idx) = ?;  % 将第n列赋值到对应位置
end

% 保存重建的音频
audiowrite(path_out, reconstruct, fs);  % 使用从原始音频读取的正确采样率fs
end

▼ 📄 ex7_analyze.m（性能分析）

function ratio = ex7_analyze(reconstruct, show_plot)
%EX7_ANALYZE 分析压缩效果并可视化对比

global path;

% 读取原始音频
[x_original, fs] = audioread(path);

% 确保长度一致（取较短的长度）
min_length = ?;  % 使用min函数
x_original = x_original(1:min_length);
reconstruct = reconstruct(1:min_length);

% 计算均方误差MSE
mse = ?;  % 使用mean函数

% 计算压缩比（简化计算，实际应考虑存储的系数数量）
% 假设原始为16位，压缩后为6位量化
ratio = 16 / 6;

% 显示结果
fprintf('压缩比: %.2f:1\n', ratio);
fprintf('均方误差MSE: %.6f\n', mse);

% 如果需要显示图形
if show_plot
    figure('Name', '原始信号 vs 重建信号', 'NumberTitle', 'off');
    
    % 绘制时域对比
    subplot(2,1,1);
    plot_range = 1:min(5000, min_length);  % 只显示前5000个点
    plot(plot_range, x_original(plot_range), 'b', 'LineWidth', 1);
    hold on;
    plot(plot_range, reconstruct(plot_range), 'r--', 'LineWidth', 1);
    legend('原始信号', '重建信号');
    title('时域波形对比');
    xlabel('采样点');
    ylabel('幅度');
    grid on;
    
    % 绘制误差
    subplot(2,1,2);
    error_signal = x_original - reconstruct;
    plot(plot_range, error_signal(plot_range), 'g', 'LineWidth', 1);
    title('重建误差');
    xlabel('采样点');
    ylabel('误差幅度');
    grid on;
end
end

📝 DSP 项目1 - 基于DCT的音频压缩实验代码模板