function [resenje]=PerioSplineInterp(cvorovi,vrednosti,vektor)

% PerioSplineInterp interpolacija funkcije periodicnim kubnim splajnom
% Tablica interpolacije mora da sadrzi bar tri cvora, svi cvorovi moraju biti razliciti,
% u neopadajucem poretku i 
% mora biti jednaka (uslov periodicnosti)
%
% Broj argumenata je unapred zadat
% Argumenti funkcije su :
%   cvorovi (1..n)      - vektor argumenata
%   vrednosti (1..n)    - vektor vrednosti funkcije
%   vektor(1..k)        - vektor za koji zelimo da se izracunaju vrednosti interpolacionog splajna
% Funkcija vraca:
%   vektor vrednosti kubnog periodicnog interpolacionog splajna u standardnom MATLAB-ovom zapisu 

% Autor : Miodrag Jankovic, MR 92/142
% Verzija : 1.0
% Datum : 01.12.2002.

%Provera uslova zadatka

if (~exist('vrednosti')) | (~exist('cvorovi'))
    error ('nisu uneti neophodni podaci. Pravilan poziv je: PeriodicalSplineInterpolation(cvorovi,vrednosti,vektor)')
end

[i,k]= size(cvorovi);
if k~=1
    error('skup argumenata mora biti vektor')
end

if (diff(cvorovi)<0)
    error('vrednosti argumenata nisu u neopadajucem poretku')
end

[i,k]= size(vrednosti);
if k~=1
    error('skup vrednosti funkcije mora biti vektor')
end

if length(cvorovi)<3
    error('Potrebno je da postoje bar dva cvora interpolacije')
end

%provera ispravnog broja argumenata i vrednosti funkcije
if length(cvorovi)~=length(vrednosti)
    error('Broj argumenata funkcije se ne poklapa sa brojem zadatih vrednosti funkcije')
end

% provera razlicitost cvorova interpolacije
if (~isequal(cvorovi,unique(cvorovi)))
    error('Svi cvorovi interpolacije bi trebalo da budu razliciti')
end

% provera uslova periodicnosti
if vrednosti(1)~=vrednosti(length(vrednosti))
    error ('Zadata funkcija nije periodicna')
end


resenje=[];
k=size(cvorovi);      	%odredjivanje dimenzije ulazne matrice
k=k(1);

h=diff(cvorovi);
dh=diff(vrednosti);     %zbog skracenog pisanja nadalje
h(k)=h(1);              %jer je splajn periodican


%izracunavanje koeficijenata h, ni, mi i lambda
for i=1:k-2
	ni(i)=h(i+1)/(h(i+1)+h(i));
	mi(i)=1-ni(i);
	lambda(i)=6*(dh(i+1)/h(i+1)-dh(i)/h(i))/(h(i)+h(i+1));
end
ni(k-1)=h(1)/(h(k)+h(1));
mi(k-1)=1-ni(k-1);
lambda(k-1)=6*(dh(1)/h(1)-dh(k-1)/h(k-1))/(h(k-1)+h(1));


%kreiranje trodijagonalne matrice sistema i resavanje sistema
ni=sparse(1:k-1,1:k-1,ni);
mi=sparse(1:k-1,1:k-1,mi);
lambda;
mi=[mi(:,2:k-1) mi(:,1)];   %pomeranje u desno
ni=[ni(:,k-1) ni(:,1:k-2)]; %pomeranje u levo
sistem=sparse(1:k-1,1:k-1,2)+ni+mi    %matrica sistema
M=sistem\lambda'
M=[M; M]     
%nadovezujem zato sto ne mogu da koristim M(0), pa cu da primenim cinjenicu da je zbog periodicnosti M(0)=M(n)


%racunanje koeficijenata splajna
for i=1:k-1
    alfa(i)=vrednosti(i);
    beta(i)=dh(i)/h(i+1)-h(i+1)/6*(2*M(i+k-2)+M(i+k-1));
    gama(i)=0.5*M(i+k-2);
    delta(i)=(M(i+k-1)-M(i+k-2))/6*h(i+1);
end

koef=[alfa' beta' gama' delta']
koef=fliplr(koef);
test=[cvorovi(1:k-1) cvorovi(2:k)];

for i=1:length(vektor)
resenje =[resenje; vektor(i) polyval(koef(find (test(:,1)<=vektor(i) & test(:,2)>vektor(i)),:),vektor(i)-cvorovi(find (test(:,1)<=vektor(i) & test(:,2)>vektor(i))))];
end