Difference between revisions of "JPEG компресия"
From Ilianko
| Line 50: | Line 50: | ||
*Честоти с малки коефициенти (малка аплитуда) остават незабелижими в цялостното изображение | *Честоти с малки коефициенти (малка аплитуда) остават незабелижими в цялостното изображение | ||
| − | + | <!--За да премахмен информацията, която не ни интересува се използва квантуване. Всеки DCT коефициент се разделя на точно определен коефициент от квантуващата матрица. -->Квантуващата матрица е генирирана в съотвествие качеството, което желаем. При квантуването коефициентите с малка стойност се нулират нула. | |
| − | + | Квантуваща матрица за яркостния сигнал при качество 50% | |
| − | |||
| − | |||
| − | За да премахмен информацията, която не ни интересува се използва квантуване. Всеки DCT коефициент се разделя на точно определен коефициент от квантуващата матрица. Квантуващата | ||
| − | |||
jpgQstepsY = [ 16 11 10 16 24 40 51 61; | jpgQstepsY = [ 16 11 10 16 24 40 51 61; | ||
12 12 14 19 26 58 60 55; | 12 12 14 19 26 58 60 55; | ||
| Line 66: | Line 62: | ||
temp = floor(Yy./(Qscale*QstepsY) + 0.5); | temp = floor(Yy./(Qscale*QstepsY) + 0.5); | ||
| − | |||
Revision as of 12:41, 30 June 2014
Компресия със загуби:
Contents
Теория
Пространствена област
Цветното изображението най-често се представя като съвкупност от три двумерни матрици. Всеки елемент от матрица е стойността за съответния пиксел на един от условно наречените основни цветове Червено, Зелено, Синьо. Тези цветове отговарят на спектралната чувствителност на окото.
Цветно (RGB) изображени представено като три измерна матрица
Визуализиране на отделните матрици (цветни канали)
Понякога при обработка на избраженията е по-удачно използването на други цветови схеми (YCrCb, HSV). В случая при компресиране се използва най-често YCrCb:
- Y - яркостен сигнал
- Cr, Cb - хроматични съставящи
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');
Преминаване в честотна област
Изображението представено в пространствената област е със сравнително висока ентропия. (няма връзка между стойностите на съседни пиксели). С помощта на двумерна DCT трансформация изображението може да се изследва в честотната област.
При преминаване в честотната област изображението се разделя на блокове от 8х8 пиксела. Всеки блок се обработва поотделно.
N=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);
...
Коефициентите на DCT трансформацията са с висока сходимост, а повечето от тях са пренебрежимо малки. В замисимост от желаното качество се избират тегловни (квантуващи) коефициенти, на които се разделят DCT коефициентите.
След преминаване в честотната област лесно може да се премахне информацията за честотите, към които човешкото възприятие е по-малко чувствително.
- Високите честоти отговарят на резки преходи, контрастни линии, дребни обекти.
- Честоти с малки коефициенти (малка аплитуда) остават незабелижими в цялостното изображение
Квантуващата матрица е генирирана в съотвествие качеството, което желаем. При квантуването коефициентите с малка стойност се нулират нула. Квантуваща матрица за яркостния сигнал при качество 50%
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 - дискретна косиносова трансформация
