%
% Adapted from a code written by Andrew Christlieb.
%
  function [tree] = m671b_build_tree(particles)
%
% input - particles : a row vector containing the particle locations
% output - tree: an array of nodes comprising the tree
%
% build root node
%
  x = [];
  tree(1) = struct( 'members' , x , 'interval' , x , 'children' , x );
  N = length(particles);    % assume N >= 3
  tree(1).members = 1:N;    % indices of particles belonging to node
  tree(1).interval = [0,1]; % interval defining the node
  tree(1).children = [2,3]; % indices in tree array of children
%
% build child nodes
%
  child(1) = struct( 'members' , x , 'interval' , x , 'children' , x );
  child(2) = struct( 'members' , x , 'interval' , x , 'children' , x );
  i = 0;
  node_count = 1;
%
  while i < length(tree)
  i = i + 1;
  midpoint = ( tree(i).interval(1) + tree(i).interval(2) )/2;
  child(1).interval = [tree(i).interval(1),midpoint];
  child(2).interval = [midpoint,tree(i).interval(2)];
  child(1).members = []; child(2).members = []; % clear temp arrays
  n = length(tree(i).members);
  if n >= 2
    count(1) = 1; count(2) = 1;
    for j = 1:n
      index = tree(i).members(j);
      if particles(index) <= midpoint
        child(1).members(count(1)) = index; count(1) = count(1) + 1;
      else; 
        child(2).members(count(2)) = index; count(2) = count(2) + 1; 
      end
    end
    for j = 1:2
    if count(j) >= 2;
      node_count = node_count + 1;
      tree(i).children(j) = node_count;
      tree(node_count) = child(j);
    end
    end
  end
  end