function [y,xWTwien]=wiener(x,xp,oddo,wfiltr,PocetUrovni)
% wienerovska vlnkova filtrace
% x ... vstupni signal
% xp ... signal pilotniho odhadu
% wfiltr ... banka filtru pro wienerovskou filtraci, napr. 'bior2.2'
% PocetUrovni ... pocet stupnu rozkladu 
% sumod, sumdo ... indexy pro odhad urovne sumu (od, do)
%
% vystupy:
% y ... filtrovany signal
% xWTwien ... vektor DTWT vstupu po wienerovske korekci
%         pasma pro PocetUrovni=4: ... 4. 3. 2. 1. a 0. (0. je vystup DP)


sumod=oddo(1);
sumdo=oddo(2);


% ********** Wieneruv filtr s pilotnim odhadem uzitecnych koeficientu
% **********
% vstup: x
% vstup pilotniho odhadu: xp
% postup: xDWT=DTWT(x)
%         xpDWT=DTWT(xp)

[xDWT, index, tau1] = NonDwt(x,wfiltr,PocetUrovni,1); % transformace vstupu
[xpDWT, index, tau1] = NonDwt(xp,wfiltr,PocetUrovni,1);  % transformace pilotniho signalu


% ---------------------------------------------------------------------------
% filtrace Wienerovym filtrem
% odhad smerodatne odchylky sumu uz byl proveden ***** zde upravuji *****

pasem = PocetUrovni+1;
z = 1;
for k = 1:pasem
    z = z+index(k);
    zac(k) = z;
end

delka=sumdo-sumod+1;

for k = PocetUrovni : -1 : 1
    stred=sum(xDWT(zac(k)+sumod:zac(k)+sumdo))/delka;
    sigma=sqrt(sum((xDWT(zac(k)+sumod:zac(k)+sumdo)-stred).^2)/delka);
    for i=zac(k):zac(k)+index(k+1)-1
        xDWT(i)=xDWT(i)*(xpDWT(i)*xpDWT(i)/(xpDWT(i)*xpDWT(i)+sigma*sigma));
    end
end

xWTwien = xDWT;

% -----------------------------------------------------------------------
% IDWT, korekce zpozdeni, kresby

[xIDWT, tau2] = NonIDwt(xWTwien,index,wfiltr,PocetUrovni);

xIDWT = CancelDelay(xIDWT,wfiltr,PocetUrovni,tau1,tau2);

y = xIDWT;



