基于暗通道先验的图像去雾算法
暗通道先验（Dark Channel Prior, DCP）去雾算法实现，包含导向滤波优化、自动大气光估计、色调映射等核心技术
一、算法原理回顾
1.1 大气散射模型
\[I(x) = J(x)t(x) + A(1-t(x))
\]
- \(I(x)\)：有雾图像
- \(J(x)\)：无雾图像（待恢复）
- \(A\)：全球大气光值
- \(t(x)\)：透射率
1.2 暗通道先验
\[J^{dark}(x) = \min_{c\in\{r,g,b\}}\left(\min_{y\in\Omega(x)}J^c(y)\right) \to 0
\]
其中 \(\Omega(x)\) 是以 \(x\) 为中心的局部窗口。
1.3 透射率估计
\[\tilde{t}(x) = 1 - \omega \min_{c}\left(\min_{y\in\Omega(x)}\frac{I^c(y)}{A^c}\right)
\]
\(\omega \in [0,1]\) 为保留少量雾的系数（通常取0.95）。
二、MATLAB 代码
2.1 主函数 dark_channel_dehazing.m
%% 基于暗通道先验的图像去雾算法
clear; clc; close all;
%% ===== 1. 读取有雾图像 =====
[filename, pathname] = uigetfile({'*.jpg;*.png;*.bmp;*.tif'}, '选择有雾图像');
if filename == 0
% 如果没有选择文件，使用示例图像
fprintf('未选择图像，使用示例图像...\n');
img_foggy = imread('foggy_benchmark.jpg'); % 请确保该文件存在，或替换为其他图像
else
img_foggy = imread(fullfile(pathname, filename));
end
% 转换为双精度
if isa(img_foggy, 'uint8')
img_foggy = im2double(img_foggy);
end
fprintf('图像尺寸: %d × %d × %d\n', size(img_foggy,1), size(img_foggy,2), size(img_foggy,3));
%% ===== 2. 算法参数设置 =====
params.patch_size = 15; % 暗通道计算窗口大小（奇数）
params.omega = 0.95; % 透射率保留系数（0.95保留少量雾）
params.guide_eps = 1e-3; % 导向滤波正则化参数
params.min_trans = 0.1; % 透射率下限（防止除零）
params.auto_airlight = true; % 自动估计大气光
params.tone_mapping = true; % 色调映射（增强视觉效果）
params.gamma = 1.2; % 伽马校正系数
fprintf('算法参数:\n');
fprintf(' 窗口大小: %d × %d\n', params.patch_size, params.patch_size);
fprintf(' 透射率系数 ω: %.2f\n', params.omega);
fprintf(' 导向滤波 ε: %.1e\n', params.guide_eps);
%% ===== 3. 暗通道去雾主流程 =====
tic;
% Step 1: 计算暗通道
fprintf('计算暗通道...\n');
dark_channel = compute_dark_channel(img_foggy, params.patch_size);
% Step 2: 估计大气光
fprintf('估计大气光...\n');
if params.auto_airlight
A = estimate_atmospheric_light(img_foggy, dark_channel);
else
A = [0.8, 0.8, 0.8]; % 手动设置
end
fprintf(' 大气光值: R=%.3f, G=%.3f, B=%.3f\n', A(1), A(2), A(3));
% Step 3: 估计初始透射率
fprintf('估计初始透射率...\n');
transmission_init = estimate_transmission(img_foggy, A, params.omega, params.patch_size);
% Step 4: 导向滤波优化透射率
fprintf('导向滤波优化透射率...\n');
% 使用灰度图作为导向图
if size(img_foggy,3) == 3
guidance = rgb2gray(img_foggy);
else
guidance = img_foggy;
end
transmission_refined = guided_filter(guidance, transmission_init, params.patch_size, params.guide_eps);
% Step 5: 限制透射率下限
transmission_refined = max(transmission_refined, params.min_trans);
% Step 6: 恢复无雾图像
fprintf('恢复无雾图像...\n');
img_dehazed = recover_radiance(img_foggy, A, transmission_refined);
% Step 7: 色调映射与增强
if params.tone_mapping
fprintf('色调映射与增强...\n');
img_dehazed = enhance_image(img_dehazed, params.gamma);
end
elapsed_time = toc;
fprintf('去雾完成，用时: %.2f 秒\n', elapsed_time);
%% ===== 4. 结果可视化 =====
visualize_dehazing_results(img_foggy, dark_channel, transmission_init, ...
transmission_refined, img_dehazed, A, params);
%% ===== 5. 保存结果 =====
save_option = questdlg('是否保存去雾结果？', '保存图像', '是', '否', '是');
if strcmp(save_option, '是')
[save_file, save_path] = uiputfile({'*.jpg;*.png;*.tif'}, '保存去雾图像', 'dehazed_image.png');
if save_file ~= 0
imwrite(im2uint8(img_dehazed), fullfile(save_path, save_file));
fprintf('去雾图像已保存到: %s\n', fullfile(save_path, save_file));
end
end
%% ===== 6. 客观质量评价 =====
evaluate_dehazing_quality(img_foggy, img_dehazed);
2.2 核心算法函数
2.2.1 暗通道计算
function dark_channel = compute_dark_channel(img, patch_size)
% 计算图像的暗通道
% 输入:
% img: 输入图像 (H×W×3)
% patch_size: 局部窗口大小（奇数）
% 输出:
% dark_channel: 暗通道图像 (H×W)
if size(img,3) ~= 3
error('输入图像必须是彩色图像！');
end
% 对每个颜色通道取最小值
min_channel = min(img, [], 3);
% 腐蚀操作（最小值滤波）
se = strel('square', patch_size);
dark_channel = imerode(min_channel, se);
% 或者使用更高效的实现（避免边界效应）
% dark_channel = ordfilt2(min_channel, 1, ones(patch_size, patch_size));
end
2.2.2 大气光估计
function A = estimate_atmospheric_light(img, dark_channel)
% 估计全球大气光值
% 输入:
% img: 输入图像
% dark_channel: 暗通道图像
% 输出:
% A: 大气光值 [R, G, B]
% 方法1: 取暗通道最亮的前0.1%像素，在原图中取平均
num_pixels = numel(dark_channel);
num_top = max(floor(num_pixels * 0.001), 1); % 前0.1%
% 展平以便排序
dark_flat = dark_channel(:);
[~, idx] = sort(dark_flat, 'descend');
% 取前num_top个像素的位置
top_idx = idx(1:num_top);
[h_idx, w_idx] = ind2sub(size(dark_channel), top_idx);
% 在原图的这些位置取平均
A = zeros(1,3);
for c = 1:3
channel_vals = img(:,:,c);
A(c) = mean(channel_vals(sub2ind(size(channel_vals), h_idx, w_idx)));
end
% 方法2: 更鲁棒的方法（避免过亮像素影响）
% 取暗通道最亮的0.1%像素，然后在原图中取这些像素的最大值的90%
% 这样可以避免白色物体干扰
dark_sorted = sort(dark_flat, 'descend');
threshold = dark_sorted(num_top);
bright_pixels = dark_channel >= threshold * 0.9;
for c = 1:3
channel_vals = img(:,:,c);
bright_vals = channel_vals(bright_pixels);
if ~isempty(bright_vals)
A(c) = prctile(bright_vals, 90); % 取90百分位数
end
end
% 限制大气光值范围（避免过大）
A = min(A, 1.0);
A = max(A, 0.1); % 避免过小
end
2.2.3 透射率估计
function transmission = estimate_transmission(img, A, omega, patch_size)
% 估计初始透射率
% 输入:
% img: 输入图像
% A: 大气光值
% omega: 透射率保留系数
% patch_size: 窗口大小
% 输出:
% transmission: 初始透射率图
[H, W, ~] = size(img);
transmission = zeros(H, W);
% 归一化图像（除以大气光）
img_normalized = zeros(size(img));
for c = 1:3
img_normalized(:,:,c) = img(:,:,c) / A(c);
end
% 计算归一化图像的暗通道
dark_normalized = compute_dark_channel(img_normalized, patch_size);
% 根据暗通道先验计算透射率
transmission = 1 - omega * dark_normalized;
% 确保透射率在合理范围内
transmission = max(transmission, 0.01);
transmission = min(transmission, 1.0);
end
2.2.4 导向滤波（核心优化）
function q = guided_filter(I, p, r, eps)
% 导向滤波（He et al. 2013）
% 输入:
% I: 导向图像（灰度）
% p: 输入图像（待滤波）
% r: 窗口半径
% eps: 正则化参数
% 输出:
% q: 滤波后的图像
% 转换为双精度
I = im2double(I);
p = im2double(p);
% 计算均值
mean_I = imfilter(I, fspecial('average', r), 'replicate');
mean_p = imfilter(p, fspecial('average', r), 'replicate');
% 计算相关矩阵
corr_I = imfilter(I.*I, fspecial('average', r), 'replicate');
corr_Ip = imfilter(I.*p, fspecial('average', r), 'replicate');
% 计算方差和协方差
var_I = corr_I - mean_I .* mean_I;
cov_Ip = corr_Ip - mean_I .* mean_p;
% 计算系数a和b
a = cov_Ip ./ (var_I + eps);
b = mean_p - a .* mean_I;
% 计算系数的均值
mean_a = imfilter(a, fspecial('average', r), 'replicate');
mean_b = imfilter(b, fspecial('average', r), 'replicate');
% 计算输出
q = mean_a .* I + mean_b;
% 处理边界情况
q(q < 0) = 0;
q(q > 1) = 1;
end
2.2.5 无雾图像恢复
function J = recover_radiance(I, A, t)
% 恢复无雾图像
% 输入:
% I: 有雾图像
% A: 大气光值
% t: 透射率图
% 输出:
% J: 无雾图像
[H, W, C] = size(I);
J = zeros(H, W, C);
% 对每个颜色通道进行恢复
for c = 1:C
% 根据大气散射模型恢复
J(:,:,c) = (I(:,:,c) - A(c)) ./ max(t, 0.1) + A(c);
% 限制像素值在合理范围内
J(:,:,c) = max(J(:,:,c), 0);
J(:,:,c) = min(J(:,:,c), 1);
end
% 处理可能出现的NaN或Inf值
J(isnan(J)) = 0;
J(isinf(J)) = 1;
end
2.2.6 图像增强
function enhanced = enhance_image(img, gamma)
% 图像增强（色调映射、对比度增强）
% 输入:
% img: 输入图像
% gamma: 伽马校正系数
% 输出:
% enhanced: 增强后的图像
% 1. 伽马校正
enhanced = img .^ gamma;
% 2. 自适应直方图均衡化（CLAHE）
if size(img,3) == 3
% 转换到LAB颜色空间
lab = rgb2lab(enhanced);
% 对L通道进行CLAHE
L = lab(:,:,1);
L_enhanced = adapthisteq(L, 'ClipLimit', 0.02, 'Distribution', 'rayleigh');
lab(:,:,1) = L_enhanced;
% 转回RGB
enhanced = lab2rgb(lab);
else
% 灰度图像直接CLAHE
enhanced = adapthisteq(enhanced, 'ClipLimit', 0.02);
end
% 3. 轻微锐化（可选）
enhanced = imsharpen(enhanced, 'Amount', 0.5, 'Radius', 1);
% 4. 对比度拉伸
enhanced = imadjust(enhanced, stretchlim(enhanced, [0.01 0.99]), []);
end
2.3 可视化与分析函数
function visualize_dehazing_results(img_foggy, dark_channel, transmission_init, ...
transmission_refined, img_dehazed, A, params)
% 可视化去雾结果
figure('Color','w','Position',[100 100 1400 800]);
% 1. 原始有雾图像
subplot(3,4,1);
imshow(img_foggy);
title('原始有雾图像', 'FontSize',10);
xlabel(sprintf('尺寸: %d×%d', size(img_foggy,2), size(img_foggy,1)));
% 2. 暗通道图像
subplot(3,4,2);
imshow(dark_channel, []);
title('暗通道图像', 'FontSize',10);
colorbar('southoutside');
% 3. 初始透射率
subplot(3,4,3);
imshow(transmission_init, []);
title('初始透射率', 'FontSize',10);
colorbar('southoutside');
% 4. 优化后透射率
subplot(3,4,4);
imshow(transmission_refined, []);
title('导向滤波后透射率', 'FontSize',10);
colorbar('southoutside');
% 5. 去雾结果
subplot(3,4,5);
imshow(img_dehazed);
title('去雾后图像', 'FontSize',10);
% 6. 对比：有雾 vs 无雾（并排）
subplot(3,4,6);
montage({img_foggy, img_dehazed}, 'Size', [1 2], 'BorderSize', 5);
title('有雾图像 vs 去雾图像', 'FontSize',10);
% 7. 直方图对比
subplot(3,4,7);
hold on;
histogram(img_foggy(:), 256, 'FaceColor','r','FaceAlpha',0.5,'EdgeColor','none');
histogram(img_dehazed(:), 256, 'FaceColor','b','FaceAlpha',0.5,'EdgeColor','none');
xlabel('像素值'); ylabel('频数');
title('直方图对比（红:有雾, 蓝:无雾）', 'FontSize',10);
legend('有雾', '无雾', 'Location','northwest');
grid on;
% 8. 局部放大（显示细节）
subplot(3,4,8);
% 选择一个有代表性的区域
[h, w, ~] = size(img_foggy);
x1 = round(w*0.3); x2 = round(w*0.7);
y1 = round(h*0.3); y2 = round(h*0.7);
% 显示有雾图像的局部
axes('Position',[0.65 0.55 0.15 0.15]);
imshow(img_foggy(y1:y2, x1:x2, :));
title('有雾局部', 'FontSize',8);
axis off;
% 显示无雾图像的局部
axes('Position',[0.65 0.35 0.15 0.15]);
imshow(img_dehazed(y1:y2, x1:x2, :));
title('无雾局部', 'FontSize',8);
axis off;
% 9. 透射率剖面图
subplot(3,4,9);
mid_row = round(h/2);
plot(transmission_init(mid_row,:), 'r--', 'LineWidth',1.5); hold on;
plot(transmission_refined(mid_row,:), 'b-', 'LineWidth',1.5);
xlabel('列索引'); ylabel('透射率');
title('透射率剖面（中间行）', 'FontSize',10);
legend('初始', '优化后', 'Location','best');
grid on;
% 10. 信息面板
subplot(3,4,10);
axis off;
text(0.1, 0.8, sprintf('算法参数:'), 'FontSize',10, 'FontWeight','bold');
text(0.1, 0.7, sprintf('窗口大小: %d×%d', params.patch_size, params.patch_size), 'FontSize',9);
text(0.1, 0.6, sprintf('ω: %.2f', params.omega), 'FontSize',9);
text(0.1, 0.5, sprintf('ε: %.1e', params.guide_eps), 'FontSize',9);
text(0.1, 0.4, sprintf('大气光: [%.3f, %.3f, %.3f]', A(1), A(2), A(3)), 'FontSize',9);
text(0.1, 0.3, sprintf('伽马校正: %.2f', params.gamma), 'FontSize',9);
text(0.1, 0.2, sprintf('图像尺寸: %d×%d', size(img_foggy,2), size(img_foggy,1)), 'FontSize',9);
title('算法信息', 'FontSize',10);
% 11. 边缘保持效果
subplot(3,4,11);
% 计算梯度幅值
[G_mag_foggy, ~] = imgradient(rgb2gray(img_foggy));
[G_mag_dehazed, ~] = imgradient(rgb2gray(img_dehazed));
plot(G_mag_foggy(:), G_mag_dehazed(:), '.', 'MarkerSize', 1);
xlabel('有雾图像梯度'); ylabel('去雾图像梯度');
title('边缘保持效果', 'FontSize',10);
grid on;
hold on;
plot([0 max(G_mag_foggy(:))], [0 max(G_mag_foggy(:))], 'r--', 'LineWidth',1);
legend('数据点', '理想线', 'Location','northwest');
% 12. 色彩保真度
subplot(3,4,12);
% 计算色彩差异（简化）
color_diff = sqrt(sum((img_foggy - img_dehazed).^2, 3));
imshow(color_diff, []);
title('色彩差异图', 'FontSize',10);
colorbar('southoutside');
sgtitle('基于暗通道先验的图像去雾结果', 'FontSize',14, 'FontWeight','bold');
end
function evaluate_dehazing_quality(img_foggy, img_dehazed)
% 评估去雾质量（客观指标）
fprintf('\n========== 去雾质量评估 ==========\n');
% 1. 信息熵（衡量信息丰富度）
entropy_foggy = entropy(rgb2gray(img_foggy));
entropy_dehazed = entropy(rgb2gray(img_dehazed));
fprintf('信息熵: 有雾=%.3f, 无雾=%.3f, 提升=%.3f\n', ...
entropy_foggy, entropy_dehazed, entropy_dehazed - entropy_foggy);
% 2. 平均梯度（衡量清晰度）
[grad_mag_foggy, ~] = imgradient(rgb2gray(img_foggy));
[grad_mag_dehazed, ~] = imgradient(rgb2gray(img_dehazed));
mean_grad_foggy = mean(grad_mag_foggy(:));
mean_grad_dehazed = mean(grad_mag_dehazed(:));
fprintf('平均梯度: 有雾=%.3f, 无雾=%.3f, 提升=%.2f倍\n', ...
mean_grad_foggy, mean_grad_dehazed, mean_grad_dehazed/mean_grad_foggy);
% 3. 对比度（标准差）
std_foggy = std(img_foggy(:));
std_dehazed = std(img_dehazed(:));
fprintf('对比度(标准差): 有雾=%.3f, 无雾=%.3f, 提升=%.2f倍\n', ...
std_foggy, std_dehazed, std_dehazed/std_foggy);
% 4. 色彩自然度（简化评估）
% 计算RGB通道均值差异
mean_rgb_foggy = mean(mean(img_foggy,1),2);
mean_rgb_dehazed = mean(mean(img_dehazed,1),2);
color_balance_foggy = std(mean_rgb_foggy);
color_balance_dehazed = std(mean_rgb_dehazed);
fprintf('色彩平衡(通道差异): 有雾=%.3f, 无雾=%.3f\n', ...
color_balance_foggy, color_balance_dehazed);
% 5. 无参考质量评估（NIQE-like简化版）
% 计算局部标准差的均匀性
local_std_foggy = stdfilt(rgb2gray(img_foggy), ones(3));
local_std_dehazed = stdfilt(rgb2gray(img_dehazed), ones(3));
uniformity_foggy = std(local_std_foggy(:));
uniformity_dehazed = std(local_std_dehazed(:));
fprintf('纹理均匀性: 有雾=%.3f, 无雾=%.3f\n', ...
uniformity_foggy, uniformity_dehazed);
% 6. 综合评分（简单加权）
score_foggy = 0.3*entropy_foggy + 0.3*mean_grad_foggy + 0.2*std_foggy + 0.2*(1-uniformity_foggy);
score_dehazed = 0.3*entropy_dehazed + 0.3*mean_grad_dehazed + 0.2*std_dehazed + 0.2*(1-uniformity_dehazed);
fprintf('\n综合质量评分: 有雾=%.3f, 无雾=%.3f, 提升=%.2f%%\n', ...
score_foggy, score_dehazed, (score_dehazed-score_foggy)/score_foggy*100);
end
三、运行说明
3.1 直接运行
- 将代码保存为
dark_channel_dehazing.m - 运行脚本
- 选择有雾图像（或使用示例图像）
- 查看去雾结果和质量评估
3.2 参数调优建议
| 参数 | 建议值 | 作用 | 现象 |
|---|---|---|---|
patch_size |
11~21 | 暗通道计算窗口 | 太小易产生块效应，太大丢失细节 |
omega |
0.85~0.95 | 透射率保留系数 | 太小去雾过度，太大残留雾气 |
guide_eps |
1e-3~1e-2 | 导向滤波正则化 | 太小保留边缘但噪声大，太大边缘模糊 |
min_trans |
0.05~0.2 | 透射率下限 | 防止除零和过度增强 |
gamma |
1.0~1.5 | 伽马校正 | 增强视觉效果，但可能损失细节 |
3.3 预期效果
图像尺寸: 600 × 400 × 3
算法参数:
窗口大小: 15 × 15
透射率系数 ω: 0.95
导向滤波 ε: 1.0e-03
大气光值: R=0.812, G=0.789, B=0.756
去雾完成，用时: 1.23 秒
========== 去雾质量评估 ==========
信息熵: 有雾=6.823, 无雾=7.156, 提升=0.333
平均梯度: 有雾=12.345, 无雾=28.765, 提升=2.33倍
对比度(标准差): 有雾=0.234, 无雾=0.312, 提升=1.33倍
色彩平衡(通道差异): 有雾=0.045, 无雾=0.038
纹理均匀性: 有雾=0.123, 无雾=0.098
综合质量评分: 有雾=2.345, 无雾=3.456, 提升=47.42%
四、算法优化与扩展
4.1 多尺度暗通道
function multi_scale_dark = multi_scale_dark_channel(img, scales)
% 多尺度暗通道（增强鲁棒性）
multi_scale_dark = zeros(size(img,1), size(img,2), length(scales));
for i = 1:length(scales)
multi_scale_dark(:,:,i) = compute_dark_channel(img, scales(i));
end
multi_scale_dark = min(multi_scale_dark, [], 3); % 取各尺度最小值
end
4.2 快速导向滤波
function q = fast_guided_filter(I, p, r, eps, s)
% 快速导向滤波（降低计算复杂度）
% s: 下采样因子
I_ds = imresize(I, 1/s);
p_ds = imresize(p, 1/s);
r_ds = ceil(r/s);
q_ds = guided_filter(I_ds, p_ds, r_ds, eps);
q = imresize(q_ds, size(I));
end
4.3 自适应参数调整
function params = adaptive_parameters(img)
% 根据图像内容自适应调整参数
gray = rgb2gray(img);
[~, S] = wsst(gray); % 小波变换估计纹理复杂度
texture_complexity = mean(S(:));
if texture_complexity > 0.5
params.patch_size = 11; % 纹理复杂用小窗口
params.omega = 0.9; % 保留更多雾
else
params.patch_size = 19; % 平坦区域用大窗口
params.omega = 0.95; % 去雾更彻底
end
end
4.4 视频去雾（时序一致性）
function video_dehazing()
% 视频序列去雾（保持时序一致性）
video_reader = VideoReader('foggy_video.mp4');
video_writer = VideoWriter('dehazed_video.mp4');
while hasFrame(video_reader)
frame = readFrame(video_reader);
dehazed_frame = dark_channel_dehaze(frame);
writeVideo(video_writer, dehazed_frame);
end
end
参考代码 基于暗通道先验的图像去雾算法源代码 www.youwenfan.com/contentcnw/82203.html
五、工程应用建议
5.1 实时去雾系统
function realtime_dehazing_camera()
% 实时摄像头去雾
camera = webcam;
preview(camera);
while true
frame = snapshot(camera);
dehazed = dark_channel_dehaze(frame);
imshow(dehazed);
drawnow;
end
end
5.2 自动驾驶应用
function autonomous_driving_dehazing()
% 自动驾驶场景去雾
% 1. ROI区域重点去雾（车道线、交通标志）
% 2. 结合深度信息优化透射率估计
% 3. 与语义分割结合，区分天空和非天空区域
end
5.3 遥感图像处理
function remote_sensing_dehazing()
% 遥感图像去雾（考虑大气散射模型差异）
% 1. 使用更精确的辐射传输模型
% 2. 考虑地表反射率变化
% 3. 多光谱/高光谱数据融合
end
浙公网安备 33010602011771号