function [hh,Ht,H] = mxgethhn(P,m,n)
hh = zeros(n);
for i=1:n
  for j=i+1:n
    h = mxetaij(P,m,n,i,j);
    hh(i,j) = h;
  end
end
Ht = 1+sum(hh,2)';
H = max(Ht);
hh0 = gethhn(P,m,n);
err = totdif(hh,hh0)
return


function [hh,Ht,H] = gethhn(P,m,n)
hh = zeros(n);
for i=1:n
  for j=i+1:n
    h = etaij(P,m,n,i,j);
    hh(i,j) = h;
  end
end
Ht = 1+sum(hh,2)';
H = max(Ht);
return


function [h,xim] = etaij(P,m,n,i,j)
global A B C mobs mhid hhat

h = 0;
xim = 0;
dims = repmat([m],[1 i-1]);
for xi = 1:m^(i-1)
  X = ind2coord(xi,dims);
  hX = etaijX0(P,m,n,i,j,X);
  if hX > h
    h = hX;
    xim = xi;
  end
end
return



function h = etaijX0(P,m,n,i,j,X)
h = 0;
for y1 = 1:m
  for y2 = y1+1:m
    P0 = getcondP(P,m,n,[j:n],[1:i],[X y1]);
    P1 = getcondP(P,m,n,[j:n],[1:i],[X y2]);
    h12 = .5*sum(abs(P0-P1));
    h = max(h,h12);
  end
end
return