基于DCT的彩色图像压缩技术详解

基于离散余弦变换（DCT）的图像压缩是JPEG等标准的核心技术，特别适合处理彩色图像。

DCT压缩核心原理

1. 基本流程

flowchart TD A[原始彩色图像] --> B[颜色空间转换] B --> C[分块处理8×8] C --> D[DCT变换] D --> E[量化] E --> F[熵编码] F --> G[压缩数据] G --> H[熵解码] H --> I[反量化] I --> J[逆DCT变换] J --> K[块合并] K --> L[颜色空间逆转换] L --> M[重建图像]

完整MATLAB实现

1. 主函数 - DCT彩色图像压缩

function [compressed_data, reconstructed_img, compression_info] = dct_color_compress(img, quality_factor)
    % 基于DCT的彩色图像压缩
    % 输入:
    %   img - 输入彩色图像
    %   quality_factor - 质量因子 (1-100)
    % 输出:
    %   compressed_data - 压缩后的数据
    %   reconstructed_img - 重建图像
    %   compression_info - 压缩信息
    
    % 参数验证
    if nargin < 2
        quality_factor = 50;
    end
    quality_factor = max(1, min(100, quality_factor));
    
    % 转换为double类型
    original_img = im2double(img);
    
    % 颜色空间转换: RGB -> YCbCr
    ycbcr_img = rgb2ycbcr(original_img);
    
    % 提取Y, Cb, Cr分量
    Y = ycbcr_img(:, :, 1);
    Cb = ycbcr_img(:, :, 2);
    Cr = ycbcr_img(:, :, 3);
    
    % 色度下采样 (4:2:0)
    [Cb_subsampled, Cr_subsampled] = chroma_subsample(Cb, Cr);
    
    % 对每个分量进行DCT压缩
    fprintf('正在压缩Y分量...\n');
    [compressed_Y, info_Y] = dct_compress_channel(Y, quality_factor, 'luminance');
    
    fprintf('正在压缩Cb分量...\n');
    [compressed_Cb, info_Cb] = dct_compress_channel(Cb_subsampled, quality_factor, 'chrominance');
    
    fprintf('正在压缩Cr分量...\n');
    [compressed_Cr, info_Cr] = dct_compress_channel(Cr_subsampled, quality_factor, 'chrominance');
    
    % 组合压缩数据
    compressed_data.Y = compressed_Y;
    compressed_data.Cb = compressed_Cb;
    compressed_data.Cr = compressed_Cr;
    compressed_data.original_size = size(img);
    compressed_data.quality_factor = quality_factor;
    
    % 重建图像
    fprintf('重建图像...\n');
    reconstructed_img = dct_color_decompress(compressed_data);
    
    % 计算压缩信息
    compression_info = calculate_compression_info(original_img, compressed_data, info_Y, info_Cb, info_Cr);
    
    fprintf('压缩完成!\n');
end

2. 色度下采样与上采样

function [Cb_subsampled, Cr_subsampled] = chroma_subsample(Cb, Cr)
    % 色度下采样 (4:2:0)
    % 在水平和垂直方向都进行2:1下采样
    
    % 使用抗混叠滤波
    h = fspecial('gaussian', [3 3], 0.5);
    Cb_filtered = imfilter(Cb, h, 'replicate');
    Cr_filtered = imfilter(Cr, h, 'replicate');
    
    % 下采样
    Cb_subsampled = Cb_filtered(1:2:end, 1:2:end);
    Cr_subsampled = Cr_filtered(1:2:end, 1:2:end);
end

function [Cb_upsampled, Cr_upsampled] = chroma_upsample(Cb_subsampled, Cr_subsampled, target_size)
    % 色度上采样
    % 使用双线性插值恢复原始尺寸
    
    Cb_upsampled = imresize(Cb_subsampled, target_size, 'bilinear');
    Cr_upsampled = imresize(Cr_subsampled, target_size, 'bilinear');
end

3. 单通道DCT压缩

function [compressed_channel, info] = dct_compress_channel(channel, quality_factor, channel_type)
    % 对单个通道进行DCT压缩
    % 输入:
    %   channel - 输入通道数据
    %   quality_factor - 质量因子
    %   channel_type - 'luminance' 或 'chrominance'
    
    [height, width] = size(channel);
    block_size = 8;
    
    % 确保图像尺寸是8的倍数
    padded_height = ceil(height / block_size) * block_size;
    padded_width = ceil(width / block_size) * block_size;
    channel_padded = padarray(channel, [padded_height-height, padded_width-width], 'replicate', 'post');
    
    % 获取量化表
    quant_table = get_quantization_table(quality_factor, channel_type);
    
    % 分块处理
    compressed_blocks = cell(padded_height/block_size, padded_width/block_size);
    dc_coefficients = [];
    ac_coefficients = [];
    
    for i = 1:block_size:padded_height
        for j = 1:block_size:padded_width
            % 提取8x8块
            block = channel_padded(i:i+block_size-1, j:j+block_size-1);
            
            % 电平偏移 (0-255 -> -128 to 127)
            block_shifted = block - 128;
            
            % DCT变换
            dct_coeffs = dct2(block_shifted);
            
            % 量化
            quantized_coeffs = round(dct_coeffs ./ quant_table);
            
            % 保存压缩后的块
            compressed_blocks{(i-1)/block_size+1, (j-1)/block_size+1} = quantized_coeffs;
            
            % 提取DC和AC系数用于熵编码分析
            dc_coefficients = [dc_coefficients; quantized_coeffs(1,1)];
            ac_block = quantized_coeffs(:);
            ac_block(1) = []; % 移除DC系数
            ac_coefficients = [ac_coefficients; ac_block];
        end
    end
    
    % 组织压缩数据
    compressed_channel.blocks = compressed_blocks;
    compressed_channel.quant_table = quant_table;
    compressed_channel.original_size = size(channel);
    compressed_channel.padded_size = [padded_height, padded_width];
    compressed_channel.block_size = block_size;
    
    % 统计信息
    info.dc_coefficients = dc_coefficients;
    info.ac_coefficients = ac_coefficients;
    info.non_zero_ac = sum(ac_coefficients ~= 0);
    info.total_ac = length(ac_coefficients);
end

4. 量化表生成

function quant_table = get_quantization_table(quality_factor, channel_type)
    % 生成量化表
    % 基于JPEG标准量化表
    
    % 标准亮度量化表
    luminance_quant = [16  11  10  16  24  40  51  61;
                       12  12  14  19  26  58  60  55;
                       14  13  16  24  40  57  69  56;
                       14  17  22  29  51  87  80  62;
                       18  22  37  56  68  109 103 77;
                       24  35  55  64  81  104 113 92;
                       49  64  78  87  103 121 120 101;
                       72  92  95  98  112 100 103 99];
    
    % 标准色度量化表
    chrominance_quant = [17  18  24  47  99  99  99  99;
                         18  21  26  66  99  99  99  99;
                         24  26  56  99  99  99  99  99;
                         47  66  99  99  99  99  99  99;
                         99  99  99  99  99  99  99  99;
                         99  99  99  99  99  99  99  99;
                         99  99  99  99  99  99  99  99;
                         99  99  99  99  99  99  99  99];
    
    % 选择量化表
    if strcmpi(channel_type, 'luminance')
        base_quant = luminance_quant;
    else
        base_quant = chrominance_quant;
    end
    
    % 根据质量因子调整量化表
    if quality_factor < 50
        scale_factor = 50 / quality_factor;
    else
        scale_factor = (100 - quality_factor) / 50;
    end
    
    quant_table = max(1, min(255, round(base_quant * scale_factor)));
end

5. 解压缩函数

function reconstructed_img = dct_color_decompress(compressed_data)
    % DCT彩色图像解压缩
    
    % 提取压缩数据
    compressed_Y = compressed_data.Y;
    compressed_Cb = compressed_data.Cb;
    compressed_Cr = compressed_data.Cr;
    original_size = compressed_data.original_size;
    
    % 解压缩各分量
    fprintf('解压缩Y分量...\n');
    Y_reconstructed = dct_decompress_channel(compressed_Y);
    
    fprintf('解压缩Cb分量...\n');
    Cb_subsampled = dct_decompress_channel(compressed_Cb);
    
    fprintf('解压缩Cr分量...\n');
    Cr_subsampled = dct_decompress_channel(compressed_Cr);
    
    % 色度上采样
    [Cb_upsampled, Cr_upsampled] = chroma_upsample(Cb_subsampled, Cr_subsampled, original_size(1:2));
    
    % 组合YCbCr图像
    ycbcr_reconstructed = cat(3, Y_reconstructed, Cb_upsampled, Cr_upsampled);
    
    % 转换回RGB
    reconstructed_img = ycbcr2rgb(ycbcr_reconstructed);
    
    % 裁剪到原始尺寸
    reconstructed_img = reconstructed_img(1:original_size(1), 1:original_size(2), :);
end

function channel_reconstructed = dct_decompress_channel(compressed_channel)
    % 单通道DCT解压缩
    
    blocks = compressed_channel.blocks;
    quant_table = compressed_channel.quant_table;
    block_size = compressed_channel.block_size;
    padded_size = compressed_channel.padded_size;
    
    % 重建图像
    channel_reconstructed = zeros(padded_size);
    
    for i = 1:block_size:padded_size(1)
        for j = 1:block_size:padded_size(2)
            % 获取量化后的系数
            quantized_coeffs = blocks{(i-1)/block_size+1, (j-1)/block_size+1};
            
            % 反量化
            dct_coeffs = quantized_coeffs .* quant_table;
            
            % 逆DCT变换
            block_reconstructed = idct2(dct_coeffs);
            
            % 电平偏移恢复
            block_reconstructed = block_reconstructed + 128;
            
            % 放回重建图像
            channel_reconstructed(i:i+block_size-1, j:j+block_size-1) = block_reconstructed;
        end
    end
    
    % 裁剪到原始尺寸
    original_size = compressed_channel.original_size;
    channel_reconstructed = channel_reconstructed(1:original_size(1), 1:original_size(2));
end

6. 压缩性能分析

function compression_info = calculate_compression_info(original_img, compressed_data, info_Y, info_Cb, info_Cr)
    % 计算压缩性能指标
    
    % 原始数据量 (字节)
    original_size_bytes = numel(original_img) * 8; % 假设8位/像素
    
    % 估算压缩后数据量
    % 这里简化计算，实际应该包括量化表、霍夫曼表等开销
    compressed_size_estimate = 0;
    
    % Y分量
    num_blocks_Y = numel(compressed_data.Y.blocks);
    non_zero_coeffs_Y = info_Y.non_zero_ac + num_blocks_Y; % 包括DC系数
    
    % Cb分量
    num_blocks_Cb = numel(compressed_data.Cb.blocks);
    non_zero_coeffs_Cb = info_Cb.non_zero_ac + num_blocks_Cb;
    
    % Cr分量
    num_blocks_Cr = numel(compressed_data.Cr.blocks);
    non_zero_coeffs_Cr = info_Cr.non_zero_ac + num_blocks_Cr;
    
    % 假设每个非零系数平均需要12位 (包括游程和值)
    compressed_size_estimate = (non_zero_coeffs_Y + non_zero_coeffs_Cb + non_zero_coeffs_Cr) * 1.5;
    
    % 压缩比
    compression_ratio = original_size_bytes / compressed_size_estimate;
    
    % 保存压缩信息
    compression_info.original_size_bytes = original_size_bytes;
    compression_info.compressed_size_estimate = compressed_size_estimate;
    compression_info.compression_ratio = compression_ratio;
    compression_info.bit_rate = compressed_size_estimate * 8 / (numel(original_img) / 3); % bpp
    
    % 系数统计
    compression_info.non_zero_coeffs.total = non_zero_coeffs_Y + non_zero_coeffs_Cb + non_zero_coeffs_Cr;
    compression_info.non_zero_coeffs.Y = non_zero_coeffs_Y;
    compression_info.non_zero_coeffs.Cb = non_zero_coeffs_Cb;
    compression_info.non_zero_coeffs.Cr = non_zero_coeffs_Cr;
    
    fprintf('压缩性能分析:\n');
    fprintf('原始数据量: %.2f KB\n', original_size_bytes/1024);
    fprintf('压缩后估算: %.2f KB\n', compressed_size_estimate/1024);
    fprintf('压缩比: %.2f:1\n', compression_ratio);
    fprintf('比特率: %.2f bpp\n', compression_info.bit_rate);
end

可视化与分析工具

1. DCT系数分析

function analyze_dct_coefficients(compressed_data)
    % 分析DCT系数分布
    
    figure('Position', [100, 100, 1200, 800]);
    
    % Y分量系数分析
    subplot(2,3,1);
    dc_y = [];
    for i = 1:numel(compressed_data.Y.blocks)
        block = compressed_data.Y.blocks{i};
        dc_y = [dc_y; block(1,1)];
    end
    hist(dc_y, 50);
    title('Y分量DC系数分布');
    xlabel('DC系数值');
    ylabel('频数');
    
    subplot(2,3,2);
    ac_y = [];
    for i = 1:numel(compressed_data.Y.blocks)
        block = compressed_data.Y.blocks{i};
        ac_block = block(:);
        ac_block(1) = []; % 移除DC
        ac_y = [ac_y; ac_block];
    end
    hist(ac_y, 100);
    title('Y分量AC系数分布');
    xlabel('AC系数值');
    ylabel('频数');
    
    % 显示量化表
    subplot(2,3,3);
    imagesc(compressed_data.Y.quant_table);
    colorbar;
    title('Y分量量化表');
    axis image;
    
    % 显示一个示例块的DCT系数
    subplot(2,3,4);
    example_block = compressed_data.Y.blocks{1,1};
    imagesc(log(abs(example_block) + 1));
    colorbar;
    title('示例块DCT系数 (对数尺度)');
    axis image;
    
    % 零系数比例
    subplot(2,3,5);
    zero_ratio_y = sum(ac_y == 0) / length(ac_y);
    zero_ratio_cb = sum([compressed_data.Cb.blocks{:}](:) == 0) / numel([compressed_data.Cb.blocks{:}]);
    zero_ratio_cr = sum([compressed_data.Cr.blocks{:}](:) == 0) / numel([compressed_data.Cr.blocks{:}]);
    
    bar([zero_ratio_y, zero_ratio_cb, zero_ratio_cr]);
    set(gca, 'XTickLabel', {'Y', 'Cb', 'Cr'});
    title('零系数比例');
    ylabel('比例');
    
    % 系数能量分布
    subplot(2,3,6);
    energy_y = accum_dct_energy(compressed_data.Y.blocks);
    plot(1:64, energy_y, 'b-', 'LineWidth', 2);
    title('DCT系数累积能量分布');
    xlabel('Zig-zag扫描顺序');
    ylabel('累积能量比例');
    grid on;
end

function energy_accum = accum_dct_energy(blocks)
    % 计算DCT系数累积能量分布
    
    total_energy = 0;
    energy_by_position = zeros(64, 1);
    
    for i = 1:numel(blocks)
        block = blocks{i};
        zigzag_coeffs = zigzag_scan(block);
        energy = zigzag_coeffs.^2;
        total_energy = total_energy + sum(energy);
        energy_by_position = energy_by_position + energy;
    end
    
    % 按能量排序
    [~, sort_idx] = sort(energy_by_position, 'descend');
    
    % 计算累积能量
    energy_accum = zeros(64, 1);
    for i = 1:64
        energy_accum(i) = sum(energy_by_position(sort_idx(1:i))) / total_energy;
    end
end

function zigzag = zigzag_scan(block)
    % Zig-zag扫描
    zigzag = zeros(64, 1);
    index = 1;
    
    for sum = 0:2:14
        if rem(sum, 2) == 0
            for i = max(1, sum-7):min(8, sum+1)
                j = sum - i + 2;
                if j >= 1 && j <= 8
                    zigzag(index) = block(i, j);
                    index = index + 1;
                end
            end
        else
            for i = max(1, sum-7):min(8, sum+1)
                j = sum - i + 2;
                if j >= 1 && j <= 8
                    zigzag(index) = block(j, i);
                    index = index + 1;
                end
            end
        end
    end
end

2. 主演示函数

function dct_compression_demo()
    % DCT压缩演示函数
    
    % 读取测试图像
    try
        img = imread('peppers.png');
    catch
        % 如果找不到图像，创建一个测试图像
        fprintf('创建测试图像...\n');
        img = uint8(255 * rand(256, 256, 3));
    end
    
    % 测试不同质量因子
    quality_factors = [10, 25, 50, 75, 90];
    
    figure('Position', [50, 50, 1400, 900]);
    
    for i = 1:length(quality_factors)
        qf = quality_factors(i);
        
        % 压缩和解压缩
        [compressed_data, reconstructed_img, compression_info] = dct_color_compress(img, qf);
        
        % 计算PSNR
        psnr_val = psnr(im2double(reconstructed_img), im2double(img));
        
        % 显示结果
        subplot(2, 3, i);
        imshow(reconstructed_img);
        title(sprintf('质量因子: %d\nPSNR: %.2f dB, 压缩比: %.2f:1', ...
            qf, psnr_val, compression_info.compression_ratio));
        
        fprintf('质量因子 %d: PSNR = %.2f dB, 压缩比 = %.2f:1\n', ...
            qf, psnr_val, compression_info.compression_ratio);
    end
    
    % 显示原始图像
    subplot(2, 3, 6);
    imshow(img);
    title('原始图像');
    
    % 详细分析中等质量的情况
    fprintf('\n进行详细分析...\n');
    [compressed_data, reconstructed_img, compression_info] = dct_color_compress(img, 50);
    
    % 分析DCT系数
    figure('Position', [100, 100, 1200, 800]);
    analyze_dct_coefficients(compressed_data);
    
    % 显示压缩前后对比
    figure('Position', [150, 150, 1000, 400]);
    subplot(1,2,1);
    imshow(img);
    title('原始图像');
    
    subplot(1,2,2);
    imshow(reconstructed_img);
    title('压缩后图像 (质量因子50)');
    
    % 计算并显示差异
    diff_img = im2double(img) - reconstructed_img;
    diff_img = (diff_img - min(diff_img(:))) / (max(diff_img(:)) - min(diff_img(:)));
    
    figure('Position', [200, 200, 600, 400]);
    imshow(diff_img);
    title('压缩误差 (归一化显示)');
    colorbar;
end

高级特性扩展

1. 自适应量化

function adaptive_quant_table = get_adaptive_quantization(block, base_quant_table, activity)
    % 基于块活动性的自适应量化
    % activity: 块的活动性度量 (如方差)
    
    if activity > 0.1
        % 高活动性区域，使用较细的量化
        scale_factor = 0.7;
    else
        % 平坦区域，使用较粗的量化
        scale_factor = 1.3;
    end
    
    adaptive_quant_table = max(1, round(base_quant_table * scale_factor));
end

2. 渐进式压缩

function progressive_data = progressive_compress(img, quality_levels)
    % 渐进式DCT压缩
    % quality_levels: 质量级别数组，如 [10, 30, 50, 70, 90]
    
    progressive_data = cell(length(quality_levels), 1);
    
    for i = 1:length(quality_levels)
        fprintf('压缩级别 %d/%d (质量因子: %d)...\n', i, length(quality_levels), quality_levels(i));
        [compressed_data, ~, compression_info] = dct_color_compress(img, quality_levels(i));
        progressive_data{i} = compressed_data;
    end
end

参考代码 基于DCT的彩色图像压缩 www.3dddown.com/cna/82678.html

这个完整的DCT彩色图像压缩系统提供了从基本压缩到高级分析的完整功能。