sclass PtTree {
  Node root = new Leaf;
  
  abstract class Node {
    Node parent;
    
    abstract bool add(Pt p);
    abstract void collectPointsIn(Rect r, L<Pt> out);
    
    abstract void replaceChild(Node from, Node to);
    
    void addAll(Iterable<Pt> l) {
      fOr (Pt p : l) add(p);
    }
  }
  
  class Leaf extends Node {
    Cl<Pt> points;
    
    *(Node *parent) {}
    
    bool add(Pt p) {
      if (contains(points, p)) false;
      if (l(points) >= maxPointsPerNode)
        ret split().add(p);
      if (points == null) points = new L;
      points.add(p);
      true;
    }
    
    void collectPointsIn(Rect r, L<Pt> out) {
      addAll(out, points);
    }
    
    Split split() {
      int n = l(points), half = n/2;
      
      IntPair xSplit = bestXSplitPoint();
      IntPair ySplit = bestXSplitPoint();
      
      Split split;
      IPred<Pt> pred;
      if (abs(xSplit.b-half) < abs(ySplit.b-half)) {
        split = new SplitX;
        split.splitPoint = xSplit.a;
        pred = p -> p.x >= split.splitPoint;
      } else {
        split = new SplitY;
        split.splitPoint = ySplit.a;
        pred = p -> p.y >= split.splitPoint;
      }
      
      Leaf a = new Leaf(split);
      a.addAll(antiFilter pred(points));
      split.a = a;
      Leaf b = new Leaf(split);
      b.addAll(filter pred(points));
      split.b = b;
      ret replaceMeWith(split);
    }
    
    Node replaceMeWith(Node node) {
      node.parent = parent;
      if (parent != null) parent.replaceChild(this, node);
      else if (root == this) root = node;
      ret node;
    }

    // (splitPoint, count in upper set)
    IntPair bestXSplitPoint() {
      L<Pt> lx = sortedBy(points, p -> p.x);
      
      int n = l(points), half = n/2;
      int i = 0;
      int splitPoint = Int.MIN_VALUE;
      while true {
        int lastSplitpoint = splitPoint;
        splitPoint = lx.get(i++).x+1;
        int lastI = i;
        while (i < n && lx.get(i).x < ySplitPoint) ++i;
        if (i >= half)
          ret abs(lastI-half) < abs(i-half) ? pair(lastSplitPoint, lastI) : pair(splitPoint, i);
      }
    }
    
    // (splitPoint, count in upper set)
    IntPair bestYSplitPoint() {
      L<Pt> ly = sortedBy(points, p -> p.y);
      
      int n = l(points), half = n/2;
      int i = 0;
      int splitPoint = Int.MIN_VALUE;
      while true {
        int lastSplitpoint = splitPoint;
        splitPoint = ly.get(i++).y+1;
        int lastI = i;
        while (i < n && ly.get(i).y < ySplitPoint) ++i;
        if (i >= half)
          ret abs(lastI-half) < abs(i-half) ? pair(lastSplitPoint, lastI) : pair(splitPoint, i);
      }
    }
    
    void replaceChild(Node from, Node to) { unimplemented(); } // doesn't need
  }
  
  abstract class Split extends Node {
    int splitPoint;
    Node a, b;
    
    void collectPointsIn(Rect r, L<Pt> out) {
      a.collectPointsIn(r, out);
      b.collectPointsIn(r, out);
    }
    
    void replaceChild(Node from, Node to) {
      if (a == from) a = to;
      if (b == from) b = to;
    } 
  }
  
  class SplitH extends Split {
    bool add(Pt p) {
      ret (p.x >= splitPoint ? b : a).add(p);
    }
  }
  
  class SplitV extends Split {
    bool add(Pt p) {
      ret (p.y >= splitPoint ? b : a).add(p);
    }
  }
  
  int maxPointsPerNode = 4;
  
  bool add(Pt p) {
    ret root.add(p);
  }
  
  L<Pt> pointsIn(Rect r) {
    new L<Pt> out;
    root.collectPointsIn(r, out);
    ret out;
  }
}