function [Ef,K] = condexp(P,F,Xt,m,n)
% computes the conditional expectation
% E[f(X)|Xt]
% a numerically stable routine

t = length(Xt);
dimt = repmat([m],[1 n-t]);
dims = repmat([m],[1 n]);
%Snt =  list_all_ind(dimt);
%Sn  =  list_all_ind(dims);

pXt = getPXt(Xt,P,m,n);
Ef = 0;
K = zeros(size(P));
for xint = 1:m^(n-t)
  xtn = ind2coord(xint,dimt);
  x = [Xt xtn];
  %xprnt = x-1
  xins = coord2ind(x,dims);
  f = F(xins);
  p = P(xins)/(pXt+realmin);
  K(xins) = p;
  Ef = Ef + f*p;
end
return



function p = getPXt(Xt,P,m,n)
t = length(Xt);
if t==0
  p = 1;
  return
end
dimt = repmat([m],[1 n-t]);
dims = repmat([m],[1 n]);
p = 0;
for xint = 1:m^(n-t)
  xtn = ind2coord(xint,dimt);
  x = [Xt xtn];
  xins = coord2ind(x,dims);
  p = p + P(xins);
end
return

% BAD, saved only as an illustration of the 
% indexation discrepancy between list_all_ind
% and ind2coord
function p = getPXt00(Xt,Sn,P,m,n)
t = length(Xt);
if t==0
  p = 1;
  return
end
inds = 1:m^n;
for i=1:t
  indi = find(Sn(:,i)==Xt(i));
  inds = intersect(inds,indi);
end
p = sum(P(inds));
return