function E = maxdev(P,F,m,n)
% computes max_{f\in F} Pf - P_n f
% where f: X -> R
% each column of F is such an f

P1 = getmargX(P,m,n,1); % process is stationary

E = 0;
for xi=1:m^n
  x = i2c(xi,m,n);
  PF = P1' * F;
  
  Fx = getFX(F,x);
  
  Z = max(abs(PF-Fx'));
  E = E + Z * P(xi);
  
end

return

function Fx = getFX(calF,x)
[m,F] = size(calF);
n = length(x);
FX = zeros(F,n);
for f=1:F
  for i=1:n
    FX(f,i) = calF(x(i),f);
  end
end
Fx = mean(FX,2);
return


function c=i2c(i,m,n)
dims = repmat(m,[1 n]);
c = ind2coord(i,dims);
return