]> git.donarmstrong.com Git - lilypond.git/blobdiff - lily/skyline.cc
Imported Upstream version 2.16.0
[lilypond.git] / lily / skyline.cc
index b729e5312b15a56471f9f7b58ca542e0db947281..0250fc07f4743ba44e7bc1ac8af23ed6057aa435 100644 (file)
@@ -1,12 +1,25 @@
-/* skyline.cc -- implement the Skyline class
+/*
+  This file is part of LilyPond, the GNU music typesetter.
+
+  Copyright (C) 2006--2012 Joe Neeman <joeneeman@gmail.com>
+
+  LilyPond is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 3 of the License, or
+  (at your option) any later version.
 
-   source file of the GNU LilyPond music typesetter
-   (c) 2006--2007 Joe Neeman <joeneeman@gmail.com>
+  LilyPond is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with LilyPond.  If not, see <http://www.gnu.org/licenses/>.
 */
 
 #include "skyline.hh"
 #include <deque>
+#include <cstdio>
 
 #include "ly-smobs.icc"
 
    but the distance routine does.
 */
 
-/*
-  FIXME:
-
-  * Consider to use
-
-  typedef list<Skyline_point> Skyline;
-  struct Skyline_point
-  {
-    Real x;
-    Drul_array<Real> ys; 
-  };
-
-  this is a cleaner representation, as it doesn't duplicate the X, and
-  doesn't need bogus buildings at infinity  --hwn.
-
-
-  * All the messing around with EPS is very fishy.  There are no
-  complicated numerical algorithms involved, so EPS should not be
-  necessary.
-
-  --hwn
-  
-  
- */
-
-#define EPS 1e-10
-
-static inline bool
-approx_equal (Real x, Real y)
-{
-  return abs (x - y) < EPS || (isinf (x) && isinf (y) && ((x > 0) == (y > 0)));
-}
-
-static inline bool
-approx_greater_than (Real x, Real y)
-{
-  return x > y + EPS;
-}
-
-static inline bool
-approx_less_than (Real x, Real y)
-{
-  return x < y - EPS;
-}
-
-static inline bool
-approx_less_equal (Real x, Real y)
-{
-  return x <= y + EPS;
-}
+/* If we start including very thin buildings, numerical accuracy errors can
+   arise. Therefore, we ignore all buildings that are less than epsilon wide. */
+#define EPS 1e-5
 
-static inline bool
-approx_greater_equal (Real x, Real y)
+static void
+print_buildings (list<Building> const &b)
 {
-  return x >= y - EPS;
+  for (list<Building>::const_iterator i = b.begin (); i != b.end (); i++)
+    i->print ();
 }
 
 void
 Skyline::print () const
 {
-  for (list<Building>::const_iterator i = buildings_.begin ();
-       i != buildings_.end (); i++)
-    {
-      (*i).print ();
-    }
+  print_buildings (buildings_);
 }
 
-bool
-is_legal_skyline (list<Building> const &buildings)
+void
+Skyline::print_points () const
 {
-  list<Building>::const_iterator i;
-  Real last_x = -infinity_f;
-  for (i = buildings.begin (); i != buildings.end (); i++)
-    {
-      if (i->iv_[LEFT] != last_x)
-       return false;
-      last_x = i->iv_[RIGHT];
-      if (isinf (i->iv_.length ()) && i->height_[LEFT] != i->height_[RIGHT])
-       return false;
-    }
-  return last_x == infinity_f;
+  vector<Offset> ps (to_points (X_AXIS));
+
+  for (vsize i = 0; i < ps.size (); i++)
+    printf ("(%f,%f)%s", ps[i][X_AXIS], ps[i][Y_AXIS],
+            (i % 2) == 1 ? "\n" : " ");
 }
 
 Building::Building (Real start, Real start_height, Real end_height, Real end)
-  : iv_ (start, end)
 {
-  height_[LEFT] = start_height;
-  height_[RIGHT] = end_height;
-
   if (isinf (start) || isinf (end))
     assert (start_height == end_height);
 
-  precompute ();
+  end_ = end;
+  precompute (start, start_height, end_height, end);
 }
 
 Building::Building (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
 {
+  Real start = b[horizon_axis][LEFT] - horizon_padding;
+  Real end = b[horizon_axis][RIGHT] + horizon_padding;
   Real height = sky * b[other_axis (horizon_axis)][sky];
 
-  iv_ = b[horizon_axis];
-  iv_.widen (horizon_padding + EPS);
-  height_[LEFT] = height;
-  height_[RIGHT] = height;
-
-  if (sane ())
-    precompute ();
+  end_ = end;
+  precompute (start, height, height, end);
 }
 
 void
-Building::precompute ()
+Building::precompute (Real start, Real start_height, Real end_height, Real end)
 {
-  slope_ = (height_[RIGHT] - height_[LEFT]) / (iv_.length());
-  if (height_[LEFT] == height_[RIGHT]) /* in case they're both infinity */
+  slope_ = (end_height - start_height) / (end - start);
+  if (start_height == end_height) /* if they were both infinite, we would get nan, not 0, from the prev line */
     slope_ = 0;
 
   assert (!isinf (slope_) && !isnan (slope_));
 
-  if (isinf (iv_[START]))
+  if (isinf (start))
     {
-      assert (slope_ == 0);
-      y_intercept_ = height_[LEFT];
+      assert (start_height == end_height);
+      y_intercept_ = start_height;
     }
   else
-    y_intercept_ = height_[LEFT] - slope_ * iv_[START];
+    y_intercept_ = start_height - slope_ * start;
 }
 
-Real 
+Real
 Building::height (Real x) const
 {
-  if (isinf (x))
-    return (x > 0) ? height_[RIGHT] : height_[LEFT];
-  return slope_*x + y_intercept_;
+  return isinf (x) ? y_intercept_ : slope_ * x + y_intercept_;
 }
 
 void
 Building::print () const
 {
-  printf ("X[%f,%f] -> Y[%f,%f]\n",
-         iv_[LEFT], iv_[RIGHT],
-         height_[LEFT], height_[RIGHT]);
+  printf ("%f x + %f ends at %f\n", slope_, y_intercept_, end_);
 }
 
 Real
 Building::intersection_x (Building const &other) const
 {
-  return (y_intercept_ - other.y_intercept_) / (other.slope_ - slope_);
+  Real ret = (y_intercept_ - other.y_intercept_) / (other.slope_ - slope_);
+  return isnan (ret) ? -infinity_f : ret;
 }
 
 void
 Building::leading_part (Real chop)
 {
-  assert (chop > iv_[LEFT] && chop <= iv_[RIGHT] && !approx_equal (chop, iv_[LEFT]));
-  iv_[RIGHT] = chop;
-  height_[RIGHT] = height (chop);
+  assert (chop <= end_);
+  end_ = chop;
 }
 
 Building
-Building::sloped_neighbour (Real horizon_padding, Direction d) const
+Building::sloped_neighbour (Real start, Real horizon_padding, Direction d) const
 {
-  Real left = iv_[d];
-  Real right = iv_[d] + d * horizon_padding;
-  Real left_height = height_[d];
-  Real right_height = height_[d] - horizon_padding;
+  Real x = (d == LEFT) ? start : end_;
+  Real left = x;
+  Real right = x + d * horizon_padding;
+  Real left_height = height (x);
+  Real right_height = left_height - horizon_padding;
   if (d == LEFT)
     {
       swap (left, right);
@@ -213,88 +161,74 @@ Building::sloped_neighbour (Real horizon_padding, Direction d) const
   return Building (left, left_height, right_height, right);
 }
 
-bool
-Building::sane () const
+static Real
+first_intersection (Building const &b, list<Building> *const s, Real start_x)
 {
-  return approx_less_than (iv_[LEFT], iv_[RIGHT])
-    && !isinf (height_[RIGHT])
-    && !isinf (height_[LEFT]);
-}
+  while (!s->empty () && start_x < b.end_)
+    {
+      Building c = s->front ();
+      if (c.conceals (b, start_x))
+        return start_x;
 
-static void
-skyline_trailing_part (list<Building> *sky, Real x)
-{
-  if (approx_equal (x, sky->front ().iv_[RIGHT]))
-    sky->pop_front ();
-  else
-    assert (x < sky->front ().iv_[RIGHT]);
+      Real i = b.intersection_x (c);
+      if (i > start_x && i <= b.end_ && i <= c.end_)
+        return i;
 
-  if (!sky->empty ())
-    {
-      sky->front ().iv_[LEFT] = x;
-      sky->front ().height_[LEFT] = sky->front ().height (x);
+      start_x = c.end_;
+      if (b.end_ > c.end_)
+        s->pop_front ();
     }
+  return b.end_;
 }
 
 bool
-Building::conceals_beginning (Building const &other) const
+Building::conceals (Building const &other, Real x) const
 {
-  bool w = false;
-  Real h = other.height (iv_[LEFT]);
-  if (approx_equal (height_[LEFT], h))
-    w = slope_ > other.slope_;    
-  else if (height_[LEFT] > h) 
-    w = true;
-  else 
-    w = false;
-
-  return w;
-}
+  if (slope_ == other.slope_)
+    return y_intercept_ > other.y_intercept_;
 
-bool
-Building::conceals (Building const &other) const
-{
-  assert (iv_[LEFT] <= other.iv_[LEFT]);
-  return (iv_[RIGHT] >= other.iv_[RIGHT])
-    && approx_greater_equal (height (other.iv_[LEFT]), other.height_[LEFT])
-    && approx_greater_equal (height (other.iv_[RIGHT]), other.height_[RIGHT]);
+  /* their slopes were not equal, so there is an intersection point */
+  Real i = intersection_x (other);
+  return (i <= x && slope_ > other.slope_)
+         || (i > x && slope_ < other.slope_);
 }
 
 void
 Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
-                                list<Building> *const result)
+                                 list<Building> *const result)
 {
+  if (s1->empty () || s2->empty ())
+    {
+      programming_error ("tried to merge an empty skyline");
+      return;
+    }
+
+  Real x = -infinity_f;
   while (!s1->empty ())
     {
-      if (s2->front ().conceals_beginning (s1->front ()))
-       swap (s1, s2);
+      if (s2->front ().conceals (s1->front (), x))
+        swap (s1, s2);
 
       Building b = s1->front ();
-      while (!s2->empty () && b.conceals (s2->front ()))
-       s2->pop_front ();
+      Real end = first_intersection (b, s2, x);
+
       if (s2->empty ())
-       {
-         result->push_front (b);
-         break;
-       }
-
-      /* s2 either intersects with b or it ends after b */
-      Real end = infinity_f;
-      Real s2_start_height = s2->front ().height_[LEFT];
-      Real s2_end_height = s2->front ().height_[RIGHT];
-      Real s1_start_height = b.height (s2->front ().iv_[LEFT]);
-      Real s1_end_height = b.height (s2->front ().iv_[RIGHT]);
-      if (approx_greater_than (s2_start_height, s1_start_height))
-       end = s2->front ().iv_[LEFT];
-      else if (approx_greater_than (s2_end_height, s1_end_height))
-       end = b.intersection_x (s2->front ());
-      end = min (end, b.iv_[RIGHT]);
-
-      b.leading_part (end);
-      result->push_front (b);
-
-      skyline_trailing_part (s1, end);
-      skyline_trailing_part (s2, end);
+        {
+          result->push_front (b);
+          break;
+        }
+
+      /* only include buildings wider than epsilon */
+      if (end > x + EPS)
+        {
+          b.leading_part (end);
+          result->push_front (b);
+        }
+
+      if (end >= s1->front ().end_)
+        s1->pop_front ();
+
+      x = end;
     }
   result->reverse ();
 }
@@ -305,63 +239,90 @@ empty_skyline (list<Building> *const ret)
   ret->push_front (Building (-infinity_f, -infinity_f, -infinity_f, infinity_f));
 }
 
+/*
+  Given Building 'b' with starting wall location 'start', extend each side
+  with a sloped roofline of width 'horizon_padding'; put the skyline in 'ret'
+*/
 static void
-single_skyline (Building b, Real horizon_padding, list<Building> *const ret)
-{
-  b.iv_.widen (horizon_padding);
-  
-  if (!isinf (b.iv_[RIGHT]))
-    ret->push_front (Building (b.iv_[RIGHT], -infinity_f,
-                              -infinity_f, infinity_f));
-  if (horizon_padding > 0 && !isinf (b.iv_.length ()))
-    ret->push_front (b.sloped_neighbour (horizon_padding, RIGHT));
-  
-  if (b.iv_[RIGHT] > b.iv_[LEFT])
+single_skyline (Building b, Real start, Real horizon_padding, list<Building> *const ret)
+{
+  bool sloped_neighbours = horizon_padding > 0 && !isinf (start) && !isinf (b.end_);
+  if (!isinf (b.end_))
+    ret->push_front (Building (b.end_ + horizon_padding, -infinity_f,
+                               -infinity_f, infinity_f));
+  if (sloped_neighbours)
+    ret->push_front (b.sloped_neighbour (start, horizon_padding, RIGHT));
+
+  if (b.end_ > start + EPS)
     ret->push_front (b);
 
-  if (horizon_padding > 0 && !isinf (b.iv_.length ()))
-    ret->push_front (b.sloped_neighbour (horizon_padding, LEFT));
-  if (!isinf (b.iv_[LEFT]))
+  if (sloped_neighbours)
+    ret->push_front (b.sloped_neighbour (start, horizon_padding, LEFT));
+
+  if (!isinf (start))
     ret->push_front (Building (-infinity_f, -infinity_f,
-                              -infinity_f, b.iv_[LEFT]));
+                               -infinity_f, start - horizon_padding));
 }
 
-/* remove a non-overlapping set of buildings from BUILDINGS and build a skyline
+/* remove a non-overlapping set of boxes from BOXES and build a skyline
    out of them */
 static list<Building>
-non_overlapping_skyline (list<Building> *const buildings)
+non_overlapping_skyline (list<Box> *const boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
 {
   list<Building> result;
   Real last_end = -infinity_f;
-  list<Building>::iterator i = buildings->begin ();
-  while (i != buildings->end ())
+  list<Box>::iterator i = boxes->begin ();
+  while (i != boxes->end ())
     {
-      if (approx_less_than (i->iv_[LEFT], last_end))
-       {
-         i++;
-         continue;
-       }
-
-      if (approx_greater_than (i->iv_[LEFT], last_end))
-       result.push_back (Building (last_end, -infinity_f, -infinity_f, i->iv_[LEFT]));
-      else
-       i->iv_[LEFT] = last_end;
-
-      last_end = i->iv_[RIGHT];
-      list<Building>::iterator j = i;
-      i++;
-      result.splice (result.end (), *buildings, j);
+      Interval iv = (*i)[horizon_axis];
+
+      if (iv[LEFT] - horizon_padding < last_end)
+        {
+          i++;
+          continue;
+        }
+
+      if (iv[LEFT] - horizon_padding > last_end + EPS)
+        result.push_front (Building (last_end, -infinity_f, -infinity_f, iv[LEFT] - 2 * horizon_padding));
+
+      Building b (*i, horizon_padding, horizon_axis, sky);
+      bool sloped_neighbours = horizon_padding > 0 && !isinf (iv.length ());
+      if (sloped_neighbours)
+        result.push_front (b.sloped_neighbour (iv[LEFT] - horizon_padding, horizon_padding, LEFT));
+      result.push_front (b);
+      if (sloped_neighbours)
+        result.push_front (b.sloped_neighbour (iv[LEFT] - horizon_padding, horizon_padding, RIGHT));
+
+      list<Box>::iterator j = i++;
+      boxes->erase (j);
+      last_end = result.front ().end_;
     }
   if (last_end < infinity_f)
-    result.push_back (Building (last_end, -infinity_f, -infinity_f, infinity_f));
-  assert (is_legal_skyline (result));
+    result.push_front (Building (last_end, -infinity_f, -infinity_f, infinity_f));
+  result.reverse ();
   return result;
 }
 
+class LessThanBox
+{
+  Axis a_;
+
+public:
+  LessThanBox (Axis a)
+  {
+    a_ = a;
+  }
+
+  bool operator () (Box const &b1, Box const &b2)
+  {
+    return b1[a_][LEFT] < b2[a_][LEFT];
+  }
+};
+
 list<Building>
-Skyline::internal_build_skyline (list<Building> *buildings)
+Skyline::internal_build_skyline (list<Box> *boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
 {
-  vsize size = buildings->size ();
+  vsize size = boxes->size ();
 
   if (size == 0)
     {
@@ -372,16 +333,18 @@ Skyline::internal_build_skyline (list<Building> *buildings)
   else if (size == 1)
     {
       list<Building> result;
-      single_skyline (buildings->front (), 0, &result);
+      single_skyline (Building (boxes->front (), horizon_padding, horizon_axis, sky),
+                      boxes->front ()[horizon_axis][LEFT] - horizon_padding,
+                      horizon_padding, &result);
       return result;
     }
 
   deque<list<Building> > partials;
-  buildings->sort ();
-  while (!buildings->empty ())
-    partials.push_back (non_overlapping_skyline (buildings));
+  boxes->sort (LessThanBox (horizon_axis));
+  while (!boxes->empty ())
+    partials.push_back (non_overlapping_skyline (boxes, horizon_padding, horizon_axis, sky));
 
-  /* we'd like to say while (partials->size () > 1) but that's O(n).
+  /* we'd like to say while (partials->size () > 1) but that's O (n).
      Instead, we exit in the middle of the loop */
   while (!partials.empty ())
     {
@@ -389,7 +352,7 @@ Skyline::internal_build_skyline (list<Building> *buildings)
       list<Building> one = partials.front ();
       partials.pop_front ();
       if (partials.empty ())
-       return one;
+        return one;
 
       list<Building> two = partials.front ();
       partials.pop_front ();
@@ -403,13 +366,14 @@ Skyline::internal_build_skyline (list<Building> *buildings)
 Skyline::Skyline ()
 {
   sky_ = UP;
-  empty_skyline (&buildings_);  
+  empty_skyline (&buildings_);
 }
 
 Skyline::Skyline (Skyline const &src)
 {
   sky_ = src.sky_;
+
+  /* doesn't a list's copy constructor do this? -- jneem */
   for (list<Building>::const_iterator i = src.buildings_.begin ();
        i != src.buildings_.end (); i++)
     {
@@ -423,6 +387,37 @@ Skyline::Skyline (Direction sky)
   empty_skyline (&buildings_);
 }
 
+/*
+  build padded skyline from an existing skyline with padding
+  added to it.
+*/
+
+Skyline::Skyline (Skyline const &src, Real horizon_padding, Axis /*a*/)
+{
+  /*
+     We extract boxes from the skyline, then build a new skyline from
+     the boxes.
+     A box is created for every horizontal portion of the skyline
+     Because skylines are defined positive, and then inverted if they
+     are to be down-facing, we create the new skyline in the UP
+     direction, then give it the down direction if needed.
+  */
+  Real start = -infinity_f;
+  list<Box> boxes;
+
+  // establish a baseline box
+  // FIXME: This has hardcoded logic, assuming a == X_AXIS!
+  boxes.push_back (Box (Interval (-infinity_f, infinity_f),
+                        Interval (0, 0)));
+  list<Building>::const_iterator end = src.buildings_.end ();
+  for (list<Building>::const_iterator i = src.buildings_.begin (); i != end; start = i->end_, i++)
+    if ((i->slope_ == 0) && !isinf (i->y_intercept_))
+      boxes.push_back (Box (Interval (start, i->end_),
+                            Interval (-infinity_f, i->y_intercept_)));
+  buildings_ = internal_build_skyline (&boxes, horizon_padding, X_AXIS, UP);
+  sky_ = src.sky_;
+}
+
 /*
   build skyline from a set of boxes. If horizon_padding > 0, expand all the boxes
   by that amount and add 45-degree sloped boxes to the edges of each box (of
@@ -434,32 +429,27 @@ Skyline::Skyline (Direction sky)
  */
 Skyline::Skyline (vector<Box> const &boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
 {
-  list<Building> bldgs;
+  list<Box> filtered_boxes;
   sky_ = sky;
 
+  Axis vert_axis = other_axis (horizon_axis);
   for (vsize i = 0; i < boxes.size (); i++)
     {
-      Building front (boxes[i], horizon_padding, horizon_axis, sky);
-      if (front.sane ())
-       {
-         bldgs.push_front (front);
-         if (horizon_padding > 0 && !isinf (front.iv_.length ()))
-           {
-             bldgs.push_front (front.sloped_neighbour (horizon_padding, LEFT));
-             bldgs.push_front (front.sloped_neighbour (horizon_padding, RIGHT));
-           }
-       }
+      Interval iv = boxes[i][horizon_axis];
+      iv.widen (horizon_padding);
+      if (iv.length () > EPS && !boxes[i][vert_axis].is_empty ())
+        filtered_boxes.push_front (boxes[i]);
     }
-  
-  buildings_ = internal_build_skyline (&bldgs);
-  assert (is_legal_skyline (buildings_));
+
+  buildings_ = internal_build_skyline (&filtered_boxes, horizon_padding, horizon_axis, sky);
 }
 
 Skyline::Skyline (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
 {
   sky_ = sky;
-  Building front (b, 0, horizon_axis, sky);
-  single_skyline (front, horizon_padding, &buildings_);
+  Building front (b, horizon_padding, horizon_axis, sky);
+  single_skyline (front, b[horizon_axis][LEFT] - horizon_padding,
+                  horizon_padding, &buildings_);
 }
 
 void
@@ -471,7 +461,6 @@ Skyline::merge (Skyline const &other)
   list<Building> my_bld;
   my_bld.splice (my_bld.begin (), buildings_);
   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
-  assert (is_legal_skyline (buildings_));
 }
 
 void
@@ -480,10 +469,23 @@ Skyline::insert (Box const &b, Real horizon_padding, Axis a)
   list<Building> other_bld;
   list<Building> my_bld;
 
+  if (isnan (b[other_axis (a)][LEFT])
+      || isnan (b[other_axis (a)][RIGHT]))
+    {
+      programming_error ("insane box for skyline");
+      return;
+    }
+
+  /* do the same filtering as in Skyline (vector<Box> const&, etc.) */
+  Interval iv = b[a];
+  iv.widen (horizon_padding);
+  if (iv.length () <= EPS || b[other_axis (a)].is_empty ())
+    return;
+
   my_bld.splice (my_bld.begin (), buildings_);
-  single_skyline (Building (b, 0, a, sky_), horizon_padding, &other_bld);
+  single_skyline (Building (b, horizon_padding, a, sky_), b[a][LEFT] - horizon_padding,
+                  horizon_padding, &other_bld);
   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
-  assert (is_legal_skyline (buildings_));
 }
 
 void
@@ -491,43 +493,89 @@ Skyline::raise (Real r)
 {
   list<Building>::iterator end = buildings_.end ();
   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
-    {
-      i->height_[LEFT] += sky_ * r;
-      i->height_[RIGHT] += sky_ * r;
-      i->y_intercept_ += sky_ * r;
-    }
-  assert (is_legal_skyline (buildings_));
+    i->y_intercept_ += sky_ * r;
 }
 
 void
-Skyline::shift (Real r)
+Skyline::shift (Real s)
 {
   list<Building>::iterator end = buildings_.end ();
   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
     {
-      i->iv_[LEFT] += r;
-      i->iv_[RIGHT] += r;
+      i->end_ += s;
+      i->y_intercept_ -= s * i->slope_;
     }
 }
 
 Real
-Skyline::distance (Skyline const &other) const
+Skyline::distance (Skyline const &other, Real horizon_padding) const
+{
+  Real dummy;
+  return internal_distance (other, horizon_padding, &dummy);
+}
+
+Real
+Skyline::touching_point (Skyline const &other, Real horizon_padding) const
+{
+  Real touch;
+  internal_distance (other, horizon_padding, &touch);
+  return touch;
+}
+
+Real
+Skyline::internal_distance (Skyline const &other, Real horizon_padding, Real *touch_point) const
 {
   assert (sky_ == -other.sky_);
-  list<Building>::const_iterator i = buildings_.begin ();
-  list<Building>::const_iterator j = other.buildings_.begin ();
+
+  Skyline const *padded_this = this;
+  Skyline const *padded_other = &other;
+  bool created_tmp_skylines = false;
+
+  /*
+    For systems, padding is not added at creation time.  Padding is
+    added to AxisGroup objects when outside-staff objects are added.
+    Thus, when we want to place systems with horizontal padding,
+    we do it at distance calculation time.
+  */
+  if (horizon_padding != 0.0)
+    {
+      padded_this = new Skyline (*padded_this, horizon_padding, X_AXIS);
+      padded_other = new Skyline (*padded_other, horizon_padding, X_AXIS);
+      created_tmp_skylines = true;
+    }
+
+  list<Building>::const_iterator i = padded_this->buildings_.begin ();
+  list<Building>::const_iterator j = padded_other->buildings_.begin ();
 
   Real dist = -infinity_f;
-  while (i != buildings_.end () && j != other.buildings_.end ())
+  Real start = -infinity_f;
+  Real touch = -infinity_f;
+  while (i != padded_this->buildings_.end () && j != padded_other->buildings_.end ())
     {
-      Interval iv = intersection (i->iv_, j->iv_);
-      dist = max (dist, max (i->height (iv[LEFT]) + j->height (iv[LEFT]),
-                            i->height (iv[RIGHT]) + j->height (iv[RIGHT])));
-      if (i->iv_[RIGHT] <= j->iv_[RIGHT])
-       i++;
+      Real end = min (i->end_, j->end_);
+      Real start_dist = i->height (start) + j->height (start);
+      Real end_dist = i->height (end) + j->height (end);
+      dist = max (dist, max (start_dist, end_dist));
+
+      if (end_dist == dist)
+        touch = end;
+      else if (start_dist == dist)
+        touch = start;
+
+      if (i->end_ <= j->end_)
+        i++;
       else
-       j++;
+        j++;
+      start = end;
     }
+
+  if (created_tmp_skylines)
+    {
+      delete padded_this;
+      delete padded_other;
+    }
+
+  *touch_point = touch;
   return dist;
 }
 
@@ -539,8 +587,8 @@ Skyline::height (Real airplane) const
   list<Building>::const_iterator i;
   for (i = buildings_.begin (); i != buildings_.end (); i++)
     {
-      if (i->iv_[RIGHT] >= airplane)
-       return sky_ * i->height (airplane);
+      if (i->end_ >= airplane)
+        return sky_ * i->height (airplane);
     }
 
   assert (0);
@@ -555,103 +603,67 @@ Skyline::max_height () const
   return sky_ * distance (s);
 }
 
+Real
+Skyline::max_height_position () const
+{
+  Skyline s (-sky_);
+  s.set_minimum_height (0);
+  return touching_point (s);
+}
+
 void
 Skyline::set_minimum_height (Real h)
 {
   Skyline s (sky_);
-  s.buildings_.front ().height_[LEFT] = h * sky_;
-  s.buildings_.front ().height_[RIGHT] = h * sky_;
   s.buildings_.front ().y_intercept_ = h * sky_;
   merge (s);
 }
 
-
 vector<Offset>
-Skyline::to_points () const
+Skyline::to_points (Axis horizon_axis) const
 {
   vector<Offset> out;
 
+  Real start = -infinity_f;
   for (list<Building>::const_iterator i (buildings_.begin ());
        i != buildings_.end (); i++)
     {
-      if (!isinf (i->iv_[LEFT]) && !isinf (i->height_[LEFT]))
-       out.push_back (Offset (i->iv_[LEFT], sky_ * i->height_[LEFT]));
-      if (!isinf (i->iv_[RIGHT]) && !isinf (i->height_[RIGHT]))
-       out.push_back (Offset (i->iv_[RIGHT], sky_ * i->height_[RIGHT]));
+      out.push_back (Offset (start, sky_ * i->height (start)));
+      out.push_back (Offset (i->end_, sky_ * i->height (i->end_)));
+      start = i->end_;
     }
-  return out;
-}
 
-Skyline_pair::Skyline_pair ()
-  : skylines_ (Skyline (DOWN), Skyline (UP))
-{
-}
-
-Skyline_pair::Skyline_pair (vector<Box> const &boxes, Real padding, Axis a)
-  : skylines_ (Skyline (boxes, padding, a, DOWN), Skyline (boxes, padding, a, UP))
-{
-}
-
-Skyline_pair::Skyline_pair (Box const &b, Real padding, Axis a)
-  : skylines_ (Skyline (b, padding, a, DOWN), Skyline (b, padding, a, UP))
-{
-}
-
-void
-Skyline_pair::raise (Real r)
-{
-  skylines_[UP].raise (r);
-  skylines_[DOWN].raise (r);
-}
+  if (horizon_axis == Y_AXIS)
+    for (vsize i = 0; i < out.size (); i++)
+      out[i] = out[i].swapped ();
 
-void
-Skyline_pair::shift (Real r)
-{
-  skylines_[UP].shift (r);
-  skylines_[DOWN].shift (r);
+  return out;
 }
 
-void
-Skyline_pair::insert (Box const &b, Real padding, Axis a)
+bool
+Skyline::is_empty () const
 {
-  skylines_[UP].insert (b, padding, a);
-  skylines_[DOWN].insert (b, padding, a);
+  Building b = buildings_.front ();
+  return b.end_ == infinity_f && b.y_intercept_ == -infinity_f;
 }
 
 void
-Skyline_pair::merge (Skyline_pair const &other)
-{
-  skylines_[UP].merge (other[UP]);
-  skylines_[DOWN].merge (other[DOWN]);
-}
-
-Skyline&
-Skyline_pair::operator [] (Direction d)
-{
-  return skylines_[d];
-}
-
-Skyline const&
-Skyline_pair::operator [] (Direction d) const
+Skyline::clear ()
 {
-  return skylines_[d];
+  buildings_.clear ();
+  empty_skyline (&buildings_);
 }
 
 /****************************************************************/
 
-
 IMPLEMENT_SIMPLE_SMOBS (Skyline);
 IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
 IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
 
-IMPLEMENT_SIMPLE_SMOBS (Skyline_pair);
-IMPLEMENT_TYPE_P (Skyline_pair, "ly:skyline-pair?");
-IMPLEMENT_DEFAULT_EQUAL_P (Skyline_pair);
-
 SCM
-Skyline::mark_smob (SCM)
+Skyline::mark_smob (SCM s)
 {
-  ASSERT_LIVE_IS_ALLOWED();
+  ASSERT_LIVE_IS_ALLOWED (s);
   return SCM_EOL;
 }
 
@@ -666,18 +678,60 @@ Skyline::print_smob (SCM s, SCM port, scm_print_state *)
   return 1;
 }
 
+MAKE_SCHEME_CALLBACK_WITH_OPTARGS (Skyline, get_touching_point, 3, 1, "")
 SCM
-Skyline_pair::mark_smob (SCM)
+Skyline::get_touching_point (SCM skyline_scm, SCM other_skyline_scm, SCM horizon_padding_scm)
 {
-  return SCM_EOL;
+  LY_ASSERT_SMOB (Skyline, other_skyline_scm, 1);
+
+  Real horizon_padding = 0;
+  if (horizon_padding_scm != SCM_UNDEFINED)
+    {
+      LY_ASSERT_TYPE (scm_is_number, horizon_padding_scm, 3);
+      horizon_padding = scm_to_double (horizon_padding_scm);
+    }
+
+  Skyline *skyline = Skyline::unsmob (skyline_scm);
+  Skyline *other_skyline = Skyline::unsmob (other_skyline_scm);
+  return scm_from_double (skyline->touching_point (*other_skyline, horizon_padding));
 }
 
-int
-Skyline_pair::print_smob (SCM s, SCM port, scm_print_state *)
+MAKE_SCHEME_CALLBACK_WITH_OPTARGS (Skyline, get_distance, 3, 1, "")
+SCM
+Skyline::get_distance (SCM skyline_scm, SCM other_skyline_scm, SCM horizon_padding_scm)
 {
-  Skyline_pair *r = (Skyline_pair *) SCM_CELL_WORD_1 (s);
-  (void) r;
+  LY_ASSERT_SMOB (Skyline, other_skyline_scm, 1);
 
-  scm_puts ("#<Skyline-pair>", port);
-  return 1;
+  Real horizon_padding = 0;
+  if (horizon_padding_scm != SCM_UNDEFINED)
+    {
+      LY_ASSERT_TYPE (scm_is_number, horizon_padding_scm, 3);
+      horizon_padding = scm_to_double (horizon_padding_scm);
+    }
+
+  Skyline *skyline = Skyline::unsmob (skyline_scm);
+  Skyline *other_skyline = Skyline::unsmob (other_skyline_scm);
+  return scm_from_double (skyline->distance (*other_skyline, horizon_padding));
+}
+
+MAKE_SCHEME_CALLBACK (Skyline, get_max_height, 1)
+SCM
+Skyline::get_max_height (SCM skyline_scm)
+{
+  return scm_from_double (Skyline::unsmob (skyline_scm)->max_height ());
+}
+
+MAKE_SCHEME_CALLBACK (Skyline, get_max_height_position, 1)
+SCM
+Skyline::get_max_height_position (SCM skyline_scm)
+{
+  return scm_from_double (Skyline::unsmob (skyline_scm)->max_height_position ());
+}
+
+MAKE_SCHEME_CALLBACK (Skyline, get_height, 2)
+SCM
+Skyline::get_height (SCM skyline_scm, SCM x_scm)
+{
+  Real x = robust_scm2double (x_scm, 0.0);
+  return scm_from_double (Skyline::unsmob (skyline_scm)->height (x));
 }