drawtree(parent,parent2,children,children2,t,nDesc)

function [y,te,ye] = drawtree(parent,parent2,children,children2,t,nDesc)
% Draw a family tree representation.
% Parents is a 1*n vector giving the index of the unique parent of each item (zero if no parent)
% children is a 1*n cell array containing the children of each item (empty if no children)
% t is a n*1 time vector, used as the x coordinate

t0=min(t);
t=t-t0;
n=numel(t);

% The problem is to find the y coordinate for every child to result in a
% readable layout.

% Strategy 1: every child has y > its parent. 

% Strategy 2: define a envelope that we need to stay above. 
% we want y(i) at time t(i) to be > ye(max(find(t(i)>te)))

te = [0 7];
ye = [ 0 0];

y= zeros(1,n);
y(1)=0;
done=false(1,n);
% for i=2:n
%    if done(i)
%       continue;
%    end
%    iFrontier = min(find(te > t(i)));
%    y(i) = ye(iFrontier) + 1;
%    ye = [ye(1:max(iFrontier-1,1)) y(i) ye(min(iFrontier+1,numel(ye)):end)];
%    te = [te(1:max(iFrontier-1,1)) t(i) te(min(iFrontier+1,numel(te)):end)];
%    done(i)=true;
% end

[dsort,isort] = sort(nDesc,'descend');
for i = isort
   if parent2(i)==0
      [y,te,ye] = drawtree_slave(children2, parent2, t,y, te,ye, i);
   end
end

% figure(1);
% plot(t,y,'.');
% for i=1:n;
%     if parent(i); 
%         line([t(parent(i)) t(i)],[y(parent(i)) y(i)],'color','black');
%     elseif parent2(i);
%         line([t(parent2(i)) t(i)],[y(parent2(i)) y(i)],'color','black','linestyle','--');
%     end
% end
% line(te,ye,'color','red');
% 
% figure(2);
% plot(1:n,y,'.');
% for i=1:n; if parent(i); line([parent(i) i],[y(parent(i)) y(i)],'color','black'); end; end

return





function [y,te,ye] = drawtree_slave(children, parent, t,y, te,ye, i)

maxslope = 0.01*numel(parent)/(max(t)-min(t));

%if y(i)<0
if ~parent(i)   % i has no parent
    if t(i)<min(t(y>0))
        y(i)=max(y)/2;
    else
        y(i) = interp1(te,ye,t(i))+1;
    end
else   % i is a child, and we want its parent-child link to pass above all previously plotted nodes
    % first: try simply placing it above the envelope
    y(i) = interp1(te,ye,t(i))+1;
    i1 = min(find(te>t(parent(i))));  % first envelope point to right of parent
    i2 = max(find(te<t(i)));   % last envelope point to left of node i
    % if there are no envelope points between i and parent, nothing to do. 
    % otherwise...
    if i2>=i1
        ycheck = interp1([t(parent(i)) t(i)],[y(parent(i)) y(i)],te(i1:i2));
        [dmax,imax] = max(ye(i1:i2)-ycheck);
        if dmax>0    % a envelope point is above the line joining i to its parent
            imax = i1 + imax - 1;
            slope = (ye(imax) + 1 - y(parent(i))) / (te(imax)-t(parent(i)));
            if slope>maxslope
                % retrospective change to the parent of i
                y(parent(i)) = ye(imax)+1-maxslope*(te(imax)-t(parent(i)));
                % may need to update the envelope:
                ye(te==t(parent(i))) = y(parent(i));
                slope=maxslope;
            end
%            y(i) = y(parent(i)) + slope * (t(i)-t(parent(i)));
        end
    end
 %   line([t(i) t(parent(i))],[y(i) y(parent(i))],'color','black');
end

% update the envelope
if ~isempty(children{i}) | parent(i)
    if any (te==t(i))
        i1 = find(te==t(i));
        ye(i1) = y(i);
    else
        if parent(i)
            i1 = min(find(te>t(parent(i))));  % first envelope point to right of parent
            i2 = max(find(te<t(i)));   % last envelope point to left of node i
            % if there are envelope points between i and its parent, they may need to get deleted.
            if i2>=i1
                ycheck = interp1([t(parent(i)) t(i)],[y(parent(i)) y(i)],te(i1:i2));
                ilow = find(ye(i1:i2)<ycheck) + i1 -1 ;
                te(ilow)=[]; ye(ilow)=[];
            end
        end
        % add i to the envelope list
        i1 = max(find(te<t(i)));
        ye = [ye(1:max(i1,1)) y(i) ye(min(i1+1,numel(ye)):end)];
        te = [te(1:max(i1,1)) t(i) te(min(i1+1,numel(te)):end)];
    end
end
%line(t(i),y(i),'marker','.');


% does y have children?
children1 = children{i};
while ~isempty(children1) 
   % do the first unprocessed child of i:
   [y,te,ye] = drawtree_slave(children, parent, t,y, te,ye, children1(1));
   children1(1)=[];   
end

return