close all
% (1) Propuštanje zašumljenog signala kroz sustav
%zadavanje istog signala
x = linspace(0,4*pi,1024);
ulaz_orig=diric(x,7);

%zadavanje bijelog šuma
bijeli_sum=randn(1,1024);

%ulazni zašumljeni signal
ulaz=ulaz_orig+bijeli_sum/10;
figure(1),plot(ulaz_orig),title('Ulaz')
figure(2),plot(ulaz),title('Ulaz_sum')


%wavelet transformacija(prva) ulaznog signala
a=swt(ulaz,5,'db1');
figure(5),imagesc(a);

%razliite metode raunanja optimalne razine praga
% thr1 = thselect(a, 'rigrsure');
% thr2 = thselect(a, 'heursure');
% thr3 = thselect(a, 'sqtwolog');
% thr4 = thselect(a, 'minimaxi');

%soft threshold filtriranje u wavwlet domeni
%a_tresh=wthresh(a,'s',thr);
%mi smo sami izabrali prag koji daje najbolje rezultate
a_tresh=wthresh(a,'s',0.01);

%inverzna wavelet transformacija
s1=iswt(a_tresh,'db1'); %signal na izlazu iz gornje grane
figure(3), plot(s1); title('Izlaz gornje grane')

%wavelet transformacija(druga) signala na izlazu iz gornje grane i
%ulaznog signala
a2=swt(s1,5,'sym4');
a3=swt(ulaz,5,'sym4');

%raunanje wienerovog filtra
broj=a2.^2;
naz=broj+(var(ulaz))^2;
hw=broj./naz;
theta2=a3.*hw;

%inverzna wavelet transformacija
izlaz=iswt(theta2,'sym4'); %signal na izlazu
figure(4), plot(izlaz); title('Izlaz');

SNRi = 10*log10(sum(abs(fft(ulaz_orig)).^2)/sum(abs(fft(bijeli_sum/10)).^2));
SNRo = 10*log10(sum(abs(fft(ulaz_orig)).^2)/sum(abs(fft(izlaz-ulaz_orig)).^2));
% %SNR na ulazu
% % SNR = 10*log10(Ps/Pe)
% % Ps = Sum of s squared
% % Pe = Sum of e squared
% Psi=(sum(ulaz_orig)*sum(ulaz_orig)); % setting Ps to sum of s squared
% Pei=(sum(bijeli_sum/10)*sum(bijeli_sum/10)); % setting Pe to sum of e squared
% SNRi=10*log10(Psi/Pei);
% 
% %SNR na izlazu
% Pso=(sum(ulaz_orig))*(sum(ulaz_orig)); % setting Ps to sum of s squared
% Peo=(sum(izlaz-ulaz_orig))*(sum(izlaz-ulaz_orig)); % setting Pe to sum of e squared
% SNRo=10*log10(Pso/Peo);
% 
% figure(6), plot(bijeli_sum/10);
% figure(7), plot(izlaz-ulaz_orig);

