【影像壓縮】基于matlab分形編碼影像壓縮【含Matlab原始碼 824期】

一、簡介

影像壓縮就是就是通過去除這些資料冗余來減少表示資料所需的位元數，去除多余資料，以數學的觀點來看，這一程序實際上就是將二維像素陣列變換為一個在統計上無關聯的資料集合，
影像壓縮是指以較少的位元有損或無損地表示原來的像素矩陣的技術，也稱影像編碼，
影像資料之所以能被壓縮，就是因為資料中存在著冗余，影像資料的冗余主要表現為：影像中相鄰像素間的相關性引起的空間冗余；影像序列中不同幀之間存在相關性引起的時間冗余；不同彩色平面或頻譜帶的相關性引起的頻譜冗余，

從壓縮編碼演算法原理上可以分類為：
(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

三、運行結果

四、備注

完整代碼或者代寫添加QQ 1564658423