function [l, pk] = hmm_get_probs(hparm, xx, use_diff);
% function [l, pk] = hmm_get_probs(hparm, xx, [use_diff=0]);
%
% gets raw observation state probabilities (for each state and for each
% time step).
%
%
if(nargin<3), use_diff=0; end;

%
% NAVAL UNDERSEA WARFARE CENTER
% class-specific@npt.nuwc.navy.mil
% http://www.npt.nuwc.navy.mil/csf
% Paul M. Baggenstoss
% 401-832-8240
% Do not disseminate.
%



doTrain = 0; % FIXME, should be passed in perhaps?


BAD = -9e37;
N=hparm.N;
[dim,nsamp]=size(xx);

% Compute probabilities
pk={};
l= BAD*ones(nsamp,N); % CAS: not needed, l is overwritten, zeros() good enough.
for istate=1:N,
    if(use_diff),
       gparm = hparm.pdfd(istate);
    else,
       gparm = hparm.pdf(istate);
    end;
    nmode = length(gparm.modes);
    wts = [gparm.modes.weight]';
    
    % log-PDF for each mode
    lg = lqr_evp(gparm, xx, 1);
    
    % subtract MDL correction if training
    if(doTrain), lg=lg - hparm.mdl(istate); end;
    
    % exponentiate
    if(nmode>1),
        lg = lg + log(repmat(wts(:)',nsamp,1)); % CAS: FIXME: log BEFORE repmat
        mx = max(lg,[],2);
        lg=lg-repmat(mx,1,nmode);
        g = exp(lg);
        psum = sum(g,2);
        lg = log(psum) + mx;
        p = g ./ repmat(psum,1,nmode);
    else,
        p = ones(nsamp,1);
    end;
    pk{istate}=p;
    l(:,istate)=lg;
end

