%% This file compares the output of the Software and Hardware DWT.
%% Author : Rami Zewail
%% Group  : JPEG2000 Project


% this is the output image from the ModelSim testbench,
fid=fopen('DWTresult.txt');
aa=fscanf(fid,'%2x %2x',[64 inf]);
aa=aa';
fclose(fid);
%% Reverse the DC level shifting step.
hw=aa;
for i=1:64
for j=1:64
if hw(i,j)>127
hw(i,j)=hw(i,j)-256;
end;
end;
end;
figure;imshow(mat2gray(hw));
%%%%%%

fid= fopen('testdata.txt');

a=fscanf(fid,'%2x %2x',[64 inf]);

a=a';

for i=1:64
for j=1:64
if a(i,j)>127
a(i,j)=a(i,j)-256;
end;
end;
end;

fclose(fid);
x_image=a;

 image_1d=dwt_1d(x_image);
 image_transpose=image_1d';
 image_2d=dwt_1d(image_transpose);
 image_2d=image_2d';
%%%%%%%%%%
sw_testdata=image_2d;
%%%%%%%%%%%%%
%Calculates the PSNR for Each Sub-band and the histogram

hw=(hw);
sw=(sw_testdata);

[M,N]=size(hw);
hw_A=(mat2gray(hw(1:M/2,1:N/2)));
figure;imhist(hw_A);
hw_H=(mat2gray(hw(1:M/2,N/2+1:N)));
figure;imhist(hw_H);
hw_V=(mat2gray(hw(M/2+1:M,1:N/2)));
figure;imhist(hw_V);
hw_D=(mat2gray(hw(M/2+1:M,N/2+1:N)));
 figure;imhist(hw_D);

sw_A=(mat2gray(sw(1:M/2,1:N/2)));
figure;imhist(sw_A);
sw_H=(mat2gray(sw(1:M/2,N/2+1:N)));
figure;imhist(sw_H);
Sw_V=(mat2gray(sw(M/2+1:M,1:N/2)));
figure;imhist(sw_V);
sw_D=(mat2gray(sw(M/2+1:M,N/2+1:N)));
figure;imhist(sw_D);
mse_A=0;mse_H=0;mse_V=0;mse_D=0;
for i=1:M
    for j=1:N
        mse_A=mse_A+((hw_A(i,j)-sw_A(i,j))^2);
        
        mse_H=mse_H+((hw_H(i,j)-sw_H(i,j))^2);
        
        mse_V=mse_V+((hw_V(i,j)-sw_V(i,j))^2);
        
        mse_D=mse_D+((hw_D(i,j)-sw_D(i,j))^2);
    end;
end;
mse_A=mse_A/(M*N);
mse_H=mse_H/(M*N);
mse_V=mse_V/(M*N);
mse_D=mse_D/(M*N);

%%%%%%%%%%%%%%%%%
psnr_A= 20*log10(2^8-1/sqrt(mse_A))
psnr_H= 20*log10(2^8-1/sqrt(mse_H))
psnr_V= 20*log10(2^8-1/sqrt(mse_V))
psnr_D= 20*log10(2^8-1/sqrt(mse_D))

%%%%%%%%%%%%%%%