%Author: Ashish Kapoor MIT Media Lab, 10/20/02
%Please send corrections to: ash@media.mit.edu

function pr = log_prob_obs_hmm(pi, a, b, num_states, obs)
%Function to compute probability of the observation given the parameters of HMM (P(O|Lambda))

%Memory Allocation
T = size(obs,1);
alpha_hat = zeros(T, num_states);
alpha_hat_hat = zeros(T, num_states);
c = zeros(T,1);

%**********************************************************************************
%Initialization

%Initializing alpha_hat_hat; 
%****AVOID FOR LOOPS****
%To DO:
%for i=1:num_states
%    alpha_hat_hat(1,i)=pi(i)*b(i,obs(1)+1);
%end
%Eqivalent MATLAB code:
alpha_hat_hat(1,:)=pi'.*b(:,obs(1)+1)';

c(1)=1/sum(alpha_hat_hat(1,:));

alpha_hat(1,:)=c(1)*alpha_hat_hat(1,:);

%**********************************************************************************
%Recursion
for t=2:T
    
    %Again avoiding for loop
    %for j=1:num_states
    %    alpha_hat_hat(t,j)=sum(alpha_hat(t-1,:).*a(:,j)'*b(j,obs(t)+1));
    %end
    %Equivalent line of code below
    alpha_hat_hat(t,:)=(alpha_hat(t-1,:)*a(:,:)).*b(:,obs(t)+1)';
    
    c(t)=1/sum(alpha_hat_hat(t,:));
    
    alpha_hat(t,:)=c(t)*alpha_hat_hat(t,:);
       
end

%**********************************************************************************
%Termination
pr=-sum(log(c));
