Wstęp
Rozdział traktować będzie o algorytmach porównywania z oryginałem skompresowanych obrazów i – na podstawie wyników – ocena, który obraz jest lepszej jakości. Innymi słowy – syntetyczny test jakości kompresji. Obraz wejściowy kompresowany jest do JPEG za pomocą wewnętrznej biblioteki Matlab. Konwersja z i do JPEG2000 wykonywana jest za pomocą zewnętrznej biblioteki OpenJPEG (potrzebne są pliki image_to_j2k.exe oraz j2k_to_image.exe).
Kod programu
porownaj.m
clear;[pictureone, user_canceled1] = imgetfile;width = 640;resultJPEG = [;]; resultJZIP = [;]; resultJ2K = [;];timeresultJPEG = [;]; timeresultJZIP = [;]; timeresultJ2K = [;];jpegMinQ = 0; jpegMaxQ = 100; jpegStep = 1;jp2kMinQ = 15; jp2kMaxQ = 50; jp2kStep = 1;compareMethod = 4;maxDisplayedSize = 15000;maxDifference = -36;usePSNR = false;useSize = true;terminateAfterMatch = false; if (user_canceled1) fprintf('nAnulowales wybieranie plikow do porownania!n')else %create bmp temp file tempRead = imread(pictureone); imwrite(tempRead,'temp.bmp','bmp'); %reading bmp file Aread = imread('temp.bmp'); A = Aread; %JPEG conversion for i=jpegMinQ:jpegStep:jpegMaxQ tic; imwrite(A,'temp.jpg','jpg','Quality', i); onlyJpegTime = toc; syscmd = '7za.exe a -mx9 temp.7z temp.jpg'; system(syscmd); compressionTime = toc; fileInfo = dir('temp.7z'); fileSize = fileInfo.bytes; Bread = imread('temp.jpg'); B = Bread; %if (fileSize <= maxDisplayedSize) % JPEGsaved = B; % JPEGsavedSize = fileSize; %end temp = compare(A,B,compareMethod); resultJZIP = cat(2,resultJZIP,[fileSize;temp]); timeresultJZIP = cat(2,timeresultJZIP,[fileSize;compressionTime]); fileInfo2 = dir('temp.jpg'); fileSize2 = fileInfo2.bytes; resultJPEG = cat(2,resultJPEG,[fileSize2;temp]); timeresultJPEG = cat(2,timeresultJPEG,[fileSize2;onlyJpegTime]); fprintf('nWykonano JPEG quality=%d...n',i); if ((temp >= maxDifference && usePSNR==true)||(fileSize <= maxDisplayedSize && useSize==true)) JPEGsaved = B; JPEGsavedSize = fileSize; else if (terminateAfterMatch) break; end end end %JPEG2000 conversion for i=jp2kMinQ:jp2kStep:jp2kMaxQ tic; syscmd = ['image_to_j2k.exe -i temp.bmp -o temp.jp2 -q ', int2str(i)]; system(syscmd); compressionTime = toc; fileInfo = dir('temp.jp2'); fileSize = fileInfo.bytes; syscmd = 'j2k_to_image.exe -i temp.jp2 -o temp2.bmp'; system(syscmd); B2k = imread('temp2.bmp'); %if (fileSize <= maxDisplayedSize) % JPEG2ksaved = B2k; % JPEG2ksavedSize = fileSize; %end temp = compare(A,B2k,compareMethod); resultJ2K = cat(2,resultJ2K,[fileSize;temp]); timeresultJ2K = cat(2,timeresultJ2K,[fileSize;compressionTime]); fprintf('nWykonano JPEG2000 q=%d...n',i); if ((temp >= maxDifference && usePSNR==true)||(fileSize <= maxDisplayedSize && useSize==true)) JPEG2ksaved = B2k; JPEG2ksavedSize = fileSize; else if (terminateAfterMatch) break; end end end %remove temporary files system('del temp*.*'); figure('Name','Compression quality plot'), plot(resultJPEG(1,:),resultJPEG(2,:),'mo','MarkerSize',2) hold all plot(resultJZIP(1,:),resultJZIP(2,:),'rd','MarkerSize',2) plot(resultJ2K(1,:),resultJ2K(2,:),'bs','MarkerSize',2) xlabel('JPEG size [bytes]') ylabel('Method dependent error') title(pictureone) h = legend('JPEG','JPEG+ZIP','JPEG2000'); set(h,'Interpreter','none') figure('Name','Compression time plot'), plot(timeresultJPEG(1,:),timeresultJPEG(2,:),'mo','MarkerSize',2) hold all plot(timeresultJZIP(1,:),timeresultJZIP(2,:),'rd','MarkerSize',2) plot(timeresultJ2K(1,:),timeresultJ2K(2,:),'bs','MarkerSize',2) xlabel('JPEG size [bytes]') ylabel('Time of compression') title(pictureone) h = legend('JPEG','JPEG+ZIP','JPEG2000'); set(h,'Interpreter','none') figure('Name','Source picture, JPEG and JPEG2000'), imshow([Aread,JPEGsaved,JPEG2ksaved])end
compare.m
function [out1] = compare(inA, inB, compareMethod) % simple difference; result is a sum of differences between % original pixel value and compressed one; % result for grayscale images are multiplied by 3 if (compareMethod == 1) outC = abs(inA-inB); [width height planes] = size(inA); out1 = sum(sum(sum(outC)))/(width*height); % Cartesian difference; result is a sum of Cartesian differences % between original pixel values and compressed one; % result for grayscale images are multiplied by sqrt(3) elseif (compareMethod == 2) outC = abs(inA-inB); [width height planes] = size(inA); outC = outC.^2; outD = sum(outC,3); outD = outD.^0.5; out1 = sum(sum(sum(outD)))/(width*height); % same as method 2, but with tolerance; % results are scaled by RGB width (255*sqrt(3)), % presented in log10 scale elseif (compareMethod == 3) tolerance = 3; outC = abs(inA-inB); [width height planes] = size(inA); outC = outC.^2; outD = sum(outC,3); outD = outD.^0.5; outD(outD < tolerance) = 0; outD = outD./(255*(3^0.5)); out1 = sum(sum(sum(outD)))/(width*height); out1 = 20*log10(out1); % weighted pseudoPSNR with FIR filter elseif (compareMethod == 4) outC = abs(inA-inB); outC = outC.^2; outD = sum(outC,3); outD = outD.^0.5; weightedMatrix = [1,2,1;3,7,3;1,2,1]; weightedMatrixNorm = weightedMatrix/sum(sum(weightedMatrix)); outE = filter2(weightedMatrixNorm,outD,'same'); outE = outE./(255*(3^0.5)); out1 = sum(sum(sum(outE)))/(size(outE,1)*size(outE,2)); out1 = 20*log10(out1); else out1 = 0; end end
JPEG kontra JPEG2000 - pierwsze testy
W rezultacie otrzymany zostaje wykres porównujący ze sobą dwa skompresowane obrazy, bazując na zaproponowanym algorytmie. Za pierwszy obraz testowy posłużyły "Lilie_wodne.jpg", dostępne jako jeden z przykładowych obrazów dla systemu Windows XP. Wykresy dla trzeciego algorytmu porównującego znajdują się poniżej:
Przybliżenie początkowej fazy wykresu:
Poniżej przedstawiono zbiorcze porównanie obrazów otrzymanych dla największego rozmiaru pliku poniżej 20 kB (kliknij aby powiększyć):
{modal href="images/stories/result_lilie.jpg_20k.png"}{/modal}
Na porównaniu wyraźnie widać, że oba skompresowane obrazy straciły znacznie na jakości. Ze względu jednak na charakter kompresji JPEG (bloki 8*8 pikseli) obraz w tym przypadku wydaje się być zdecydowanie mniej czytelny zarówno w miejscach o dużej szczegółowości jak i małej. JPEG2000 większość "jakości" obrazu traci w obszarach mniej szczegółowych, przez co efekt jest ten mniej zauważalny, a pozwala poprawić jakość tam, gdzie szczegółów jest więcej. Wraz ze zmniejszaniem się rozmiaru pliku (wzrost stopnia kompresji) różnica pomiędzy JPEG i JPEG2000 staje się coraz bardziej wyraźna:
{modal href="images/stories/result_lilie.jpg_15k.png"}{/modal}
Dla porównania, poniżej znajduje się fragment identycznie skompresowanego zdjęcia "Zachód słońca.jpg" (również znajduje się w przykładowych obrazach Windows XP). Rozmiar docelowy to 15 kB:
{modal href="images/stories/result_zachod.jpg_15k.png"}{/modal}
Wykres przedstawiający zależność jakości obrazu od rozmiaru skompresowanego pliku dla wybranego algorytmu porównującego i obrazu "Zachód słońca.jpg":
Kolejne testy - JPEG+PAQ8fthis2 kontra JPEG2000
Poniżej prezentuję testy porównujące kompresję JPEG2000 z połączeniem JPEG + PAQ8fthis2. Tym razem kompresji poddana została koala z przykładowych zdjęć Windows 7. Dodatkowym parametrem był mierzony czas kompresji. Użyty został czwarty algorytm porównywania obrazów, biorący pod uwagę nie tylko aktualny piksel ale i jego "sąsiadów". Pierwszy obrazek prezentuje kompresję do JPEG i JPEG2000 dla plików wynikowych o rozmiarach odpowiednio 14,5 kB i 13,3 kB. Plik źródłowy zajmował 762 kB.
{modal href="images/stories/result_koala_15k_paq8f.png"}{/modal}
{modal href="images/stories/result_koala_15k_plot_paq8f.png"}{/modal}
Z tej części eksperymentu wyraźnie widać, że nawet dla bardzo silnej, czasochłonnej i zasobożernej (tego ostatniego nie widać, ale logi Matlab-a pokazują użycie pamięci dla PAQ8f na poziomie 2 GB) opcji kompresji plików JPEG, rezultaty znacznie odbiegają od możliwości JPEG2000.
Proponowane algorytmy porównujące
Algorytm 1
outC = abs(inA-inB);[width height planes] = size(inA);out1 = sum(sum(sum(outC)))/(width*height);
Algorytm 2
outC = abs(inA-inB);[width height planes] = size(inA);outC = outC.^2;outD = sum(outC,3);outD = outD.^0.5;out1 = sum(sum(sum(outD)))/(width*height);;
Algorytm 3
tolerance = 3;outC = abs(inA-inB);[width height planes] = size(inA);outC = outC.^2;outD = sum(outC,3);outD = outD.^0.5;outD(outD < tolerance) = 0;outD = outD./(255*(3^0.5));out1 = sum(sum(sum(outD)))/(width*height);out1 = 20*log10(out1);
Algorytm 4
outC = abs(inA-inB);outC = outC.^2;outD = sum(outC,3);outD = outD.^0.5;weightedMatrix = [1,2,1;3,7,3;1,2,1];weightedMatrixNorm = weightedMatrix/sum(sum(weightedMatrix));outE = filter2(weightedMatrixNorm,outD,'same');outE = outE./(255*(3^0.5));out1 = sum(sum(sum(outE)))/(size(outE,1)*size(outE,2));out1 = 20*log10(out1);
Algorytm od poprzedniego różni się podejściem do liczenia odchyłki dla piksela. Brany jest nie tylko piksel centralny ale również piksele sąsiednie. Sąsiednim pikselom przypisane są wagi:
{tex}mathrm{weightedMatrix} = left [ begin{matrix}1 & 2 & 1\ 3 & 7 & 3\ 1 & 2 & 1end{matrix} right ]{/tex}
Po podzieleniu przez sumę wszystkich elementów otrzymamy znormalizowaną macierz wag, która wykorzystana została jako filtr FIR macierzy wejściowej. W powyższym przykładzie wagi ustawione są tak aby największa uwaga skupiała się na elemencie centralny, mniejsze na sąsiednich, przy czym uwzględniona została pozioma i pionowa zdolność rozdzielcza ludzkiego oka [do zweryfikowania].
{tex}mathrm{weightedMatrixNorm} = left [ begin{matrix}0.0476 & 0.0952 & 0.0476\0.1429 & 0.3333 & 0.1429\0.0476 & 0.0952 & 0.0476end{matrix} right ]{/tex}