function [NEWA,NEWMI,NEWSIGMA,Ptot]=reestim4(O1,O2,O3,O4,A,MI,SIGMA);
%
% Upravená verze pro 4 vzorky

Ptot=ones(4,1);
N=size(A,1);
[P1,T1]=size(O1);
[P2,T2]=size(O2);
[P3,T3]=size(O3);
[P4,T4]=size(O4);

%První
[ALFA1,Ptot(1)]=GetAlfa(O1,A,MI,SIGMA);
[BETA1]=GetBeta(O1,A,MI,SIGMA);

%Druhá
[ALFA2,Ptot(2)]=GetAlfa(O2,A,MI,SIGMA);
[BETA2]=GetBeta(O2,A,MI,SIGMA);

%Třetí
[ALFA3,Ptot(3)]=GetAlfa(O3,A,MI,SIGMA);
[BETA3]=GetBeta(O3,A,MI,SIGMA);

%Čtvrtá
[ALFA4,Ptot(4)]=GetAlfa(O4,A,MI,SIGMA);
[BETA4]=GetBeta(O4,A,MI,SIGMA);

%Get New A

NEWA=zeros(size(A));

for i=2:(N-1),
for j=2:(N-1),
up=0; down=0;
for t=1:(T1-1),
up=up+A(i,j)*normal(O1(:,t+1),MI(:,j),SIGMA(:,j))*ALFA1(i,t)*BETA1(j,t+1);
end
for t=1:(T2-1),
up=up+A(i,j)*normal(O2(:,t+1),MI(:,j),SIGMA(:,j))*ALFA2(i,t)*BETA2(j,t+1);
end
for t=1:(T3-1),
up=up+A(i,j)*normal(O3(:,t+1),MI(:,j),SIGMA(:,j))*ALFA3(i,t)*BETA3(j,t+1);
end
for t=1:(T4-1),
up=up+A(i,j)*normal(O4(:,t+1),MI(:,j),SIGMA(:,j))*ALFA4(i,t)*BETA4(j,t+1);
end

for t=1:T1,
down=down+ALFA1(i,t)*BETA1(i,t);
end
for t=1:T2,
down=down+ALFA2(i,t)*BETA2(i,t);
end
for t=1:T3,
down=down+ALFA3(i,t)*BETA3(i,t);
end
for t=1:T4,
down=down+ALFA4(i,t)*BETA4(i,t);
end
NEWA(i,j) = up / down;
end
end

for j=2:(N-1),
up=ALFA1(j,1)*BETA1(j,1);
up=up+ALFA2(j,1)*BETA2(j,1);
up=up+ALFA3(j,1)*BETA3(j,1);
up=up+ALFA4(j,1)*BETA4(j,1);
down=sum(Ptot);
NEWA(1,j) = up / down;
end

down=0;
for i=2:(N-1),
% for j=2:(N-1),
up=ALFA1(i,T1)*BETA1(i,T1) + ALFA2(i,T2)*BETA2(i,T2) + ...
ALFA3(i,T3)*BETA3(i,T3) + ALFA4(i,T4)*BETA4(i,T4);
for t=1:T1,
down=down + ALFA1(i,t)*BETA1(i,t);
end
for t=1:T2,
down=down + ALFA2(i,t)*BETA2(i,t);
end
for t=1:T3,
down=down + ALFA3(i,t)*BETA3(i,t);
end
for t=1:T4,
down=down + ALFA4(i,t)*BETA4(i,t);
end
NEWA(i,N) = up / down;
% end
end

%Get New B

L1=zeros(N,T1);
L1=ALFA1.*BETA1/Ptot(1);

L2=zeros(N,T2);
L2=ALFA2.*BETA2/Ptot(2);

L3=zeros(N,T3);
L3=ALFA3.*BETA3/Ptot(3);

L4=zeros(N,T4);
L4=ALFA4.*BETA4/Ptot(4);

NEWMI=zeros(size(MI));
NEWSIGMA=zeros(size(SIGMA));
for j=2:(N-1),
L1vec=L1(j,:); L1mat=ones(P1,1)*L1vec;
L2vec=L2(j,:); L2mat=ones(P2,1)*L2vec;
L3vec=L3(j,:); L3mat=ones(P3,1)*L3vec;
L4vec=L4(j,:); L4mat=ones(P4,1)*L4vec;

UP= sum((L1mat.*O1)')';
UP=UP+sum((L2mat.*O2)')';
UP=UP+sum((L3mat.*O3)')';
UP=UP+sum((L4mat.*O4)')';

DOWN=sum(L1vec)+sum(L2vec)+sum(L3vec)+sum(L4vec);

NEWMI(:,j) = UP / DOWN;

UP= sum((L1mat.*((O1-NEWMI(:,j)*ones(1,T1)).^2))')';
UP=UP+sum((L2mat.*((O2-NEWMI(:,j)*ones(1,T2)).^2))')';
UP=UP+sum((L3mat.*((O3-NEWMI(:,j)*ones(1,T3)).^2))')';
UP=UP+sum((L4mat.*((O4-NEWMI(:,j)*ones(1,T4)).^2))')';

NEWSIGMA(:,j) = UP / DOWN;
end