function p = probdev0(A,C,phi)
% the probability that the rv. phi(X), with X ~ Markov(A,C) 
% deviates from its mean 
% by at least delta*n
%%by 1

[m,m,t] = size(A);
if ~exist('C','var')
  C = ones(m,1)/m;
end
n = t+1;
dims = repmat([m],[1 n]);

P = zeros(m^n,1);
F = zeros(m^n,1);
for xind=1:m^n
  x = ind2coord(xind,dims);
  P(xind) = Pr(A,C,x);
  F(xind) = phi(xind);
end

Ef = P'*F;
del = abs(F-Ef);

%ii = find(del==1);
%p = sum(P(ii));

global delta
ii = find(del>=delta*n);
p = sum(P(ii));
return

D = unique(del);
pD = zeros(size(D));
for i=1:length(D)
  jj = find(del==D(i));
  pD(i) = sum(P(jj));
end
plot(D,pD,'.');

return


function p = Pr(A,C,x)
%p = 1;
p = C(x(1));
for j=2:length(x)
  p = p * A(x(j),x(j-1),j-1);
end
return