降噪处理

%% ========== 读取对齐后的数据并进行降噪处理 ==========
clear; clc; close all;% 加载对齐后的数据
fprintf('正在加载对齐后的数据...\n');
load('对齐后的数据_85点.mat', 'aligned_data', 'time_axis');% 设置参数
datasets_to_process = 4:12;  % 要处理的数据集
target_length = 85;  % 时间点数% 设置变量名称（19个变量，已去掉第6和19个）
var_names = {'变量1', '变量2', '变量3', '变量4', '变量5', '变量7', ...'变量8', '变量9', '变量10', '变量11', '变量12', '变量13', '变量14', ...'变量15', '变量16', '变量17', '变量18', '变量20', '变量21'};% 降噪参数设置
% Savitzky-Golay滤波器参数
sg_order = 3;        % 多项式阶数（推荐2-4）
sg_framelen = 11;    % 窗口长度（必须为奇数，推荐9-15）% 也可以选择其他降噪方法，取消注释以下代码：
% method = 'savitzky-golay';  % Savitzky-Golay滤波
% method = 'moving-average';  % 移动平均
% method = 'wavelet';         % 小波降噪
method = 'savitzky-golay';  % 默认使用SG滤波% 存储降噪后的数据
denoised_data = cell(12, 1);%% ========== 对每个数据集进行降噪处理 ==========
fprintf('\n开始降噪处理...\n');
fprintf('使用方法: %s\n', method);for i = datasets_to_processdata = aligned_data{i};[n, m] = size(data);denoised = zeros(n, m);% 对每个变量分别降噪for col = 1:msignal = data(:, col);switch methodcase 'savitzky-golay'% Savitzky-Golay滤波器denoised(:, col) = sgolayfilt(signal, sg_order, sg_framelen);case 'moving-average'% 移动平均滤波window_size = 7;denoised(:, col) = movmean(signal, window_size);case 'wavelet'% 小波降噪[denoised(:, col), ~] = wdenoise(signal, 'Wavelet', 'sym4', ...'DenoisingMethod', 'Bayes');endenddenoised_data{i} = denoised;fprintf('数据集 %d 降噪完成\n', i);
endfprintf('\n所有数据集降噪处理完成！\n');%% ========== 可视化对比（降噪前后） ==========
fprintf('\n开始生成对比图像...\n');% 设置颜色
colors_original = lines(9);
colors_denoised = colors_original * 0.7;  % 降噪后的颜色稍暗% 为每个变量创建对比图
for col = 1:19% 创建新图形（包含两个子图）figure('Position', [100, 100, 1400, 600]);% 子图1：降噪前subplot(1, 2, 1);hold on;plot_idx = 1;for i = datasets_to_processdata = aligned_data{i};plot(time_axis, data(:, col), 'LineWidth', 1.5, ...'Color', colors_original(plot_idx, :), ...'DisplayName', ['数据集 ' num2str(i)]);plot_idx = plot_idx + 1;endxlabel('时间步 (1-85)', 'FontSize', 11, 'FontWeight', 'bold');ylabel('数值', 'FontSize', 11, 'FontWeight', 'bold');title(['降噪前 - ' var_names{col}], 'FontSize', 13, 'FontWeight', 'bold');legend('Location', 'best', 'FontSize', 9, 'NumColumns', 2);grid on;xlim([1, target_length]);hold off;% 子图2：降噪后subplot(1, 2, 2);hold on;plot_idx = 1;for i = datasets_to_processdata = denoised_data{i};plot(time_axis, data(:, col), 'LineWidth', 1.8, ...'Color', colors_denoised(plot_idx, :), ...'DisplayName', ['数据集 ' num2str(i)], ...'LineStyle', '-');plot_idx = plot_idx + 1;endxlabel('时间步 (1-85)', 'FontSize', 11, 'FontWeight', 'bold');ylabel('数值', 'FontSize', 11, 'FontWeight', 'bold');title(['降噪后 - ' var_names{col}], 'FontSize', 13, 'FontWeight', 'bold');legend('Location', 'best', 'FontSize', 9, 'NumColumns', 2);grid on;xlim([1, target_length]);hold off;% 添加总标题sgtitle(['变量 ' num2str(col) ': ' var_names{col} ' - 降噪对比'], ...'FontSize', 15, 'FontWeight', 'bold');% 保存图像saveas(gcf, ['降噪对比_变量' num2str(col) '.fig']);saveas(gcf, ['降噪对比_变量' num2str(col) '.png']);close(gcf);fprintf('已完成变量 %d/%d 的对比图\n', col, 19);
end%% ========== 创建单独的降噪后图像（与原格式一致） ==========
fprintf('\n开始生成降噪后的单独图像...\n');colors = lines(9);
for col = 1:19figure('Position', [100, 100, 1000, 600]);hold on;plot_idx = 1;for i = datasets_to_processdata = denoised_data{i};plot(time_axis, data(:, col), 'LineWidth', 1.8, ...'Color', colors(plot_idx, :), ...'DisplayName', ['数据集 ' num2str(i)], ...'Marker', 'none');plot_idx = plot_idx + 1;endxlabel('时间步 (1-85)', 'FontSize', 12, 'FontWeight', 'bold');ylabel('数值', 'FontSize', 12, 'FontWeight', 'bold');title(['变量 ' num2str(col) ': ' var_names{col} ' (降噪后)'], ...'FontSize', 14, 'FontWeight', 'bold');legend('Location', 'best', 'FontSize', 10, 'NumColumns', 2);grid on;xlim([1, target_length]);hold off;saveas(gcf, ['降噪后_变量' num2str(col) '.fig']);saveas(gcf, ['降噪后_变量' num2str(col) '.png']);close(gcf);fprintf('已完成变量 %d/%d 的降噪后图像\n', col, 19);
end%% ========== 保存降噪后的数据为.dat和.xlsx格式 ==========
fprintf('\n开始保存降噪后的数据文件...\n');for idx = 1:length(datasets_to_process)i = datasets_to_process(idx);data_to_save = denoised_data{i};% 生成文件名（06到14）file_number = sprintf('%02d', i + 2);dat_filename = [file_number '_te.dat'];xlsx_filename = [file_number '_te.xlsx'];% 保存为.dat文件（空格分隔）dlmwrite(dat_filename, data_to_save, 'delimiter', ' ', 'precision', '%.8f');fprintf('已保存: %s\n', dat_filename);% 保存为.xlsx文件col_headers = var_names;if exist('writematrix', 'file')% MATLAB R2019a及以上writecell([col_headers; num2cell(data_to_save)], xlsx_filename);else% 旧版本MATLABxlswrite(xlsx_filename, col_headers, 1, 'A1');xlswrite(xlsx_filename, data_to_save, 1, 'A2');endfprintf('已保存: %s\n', xlsx_filename);
end%% ========== 保存降噪后的完整数据 ==========
save('降噪后的数据_85点.mat', 'denoised_data', 'time_axis', 'method', ...'sg_order', 'sg_framelen');
fprintf('\n已保存: 降噪后的数据_85点.mat\n');%% ========== 生成降噪效果评估报告 ==========
fprintf('\n========== 降噪效果评估 ==========\n');% 计算信噪比改善
snr_improvement = zeros(length(datasets_to_process), 19);
smoothness_improvement = zeros(length(datasets_to_process), 19);for idx = 1:length(datasets_to_process)i = datasets_to_process(idx);original = aligned_data{i};denoised = denoised_data{i};for col = 1:19% 计算噪声水平（使用一阶差分的标准差）noise_original = std(diff(original(:, col)));noise_denoised = std(diff(denoised(:, col)));% 噪声减少百分比snr_improvement(idx, col) = (1 - noise_denoised/noise_original) * 100;% 平滑度改善（二阶差分）smooth_original = std(diff(original(:, col), 2));smooth_denoised = std(diff(denoised(:, col), 2));smoothness_improvement(idx, col) = (1 - smooth_denoised/smooth_original) * 100;end
endfprintf('平均噪声减少: %.2f%%\n', mean(snr_improvement(:)));
fprintf('平均平滑度改善: %.2f%%\n', mean(smoothness_improvement(:)));% 创建评估热图
figure('Position', [100, 100, 1200, 500]);
subplot(1, 2, 1);
imagesc(snr_improvement');
colorbar;
xlabel('数据集编号', 'FontSize', 11);
ylabel('变量编号', 'FontSize', 11);
title('噪声减少百分比 (%)', 'FontSize', 13, 'FontWeight', 'bold');
set(gca, 'XTick', 1:9, 'XTickLabel', datasets_to_process, 'YTick', 1:19);
colormap(jet);subplot(1, 2, 2);
imagesc(smoothness_improvement');
colorbar;
xlabel('数据集编号', 'FontSize', 11);
ylabel('变量编号', 'FontSize', 11);
title('平滑度改善百分比 (%)', 'FontSize', 13, 'FontWeight', 'bold');
set(gca, 'XTick', 1:9, 'XTickLabel', datasets_to_process, 'YTick', 1:19);
colormap(jet);sgtitle('降噪效果评估', 'FontSize', 15, 'FontWeight', 'bold');
saveas(gcf, '降噪效果评估.fig');
saveas(gcf, '降噪效果评估.png');
close(gcf);fprintf('\n========== 所有处理完成！ ==========\n');
fprintf('总共生成:\n');
fprintf('  - %d 个对比图（.fig和.png）\n', 19);
fprintf('  - %d 个降噪后图像（.fig和.png）\n', 19);
fprintf('  - %d 个.dat文件\n', length(datasets_to_process));
fprintf('  - %d 个.xlsx文件\n', length(datasets_to_process));
fprintf('  - 1 个降噪后的数据.mat文件\n');
fprintf('  - 1 个降噪效果评估图\n');