图像压缩基于分形编码图像压缩matlab源码

Posted 2021-05-29 Matlab咨询QQ1575304183

一、简介

图像压缩就是就是通过去除这些数据冗余来减少表示数据所需的比特数，去除多余数据。以数学的观点来看，这一过程实际上就是将二维像素阵列变换为一个在统计上无关联的数据集合。
图像压缩是指以较少的比特有损或无损地表示原来的像素矩阵的技术，也称图像编码。
图像数据之所以能被压缩，就是因为数据中存在着冗余。图像数据的冗余主要表现为：图像中相邻像素间的相关性引起的空间冗余；图像序列中不同帧之间存在相关性引起的时间冗余；不同彩色平面或频谱带的相关性引起的频谱冗余。

从压缩编码算法原理上可以分类为：
(1)无损压缩编码种类：哈夫曼编码、算术编码、行程编码、Lempel zev编码
(2)有损压缩编码种类：
预测编码：DPCM，运动补偿
频率域方法：正交变换编码(如DCT)，子带编码
空间域方法：统计分块编码
模型方法：分形编码，模型基编码
基于重要性：滤波，子采样，比特分配，矢量量化

(3)混合编码 ?JBIG，H261，JPEG，MPEG等技术标准

二、源代码

clear;clc
tic
A = imread('Lena.tiff');%载入待编码图像
[a, b]=size(A);
blocksizeD = 16;%domain分块大小
maxblockno = 6;
blocknoD = a/blocksizeD;
blockDtotal = blocknoD^2;
% blockDtotal = 12;
D = cell(1,blockDtotal);
for t = 1:blockDtotal
    D(1, t) = {divideblock(A, blocksizeD, t, a)};
end
blocksizeR = blocksizeD/2;%range块边长设置为domain块边长的一半
blocknoR = a/blocksizeR;
blockRtotal = blocknoR^2;
R = cell(1,blockRtotal);
for t = 1:blockRtotal
    R(1, t) = {divideblock(A, blocksizeR, t, a)};
end
E = cell(1,blockDtotal);
for k = 1:blockDtotal
    for s = 1:blocksizeR
        for t = 1:blocksizeR
            E{k}(s, t) = (D{k}(2*s-1, 2*t-1)+D{k}(2*s-1, 2*t)+D{k}(2*s, 2*t-1)+D{k}(2*s, 2*t))/4;
        end
    end
end
blockUtotal = blockDtotal*8;
U = cell(1, blockUtotal);
for k = 1:blockUtotal
    switch mod(k, 8)
        case 1
            U{k} = E{(k+7)/8};
        case 2
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{(k+6)/8}(t, blocksizeR+1-s);
                end
            end
        case 3
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{(k+5)/8}(blocksizeR+1-t, s);
                end
            end
        case 4
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{(k+4)/8}(blocksizeR+1-t, blocksizeR+1-s);
                end
            end
        case 5
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{(k+3)/8}(blocksizeR-s+1, blocksizeR-t+1);
                end
            end
        case 6
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{(k+2)/8}(s, blocksizeR-t+1);
                end
            end
        case 7
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{(k+1)/8}(blocksizeR-s+1, t);
                end
            end
        case 0
            for t = 1:blocksizeR
                for s = 1:blocksizeR
                    U{k}(t, s) = E{k/8}(s, t);
                end
            end
    end
end

blockRtotalsize = blocksizeR^2;

%以下计算不同块间度量，块间度量定义为单位化某块后另一块的八种变化与其的最小度量距离
markMatrix = zeros(blockDtotal);
for k = 1:(blockDtotal-1)
    for t = (k+1):blockDtotal
        tempM = U{8*k-7}/(sum(U{8*k-7}(:))+eps);
        for s = 1:8
            tempMatrix = U{8*t-8+s};
            tempMmultiMatrix = tempM.*tempMatrix;
            tempMatrixSqure = tempMatrix.*tempMatrix;
            sMark = (blockRtotalsize*(sum(tempMmultiMatrix(:)))-sum(tempMatrix(:)))/(blockRtotalsize*sum(tempMatrixSqure(:))-(sum(tempMatrix(:)))^2+eps);
            oMark = (1-sMark*sum(tempMatrix(:)))/blockRtotalsize;
            tempHMatrix = sMark*tempMatrix+oMark-tempM;
            tempHM = tempHMatrix.*tempHMatrix;
            hMark(s) = sum(tempHM(:));
        end
        markMatrix(k,t) = min(min(min(hMark(1),hMark(2)),min(hMark(3),hMark(4))),min(min(hMark(5),hMark(6)),min(hMark(7),hMark(8))));
    end
end

[maxM,maxmark] = max(markMatrix(:));
maxmarkx = 1+mod(maxmark-1,blockDtotal);
maxmarky = 1+fix((maxmark-1)/blockDtotal);
maxIndexMatrix(1) = maxmarkx;
maxIndexMatrix(2) = maxmarky;
for k = 3:maxblockno
    maxIndexMatrix(k) = 0;
end

for p = 3:maxblockno
    for k = 1:blockDtotal
        judgeno = 0;
        for s = 1:(p-1)
            if(maxIndexMatrix(s)==k)
                judgeno = judgeno | 1;
            end
        end
        tempMaxMatrix(k) = 0;
        if (judgeno == 0)
            for t = 1:maxblockno
                for s = 1:(p-1)
                    tempMaxMatrix(k) = tempMaxMatrix(k)+markMatrix(s,k)+markMatrix(k,s);
                end
            end
        end
    end
    [C,I] = max(tempMaxMatrix);
    maxIndexMatrix(p) = I;
end
%8个索引k对应的U中的序号为8*k-7
blockWtotal = maxblockno*8;
W = cell(1, blockWtotal);
for k = 1:maxblockno
    W{8*k-7} = U{maxIndexMatrix(k)*8-7};
    W{8*k-6} = U{maxIndexMatrix(k)*8-6};
    W{8*k-5} = U{maxIndexMatrix(k)*8-5};
    W{8*k-4} = U{maxIndexMatrix(k)*8-4};
    W{8*k-3} = U{maxIndexMatrix(k)*8-3};
    W{8*k-2} = U{maxIndexMatrix(k)*8-2};
    W{8*k-1} = U{maxIndexMatrix(k)*8-1};
    W{8*k} = U{maxIndexMatrix(k)*8};
end

%R{i}为小块，D{i}为大块，U{i}为D经过变换后的小块
%对任一个R{i}，与每个U{i}进行匹配计算，选取获得最小H的压缩比和常数项，定义Radd为其中第i个元素为对应R{i}的所有元素之和
%Radd = zeros(1,blockRtotal);
%for k = 1:blockRtotal
%    for t = 1:blocksizeR
%        for s = 1:blocksizeR
%            T(1,k) = T(1,k) + R{k}(t,s);
%        end
%    end
%end
Radd = addofall(blockRtotal,R);
%Uadd为U{i}的所有元素之和构成的矩阵,P为U{i}所有元素平方和构成的矩阵
Wadd = addofall(blockWtotal,W);
Wsquareadd = multiofall(blockWtotal,W,W);
Rsquareadd = multiofall(blockRtotal,R,R);
%RmultiUadd为一个矩阵，其中的第i,j个元素对应R{i}与U{j}的乘积的和
RmultiWadd = zeros(blockRtotal,blockWtotal);
for t = 1:blockRtotal
    for s = 1:blockWtotal
        Matrixtemp = R{t}.*W{s};
        RmultiWadd(t,s) = sum(Matrixtemp(:));
    end
end
WmultiWadd = zeros(blockWtotal,blockWtotal);
for s = 1:blockWtotal
    for t = s:blockWtotal
        Matrixtemp = W{t}.*W{s};
        WmultiWadd(s,t) = sum(Matrixtemp(:));
%        UmultiUadd(t,s) = UmultiUadd(s,t);
    end
end

三、运行结果

完整代码或仿真咨询QQ1575304183