function[pocet] = sturm(P, a, b)

n = length(P);
p1 = P;
p2 = -polyder(P); %derivacia polynomu

pol = [p1; [0 p2]]; %Matica majuca za prvky koeficienty danych polynomov

[q, r] = deconv(p1, p2);%delenie polynomu so zvyskom
r = r(find(r,1):end);


while ~all(P == 0) & length(r) > 1
    p1 = p2;
    p2 = -r;
    pol = [pol; [zeros(1, n-length(r)), -r]];
    [q, r] = deconv(p1, p2);
    r = r(find(r,1):end);
end

pol = [pol; [zeros(1, n-length(r)), -r]];

Pa = []; Pb = [];
for i=1:size(pol,1)
    Pa(i) = polyval(pol(i,:), a); 
    Pb(i) = polyval(pol(i,:), b); 
end

Wa = pocetzmen(Pa);
Wb = pocetzmen(Pb);

pocet = Wb - Wa;

end
%Polynom, s ktorym sme pracovali na seminari
%P =[1 0 0 -4 1]