Difference between revisions of "JPEG компресия"

Revision as of 12:10, 30 June 2014

Компресия със загуби:

1. Цветното изображението най-често се представя като съвкупност от три двумерни матрици. Всеки елемент от матрица е представя стойността за съответния пиксел на един от условно наречените основни цветове Червено, Зелено, Синьо. Тези цветове отговарят на цветовата спектралната чувствителност на окото.

Цветно (RGB) изображени представено като три измерна матрица

Визуализиране на отделните матрици (цветни канали)

Понякога при обработка на избраженията е по удачно преминаването към други цветови схеми (YCrCb, HSV). В случая с YCrCb:

Y - яркостен сигнал
Cr, Cb - хроматични съставящи

 ${\begin{aligned}Y'&=&16&+&(65.481\cdot R'&+&128.553\cdot G'&+&24.966\cdot B')\\C_{B}&=&128&+&(-37.797\cdot R'&-&74.203\cdot G'&+&112.0\cdot B')\\C_{R}&=&128&+&(112.0\cdot R'&-&93.786\cdot G'&-&18.214\cdot B')\end{aligned}}$

A = double(rgb2ycbcr(A)); %

Поради това, че окото е по чувствително към измененията на яркоста, хроматичните съставящи могат да се предават с по-малка резолюция и да се компресират с по-ниско качество без това да окаже видимо влияние върху качеството на избражението.(4:2:2,4:2:0...)

case '4:2:0'
       Cb = imresize(A(:,:,2),[Height/2 Width/2],'cubic');
       Cr = imresize(A(:,:,3),[Height/2 Width/2],'cubic');
case '4:2:2'
       Cb = imresize(A(:,:,2),[Height Width/2],'cubic');
       Cr = imresize(A(:,:,3),[Height Width/2],'cubic');

2. Пеминаване в честотна област

Изображението представено в пространствената област е със сравнително висока ентропия. (няма връзка между стойностите на съседни пиксели).

Коефициентите на DCT трансформацията са с висока сходимост, а повечето от тях са пренебрежимо малки. В замисимост от желаното качество се избират тегловни (квантуващи) коефициенти, на които се разделят DCT коефициентите. След преминаване в честотната могат да се премахва информацията за честотите, които ще останат незабелижими за човешкото възприятие.

Високите честоти отговарят на резки преходи, контрастни линии, дребни обекти.
Честоти с малки коефициенти (малка аплитуда) остават незабелижими в цялостното изображение

При преминаване в честотната област изображението се разделя на блокове от 8х8 пиксела. Всеки блок се обработва поотделно.

for m = 1:N:Height

   for n = 1:N:Width
           
       t = y(m:m+N-1,n:n+N-1) - 128;
       Yy = dct2(t); % N x N 2D DCT of input image
           % quantize the DCT coefficients

За да премахмен информацията, която не ни интересува се използва квантуване. Всеки DCT коефициент се разделя на точно определен коефициент от квантуващата матрица. Квантуващата таблица е генирирана в съотвествие качеството, което желаем. При квантуването коефициентите с малка стойност се нулират нула.

jpgQstepsY = [ 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];

temp = floor(Yy./(Qscale*QstepsY) + 0.5);

3. Зиг-заг пренареждане

DC - компонентата се кодира диференциално (за всички блокове от по 8х8), като разликите са Hoffman кодиране
AC - коефициенти се пренареждат зигзаг образно, след това се приланга диференциално и статистическо, а накрая кодиране на хофман (или аритметична компресия), на ненулевите коефициенти, и статестическа компресия на нулевите.

референтно изображение

Лена

Matlab реализация

% JPEG Основен тип компресия

clear;
A = imread('birds.ras'); 
[Height,Width,Depth] = size(A);
N = 8; % Размер на подматрици

% Размерите трябва да са кратни на 8
if mod(Height,N) ~= 0
    Height = floor(Height/N)*N;
end
if mod(Width,N) ~= 0
    Width = floor(Width/N)*N;
end

A1 = A(1:Height,1:Width,:);
clear A
A = A1;

SamplingFormat = '4:2:0';
if Depth == 1 % черно бяло
    y = double(A);
else
    A = double(rgb2ycbcr(A)); % Цветно
    y = A(:,:,1);
    switch SamplingFormat
    case '4:2:0'
        Cb = imresize(A(:,:,2),[Height/2 Width/2],'cubic');
        Cr = imresize(A(:,:,3),[Height/2 Width/2],'cubic');
    case '4:2:2'
        Cb = imresize(A(:,:,2),[Height Width/2],'cubic');
        Cr = imresize(A(:,:,3),[Height Width/2],'cubic');
    end
end

jpgQstepsY = [  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];

QstepsY = jpgQstepsY;
Qscale = 1.5;
Yy = zeros(N,N);
xqY = zeros(Height,Width);
acBitsY = 0;
dcBitsY = 0;

if Depth > 1
    jpgQstepsC = [  17 18 24 47 66 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];
    QstepsC = jpgQstepsC;
    YCb = zeros(N,N);
    YCr = zeros(N,N);
    switch SamplingFormat
    case '4:2:0'
        xqCb = zeros(Height/2,Width/2);
        xqCr = zeros(Height/2,Width/2);
    case '4:2:2'
        xqCb = zeros(Height,Width/2);
        xqCr = zeros(Height,Width/2);
    end
    
    acBitsCb = 0;
    dcBitsCb = 0;
    acBitsCr = 0;
    dcBitsCr = 0;
end

% Изчисляване на y  компонента
for m = 1:N:Height
    for n = 1:N:Width
            
        t = y(m:m+N-1,n:n+N-1) - 128;
        Yy = dct2(t); % N x N 2D DCT of input image
            % quantize the DCT coefficients
        temp = floor(Yy./(Qscale*QstepsY) + 0.5);

        % Calculate bits for the DC difference    
        if n == 1
            DC = temp(1,1);
            dcBitsY = dcBitsY + jpgDCbits(DC,'Y');
        else
                DC = temp(1,1) - DC;
                dcBitsY = dcBitsY + jpgDCbits(DC,’Y’);
                DC = temp(1,1);
        end
            

            % Calculate the bits for the AC coefficients
        ACblkBits = jpgACbits(temp,'Y');
        acBitsY = acBitsY + ACblkBits;
        % dequantize & IDCT the DCT coefficients
        xqY(m:m+N-1,n:n+N-1)= idct2(temp .* (Qscale*QstepsY))+ 128;
    end
end
% If the input image is a color image,
% calculate the bits for the chroma components
if Depth > 1
    if strcmpi(SamplingFormat,'4:2:0')
        EndRow = Height/2;
    else
    EndRow = Height;
    end
    
    for m = 1:N:EndRow
        for n = 1:N:Width/2
            t1 = Cb(m:m+N-1,n:n+N-1) - 128;
            t2 = Cr(m:m+N-1,n:n+N-1) - 128;
            Ycb = dct2(t1); % N x N 2D DCT of Cb image
            Ycr = dct2(t2);
            temp1 = floor(Ycb./(Qscale*QstepsC) + 0.5);
            temp2 = floor(Ycr./(Qscale*QstepsC) + 0.5);
            if n == 1
                DC1 = temp1(1,1);
                DC2 = temp2(1,1);
                dcBitsCb = dcBitsCb + jpgDCbits(DC1,'C');
                dcBitsCr = dcBitsCr + jpgDCbits(DC2,'C');
            else
                DC1 = temp1(1,1) - DC1;
                DC2 = temp2(1,1) - DC2;
                dcBitsCb = dcBitsCb + jpgDCbits(DC1,'C');
                dcBitsCr = dcBitsCr + jpgDCbits(DC2,'C');
                DC1 = temp1(1,1);
                DC2 = temp2(1,1);
            end
            ACblkBits1 = jpgACbits(temp1,'C');
            ACblkBits2 = jpgACbits(temp2,'C');
            acBitsCb = acBitsCb + ACblkBits1;
            acBitsCr = acBitsCr + ACblkBits2;

            % dequantize and IDCT the coefficients
            xqCb(m:m+N-1,n:n+N-1)= idct2(temp1 .* (Qscale*QstepsC))+ 128;
            xqCr(m:m+N-1,n:n+N-1)= idct2(temp2 .* (Qscale*QstepsC))+ 128;
        end
    end
end
%
mse = std2(y-xqY);
snr = 20*log10(std2(y)/mse);
sprintf('SNR = %4.2f\n',snr)

if Depth == 1
    TotalBits = acBitsY + dcBitsY;
    figure,imshow(xqY,[])
    title(['JPG compressed ' '@ ' num2str(TotalBits/(Height*Width)) ' bpp'])
else
    
TotalBits = acBitsY + dcBitsY + dcBitsCb +...
acBitsCb + dcBitsCr + acBitsCr;
c1 = imresize(xqCb,[Height Width],'cubic');
c2 = imresize(xqCr,[Height Width],'cubic');
mseb = std2(A(:,:,2)-c1);
snrb = 20*log10(std2(A(:,:,2))/mseb);
msec = std2(A(:,:,3)-c2);
snrc = 20*log10(std2(A(:,:,3))/msec);
sprintf('SNR(Cb) = %4.2fdB\tSNR(Cr) = %4.2fdB\n',snrb,snrc)
xq(:,:,1) = xqY;
xq(:,:,2) = c1;
xq(:,:,3) = c2;
figure,imshow(ycbcr2rgb(uint8(round(xq))))
title(['JPG compressed ' '@ ' num2str(TotalBits/(Height*Width)) ' bpp'])
end
sprintf('Bit rate = %4.2f bpp\n',TotalBits/(Height*Width))

A


function y = ZigZag(x)
% y = ZigZag(x)
% returns the zigzag scan addresses of an 8 x 8 array
y(1) = x(1,1);y(2) = x(1,2);y(3) = x(2,1);y(4) = x(3,1);
y(5) = x(2,2);y(6) = x(1,3);y(7) = x(1,4);y(8) = x(2,3);
y(9) = x(3,2);y(10) = x(4,1);y(11) = x(5,1);y(12) = x(4,2);
y(13) = x(3,3);y(14) = x(2,4);y(15) = x(1,5);y(16) = x(1,6);
y(17) = x(2,5);y(18) = x(3,4);y(19) = x(4,3);y(20) = x(5,2);
y(21) = x(6,1);y(22) = x(7,1);y(23) = x(6,2);y(24) = x(5,3);
y(25) = x(4,4);y(26) = x(3,5);y(27) = x(2,6);y(28) = x(1,7);
y(29) = x(1,8);y(30) = x(2,7);y(31) = x(3,6);y(32) = x(4,5);
y(33) = x(5,4);y(34) = x(6,3);y(35) = x(7,2);y(36) = x(8,1);
y(37) = x(8,2);y(38) = x(7,3);y(39) = x(6,4);y(40) = x(5,5);
y(41) = x(4,6);y(42) = x(3,7);y(43) = x(2,8);y(44) = x(3,8);
y(45) = x(4,7);y(46) = x(5,6);y(47) = x(6,5);y(48) = x(7,4);
y(49) = x(8,3);y(50) = x(8,4);y(51) = x(7,5);y(52) = x(6,6);
y(53) = x(5,7);y(54) = x(4,8);y(55) = x(5,8);y(56) = x(6,7);
y(57) = x(7,6);y(58) = x(8,5);y(59) = x(8,6);y(60) = x(7,7);
y(61) = x(6,8);y(62) = x(7,8);y(63) = x(8,7);y(64) = x(8,8);

B


function Bits = jpgDCbits(dc,C)
% Bits = jpgDCbits(dc,C)
% Computes the exact number of bits to entropy code
% the DC differential coefficient using the JPEG default
% Huffman tables.
% Input:
% dc = predicted DC differential value
% C = ’Y’ for luma or ’C’ for chroma
% Output:
% # of bits for the DC differential
%
% CodeLengthY contains the lengths of codes for the
% luma DC differential for size categories fro 0 to 11 inclusive.
% Code length equals the length of the Huffman code for
% the size category plus the binary code for the value.
% CodeLengthC is the corresponding table for the chroma.

CodeLengthY = int16([3 4 5 5 7 8 10 12 14 16 18 20]);
CodeLengthC = int16([2 3 5 6 7 9 11 13 15 18 19 20]);
switch C
case 'Y'
    CodeLength = CodeLengthY;
case 'C'
    CodeLength = CodeLengthC;
end
if dc == 0
    Bits = double(CodeLength(1));
else
    Bits = double(CodeLength(round(log2(abs(dc))+0.5)+1));
end

C

function Bits = jpgACbits(x,C)
% y = jpgACbits(x,C)
% Computes the exact number of bits to Huffman code
% an N x N quantized DCT block
% using JPEG Huffman table for AC coefficients.
% DCT coefficients are zigzag scanned and
% the RL/Amplitude pairs are found. Then,
% using the JPEG RL/Amp table, number of bits to
% compress the block is determined.
% Input:
% x = N x N quantized DCT block
% C = ’Y’ for luma or ’C’ for chroma
% Outputs:
% Bits = total bits to encode the block
%
% RLCy is the JPEG run-length/size-category Huffman table
% for luma AC coefficients &
% RLCc is the corresponding Huffman table for the chroma component
% Each entry in the table is the code length for the RL/size
% category which equals the RL/size Huffman code plus the
% binary code of the actual non-zero coefficient amplitude.
% RLC{1}(1) is the EOB code corresponding to 0/0
% RLC{16}(1) is the Huffman code for RL/size = 16/0
%

RLCy = cell(1,15);
RLCy{1}=int16([4 3 4 6 8 10 12 14 18 25 26]);
RLCy{2} = int16([5 8 10 13 16 22 23 24 25 26]);
RLCy{3} = int16([6 10 13 20 21 22 23 24 25 26]);
RLCy{4} = int16([7 11 14 20 21 22 23 24 25 26]);
RLCy{5} = int16([7 12 19 20 21 22 23 24 25 26]);
RLCy{6} = int16([8 12 19 20 21 22 23 24 25 26]);
RLCy{7} = int16([8 13 19 20 21 22 23 24 25 26]);
RLCy{8} = int16([9 13 19 20 21 22 23 24 25 26]);
RLCy{9} = int16([9 17 19 20 21 22 23 24 25 26]);
RLCy{10} = int16([10 18 19 20 21 22 23 24 25 26]);
RLCy{11} = int16([10 18 19 20 21 22 23 24 25 26]);
RLCy{12} = int16([10 18 19 20 21 22 23 24 25 26]);
RLCy{13} = int16([11 18 19 20 21 22 23 24 25 26]);
RLCy{14} = int16([12 18 19 20 21 22 23 24 25 26]);
RLCy{15} = int16([13 18 19 20 21 22 23 24 25 26]);
RLCy{16} = int16([12 17 18 19 20 21 22 23 24 25 26]);
%
RLCc = cell(1,15);
RLCc{1}=int16([2 3 5 7 9 12 15 23 24 25 26]);
RLCc{2} = int16([5 8 11 14 20 22 23 24 25 26]);
RLCc{3} = int16([6 9 13 15 21 22 23 24 25 26]);
RLCc{4} = int16([6 9 12 15 21 22 23 24 25 26]);
RLCc{5} = int16([6 10 14 20 21 22 23 24 25 26]);
RLCc{6} = int16([7 11 17 20 21 22 23 24 25 26]);
RLCc{7} = int16([8 12 19 20 21 22 23 24 25 26]);
RLCc{8} = int16([7 12 19 20 21 22 23 24 25 26]);
RLCc{9} = int16([8 14 19 20 21 22 23 24 25 26]);
RLCc{10} = int16([9 14 19 20 21 22 23 24 25 26]);
RLCc{11} = int16([9 14 19 20 21 22 23 24 25 26]);
RLCc{12} = int16([9 14 19 20 21 22 23 24 25 26]);
RLCc{13} = int16([9 18 19 20 21 22 23 24 25 26]);
RLCc{14} = int16([11 18 19 20 21 22 23 24 25 26]);
RLCc{15} = int16([13 18 19 20 21 22 23 24 25 26]);
RLCc{16} = int16([11 17 18 19 20 21 22 23 24 25 26]);
%
switch C
case 'Y'
    RLC = RLCy;
case 'C'
    RLC = RLCc;
end
x1 = ZigZag(x); % zigzag scan the 8 x 8 DCT blocks
%
k = 2; Count = 0; Bits = double(0);
while k <= 64
    if x1(k) == 0
        Count = Count + 1;
        if k == 64
            Bits = Bits + double(RLC{1}(1));
            break;
        end
    else
        if Count == 0
            RL = Count;
            Level = round(log2(abs(x1(k)))+0.5);
            Bits = Bits + double(RLC{RL+1}(Level+1));
        elseif Count >= 1 && Count <=15
            RL = Count;
            Level = round(log2(abs(x1(k)))+0.5);
            Bits = Bits + double(RLC{RL+1}(Level));
            Count = 0;
        else
            Bits = Bits + double(RLC{16}(1));
            Count = Count - 16;
        end
    end
    k = k + 1;
end

Литература

K. S. Thyagarajan, Still Image and Video Compression with MATLAB

DCT - discrete cosine transformation - дискретна косиносова трансформация

Anonymous

Search

Difference between revisions of "JPEG компресия"

Namespaces

More

Page actions

Revision as of 12:10, 30 June 2014

Contents

референтно изображение

Matlab реализация

A

B

C

Литература

Navigation

Wiki tools

Wiki tools

Anonymous

Search

Difference between revisions of "JPEG компресия"

Revision as of 12:10, 30 June 2014

Contents

референтно изображение

Matlab реализация

A

B

C

Литература

Navigation

Wiki tools

Page tools

Categories