function vpath = viterbi(a, pi, b, o)

vpath = [];
spi = size(pi');
so = size(o');

T = so(1);
N = spi(1);
for i = 1:N
	delta(1, i) = pi(i)*b(i,o(1));
	psi(1,i) = 0;
end

for t=2:T
	for j = 1:N
		pref = 0;
		argmax = 1;
		for i = 1:N
			if((delta(t-1, i)*a(i,j)) > pref)
				pref = delta(t-1,i)*a(i,j);
				argmax = i;
			end
		end
		delta(t,j) = pref*b(j, o(t));
		psi(t, j) = argmax;
	end
end

m = 0;
for i = 1:N
	if(delta(T, i) > m)
		m = delta(T,i);
		argmax = i;
	end
end

qstar(T) = argmax;

vpath = [argmax];
for t = (T-1):-1:1
	newp = psi(t+1, qstar(t+1));
	qstar(t) = newp;
	vpath = [newp vpath];
end
