SVD、DCT图像压缩实践
文章目录
- SVD图像压缩实践
- MATLAB实现
- 效果分析
- DCT图像压缩实践
- MATLAB实现DCT图像压缩
- C++实现
SVD图像压缩实践
SVD是一种矩阵分解技术,可以将一个矩阵分解为三个矩阵的乘积,即 U、Σ 和 V。对于图像压缩,我们将图像矩阵分解为 U、Σ 和 V 三个矩阵,其中 Σ 矩阵包含了图像的主要信息,而 U 和 V 矩阵则包含了图像的细节信息。通过对 Σ 矩阵进行压缩,我们可以有效地减少图像文件的大小。
基于 SVD 的图像压缩算法主要分为以下几个步骤:
-
将图像转换为灰度图像。
-
将灰度图像转换为矩阵 A。
-
对矩阵 A 进行 SVD 分解,得到 U、Σ 和 V 三个矩阵。
-
对 Σ 矩阵进行压缩,保留前 k 个奇异值。
-
根据压缩后的 Σ 矩阵重建图像。
MATLAB实现
clc;
clear;
close all;filepath = 'cat.png';
% 1. 输入处理
if ~exist(filepath,'file')error('输入文件不存在');
end
img = imread(filepath);
if size(img, 3) == 3img = rgb2gray(img);
end
img_matrix = im2double(img);% 2. SVD分解
[U, S, V] = svd(img_matrix, 'econ');% 3. 压缩重构
k = 50;
img_compressed = U(:, 1:k) * S(1:k, 1:k) * V(:, 1:k)';original_size = numel(img_matrix);
compressed_size = k*(size(U,1)+size(V,1)+1);
fprintf('压缩率:%.2f%%\n',100*compressed_size/original_size);% 4. 输出显示
img_compressed = im2uint8(img_compressed);figure
subplot(1,2,1), imshow(img), title('原图');
subplot(1,2,2), imshow(img_compressed), title(['压缩图(k=',num2str(k),')']);
不同k值下的图像效果:
压缩率:25.70%
压缩率:51.40%
压缩率:77.10%
效果分析
基于 SVD 的图像压缩算法的压缩效果主要取决于压缩 Σ 矩阵的程度。奇异值越少,压缩比越高,但图像质量也会下降。因此,在实际应用中,需要根据具体需求权衡压缩比和图像质量之间的关系。
上面的例子展示了基于 SVD 的图像压缩算法在不同压缩比下的压缩效果。可以看出,随着压缩比的增加,图像质量逐渐下降。
DCT图像压缩实践
DCT将图像从空间域转换到频率域后,得到的DCT系数具有特定的物理意义。在频率域中,低频系数代表图像的全局特征,如主要颜色和整体亮度变化,而高频系数则对应图像的细节部分,如边缘和纹理信息。
通过观察DCT系数矩阵,可以发现大部分图像的能量都集中在左上角的低频区域。这是因为自然图像通常具有较多的低频成分,而高频成分则相对较少。因此,若想要减少图像数据量,重点应放在对高频系数的处理上。
在图像压缩中离散余弦变换(DCT)主要应用在JPEG压缩流程中,如下图所示,在此仅展示DCT变换这一步。使用DCT函数将低频的点都集中在左上角,再对右下角的高频数据丢弃。
MATLAB实现DCT图像压缩
clc;
clear;
close all;% 读取图像(使用MATLAB自带图像)
src = imread('cat.png');
if size(src,3)==3src = rgb2gray(src);
end
figure; imshow(src); title('原始图像');
set(gcf, 'Color', 'w'); % 设置figure背景为白色
set(gca, 'Color', 'none');figure, imshow(log(dct2(src)),[]), title('dct image')
set(gcf, 'Color', 'w'); % 设置figure背景为白色
colormap jet
colorbar% 转换为浮点型并归一化
src = im2double(src);
[h,w] = size(src);
big = 8;% 初始化DCT系数矩阵
srcDCT = zeros(size(src));% 分块DCT变换
for i1 = 1:big:hfor j1 = 1:big:wblock = src(i1:min(i1+big-1,h), j1:min(j1+big-1,w));srcDCT(i1:min(i1+big-1,h), j1:min(j1+big-1,w)) = dct2(block);end
end% 定义掩模矩阵(保留左上角低频系数)
mask = [1 1 1 1 1 0 0 0;1 1 1 1 0 0 0 0;1 1 1 0 0 0 0 0;1 1 0 0 0 0 0 0;1 0 0 0 0 0 0 0;0 0 0 0 0 0 0 0;0 0 0 0 0 0 0 0;0 0 0 0 0 0 0 0];% 应用掩模舍弃高频系数
for i1 = 1:big:hfor j1 = 1:big:wblock = srcDCT(i1:min(i1+big-1,h), j1:min(j1+big-1,w));block = block .* mask(1:size(block,1),1:size(block,2));srcDCT(i1:min(i1+big-1,h), j1:min(j1+big-1,w)) = block;end
end% 分块IDCT重建图像
iDct1 = zeros(size(srcDCT));
for i1 = 1:big:hfor j1 = 1:big:wblock = srcDCT(i1:min(i1+big-1,h), j1:min(j1+big-1,w));iDct1(i1:min(i1+big-1,h), j1:min(j1+big-1,w)) = idct2(block);end
end% 显示结果
figure; imshow(srcDCT,[]); title('DCT系数矩阵');
set(gcf, 'Color', 'w'); % 设置figure背景为白色
figure; imshow(iDct1); title('重建图像');
set(gcf, 'Color', 'w'); % 设置figure背景为白色% 计算PSNR
mse = mean((src(:)-iDct1(:)).^2);
psnr = 10*log10(1/mse);
fprintf('PSNR: %.2f dB\n', psnr);
C++实现
int main()
{Mat src = imread("coins.png",0);imshow("clock", src);src = Mat_<float>(src);int w = src.rows;int h = src.cols;int big;big = 8;src.convertTo(src, CV_32F, 1.0 / 255);Mat srcDCT = Mat_<float>(src);for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {dct(src(Rect(i1 * big, j1 * big, big, big)),srcDCT(Rect(i1 * big, j1 * big, big, big)));}}int mask[8][8] = {{1, 1, 1, 1, 0, 0, 0, 0},{1, 1, 1, 0, 0, 0, 0, 0},{1, 1, 0, 0, 0, 0, 0, 0},{1, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},};double s; for (int i1 = 0; i1 < (srcDCT.rows / big); i1++) {for (int j1 = 0; j1 < (srcDCT.cols / big); j1++) {for (int i = 0; i < big; i++){for (int j = 0; j < big; j++){if (mask[i % big][j % big] == 0) {srcDCT.at<int>(i1 * big + i, j1 * big + j) = 0;;}else {s = 0;}}}}}Mat iDct1 = srcDCT;for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {idct(srcDCT(Rect(i1 * big, j1 * big, big, big)), iDct1(Rect(i1 * big, j1 * big, big, big)));}}imshow("c", srcDCT);imshow("DstImage",iDct1);waitKey(0);return 0;
}
#include <opencv2/opencv.hpp>
#include <iostream>
using namespace cv;
using namespace std;double T = 100;int main()
{Mat src = imread("temp.jpg");//先读入图像,在这里读入的是RGB图像,有三个通道imshow("clock", src);//先显示图像int big;big = 8;//定义一下掩盖矩阵的大小vector<Mat> mv;split(src, mv);//通道分割Mat b = Mat_<float>(mv[0]);Mat g = Mat_<float>(mv[1]);Mat r = Mat_<float>(mv[2]);b.convertTo(b, CV_32F, 1.0 / 255);g.convertTo(g, CV_32F, 1.0 / 255);r.convertTo(r, CV_32F, 1.0 / 255);//要进行DCT变换先转换成float,我也不知道为什么Mat bDCT = Mat_<float>(b); //这里先复制矩阵,后面有用Mat gDCT = Mat_<float>(g);Mat rDCT = Mat_<float>(r);int w = b.rows;int h = b.cols;//三个通道的大小都是一样的,所以只要取一个就好了//DCT变换每个通道都做一遍for (int i1 = 0; i1 < (h / big); i1++) {//注意这里的取值和我们平常遇到的x,y不一样 一张宽度x像素、高度y像素的灰度图保存在一个y * x的矩阵中。for (int j1 = 0; j1 < (w / big); j1++) {dct(b(Rect(i1 * big, j1 * big, big, big)), bDCT(Rect(i1 * big, j1 * big, big, big)));//Rect是滑块函数,从指定位置取指定大小的矩阵,这就是跟Matlab分块处理一样}}//发现有的图像不是长宽相同的,所以取整除,而没有取到数据,跟之前相同就行了,复制的用处就来啦for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {dct(g(Rect(i1 * big, j1 * big, big, big)), gDCT(Rect(i1 * big, j1 * big, big, big)));}}for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {dct(r(Rect(i1 * big, j1 * big, big, big)), rDCT(Rect(i1 * big, j1 * big, big, big)));}}int mask[8][8] = {//定义掩盖矩阵{1, 1, 1, 1, 0, 0, 0, 0},{1, 1, 1, 0, 0, 0, 0, 0},{1, 1, 0, 0, 0, 0, 0, 0},{1, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},{0, 0, 0, 0, 0, 0, 0, 0},};//根据mask 矩阵 去掉右下角的数据,同样对每个通道做一次for (int i1 = 0; i1 < (bDCT.rows / big); i1++) {for (int j1 = 0; j1 < (bDCT.cols / big); j1++) {for (int i = 0; i < big; i++){for (int j = 0; j < big; j++){if (mask[i % big][j % big] == 0){bDCT.at<int>(i1 * big + i, j1 * big + j) = 0;//这里和mask矩阵相乘是同一个意思,同时扫描mask和DCT变换后的矩阵,当mask矩阵扫描到的值为0时,把DCT矩阵也置零,去掉高频数据}}}}}for (int i1 = 0; i1 < (gDCT.rows / big); i1++){for (int j1 = 0; j1 < (gDCT.cols / big); j1++){for (int i = 0; i < big; i++){for (int j = 0; j < big; j++){if (mask[i % big][j % big] == 0){gDCT.at<int>(i1 * big + i, j1 * big + j) = 0;;}}}}} for (int i1 = 0; i1 < (rDCT.rows / big); i1++){for (int j1 = 0; j1 < (rDCT.cols / big); j1++){for (int i = 0; i < big; i++){for (int j = 0; j < big; j++){if (mask[i % big][j % big] == 0){rDCT.at<int>(i1 * big + i, j1 * big + j) = 0;;}}}}}Mat iDctb = bDCT;Mat iDctg = gDCT;Mat iDctr = rDCT;//进行逆DCT变换,同做三次for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {idct(bDCT(Rect(i1 * big, j1 * big, big, big)), iDctb(Rect(i1 * big, j1 * big, big, big)));}}for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {idct(gDCT(Rect(i1 * big, j1 * big, big, big)), iDctg(Rect(i1 * big, j1 * big, big, big)));}}for (int i1 = 0; i1 < (h / big); i1++) {for (int j1 = 0; j1 < (w / big); j1++) {idct(rDCT(Rect(i1 * big, j1 * big, big, big)), iDctr(Rect(i1 * big, j1 * big, big, big)));}}mv[0] = iDctb;mv[1] = iDctg;mv[2] = iDctr;Mat dst;merge(mv, dst);//通道合并imshow("DstImage",dst);//显示图像waitKey(0);return 0;
}