function dpath = dpath(T,i,j)
% bounds hh(i,j) for a tree

% find j0, where
% j0 is the minimum descendant of i s.t. j0>=j


Ti = getSubt(T,i);

jj = Ti(Ti>=j);
if isempty(jj)
  dpath = [];
else
  j0 = min(jj);
  dpath = parentseq(T,j0,i);
end
return

function pseq = parentseq(T,u,r)
v = u;
el = 0;
pseq = [u];
while v~=r
  el = el+1;
  rents = find(T(:,v));
  v = rents;
  pseq = [v pseq];
end
return



function sT = getSubt(T,v)
% a list of v's descendents in T
% includes v itself
global SUBT
SUBT = [];
Pi(T,v);
sT = SUBT;
return

function Pi(T,v)
global SUBT
SUBT(end+1) = v;
vkids = find(T(v,:));
for vi=1:length(vkids)
  v1 = vkids(vi);
  Pi(T,v1);
end