]> git.donarmstrong.com Git - lilypond.git/blobdiff - lily/skyline.cc
lilypond-manuals.css: edit color scheme and some spacing
[lilypond.git] / lily / skyline.cc
index 8f62710b0947fa8e48ffcbc7e181cba4709c96de..6ef2069c39b50b0acb7273c06928c3e509cf4729 100644 (file)
@@ -1,7 +1,7 @@
 /*
   This file is part of LilyPond, the GNU music typesetter.
 
-  Copyright (C) 2006--2011 Joe Neeman <joeneeman@gmail.com>
+  Copyright (C) 2006--2015 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
 */
 
 #include "skyline.hh"
+#include "skyline-pair.hh"
 #include <deque>
 #include <cstdio>
 
-#include "ly-smobs.icc"
-
 /* A skyline is a sequence of non-overlapping buildings: something like
    this:
                    _______
    a DOWN skyline with an UP skyline, we need to flip the DOWN skyline first.
    This means that the merging routine doesn't need to be aware of direction,
    but the distance routine does.
-*/
 
-/* 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
+   From 2007 through 2012, buildings of width less than EPS were discarded,
+   citing numerical accuracy concerns.  We remember that floating point
+   comparisons of nearly-equal values can be affected by rounding error.
+   Also, some target machines use the x87 floating point unit, which provides
+   extended precision for intermediate results held in registers. On this type
+   of hardware comparisons such as
+     double c = 1.0/3.0; boolean compare = (c == 1.0/3.0)
+   could go either way because the 1.0/3.0 is allowed to be kept
+   higher precision than the variable 'c'.
+   Alert to these considerations, we now accept buildings of zero-width.
+*/
 
 static void
 print_buildings (list<Building> const &b)
@@ -87,16 +93,18 @@ Building::Building (Real start, Real start_height, Real end_height, Real end)
   if (isinf (start) || isinf (end))
     assert (start_height == end_height);
 
+  start_ = start;
   end_ = end;
   precompute (start, start_height, end_height, end);
 }
 
-Building::Building (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
+Building::Building (Box const &b, Axis horizon_axis, Direction sky)
 {
-  Real start = b[horizon_axis][LEFT] - horizon_padding;
-  Real end = b[horizon_axis][RIGHT] + horizon_padding;
+  Real start = b[horizon_axis][LEFT];
+  Real end = b[horizon_axis][RIGHT];
   Real height = sky * b[other_axis (horizon_axis)][sky];
 
+  start_ = start;
   end_ = end;
   precompute (start, height, height, end);
 }
@@ -104,9 +112,9 @@ Building::Building (Box const &b, Real horizon_padding, Axis horizon_axis, Direc
 void
 Building::precompute (Real start, Real start_height, Real end_height, Real end)
 {
-  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;
+  slope_ = 0.0; /* if they were both infinite, we would get nan, not 0, from the prev line */
+  if (start_height != end_height)
+    slope_ = (end_height - start_height) / (end - start);
 
   assert (!isinf (slope_) && !isnan (slope_));
 
@@ -115,11 +123,17 @@ Building::precompute (Real start, Real start_height, Real end_height, Real end)
       assert (start_height == end_height);
       y_intercept_ = start_height;
     }
+  else if (fabs (slope_) > 1e6)
+    // too steep to be stored in slope-intercept form, given round-off error
+    {
+      slope_ = 0.0;
+      y_intercept_ = max (start_height, end_height);
+    }
   else
     y_intercept_ = start_height - slope_ * start;
 }
 
-Real
+inline Real
 Building::height (Real x) const
 {
   return isinf (x) ? y_intercept_ : slope_ * x + y_intercept_;
@@ -128,109 +142,138 @@ Building::height (Real x) const
 void
 Building::print () const
 {
-  printf ("%f x + %f ends at %f\n", slope_, y_intercept_, end_);
+  printf ("%f x + %f from %f to %f\n", slope_, y_intercept_, start_, end_);
 }
 
-Real
+inline Real
 Building::intersection_x (Building const &other) const
 {
   Real ret = (y_intercept_ - other.y_intercept_) / (other.slope_ - slope_);
   return isnan (ret) ? -infinity_f : ret;
 }
 
-void
-Building::leading_part (Real chop)
+// Returns a shift s such that (x + s, y) intersects the roof of
+// this building.  If no such shift exists, returns infinity_f.
+Real
+Building::shift_to_intersect (Real x, Real y) const
 {
-  assert (chop <= end_);
-  end_ = chop;
+  // Solve for s: y = (x + s)*m + b
+  Real ret = (y - y_intercept_ - slope_ * x) / slope_;
+
+  if (ret >= start_ && ret <= end_ && !isinf (ret))
+    return ret;
+  return infinity_f;
 }
 
-Building
-Building::sloped_neighbour (Real start, Real horizon_padding, Direction d) const
+bool
+Building::above (Building const &other, Real x) const
 {
-  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);
-      swap (left_height, right_height);
-    }
-  return Building (left, left_height, right_height, right);
+  return (isinf (y_intercept_) || isinf (other.y_intercept_) || isinf (x))
+         ? y_intercept_ > other.y_intercept_
+         : (slope_ - other.slope_) * x + y_intercept_ > other.y_intercept_;
 }
 
-static Real
-first_intersection (Building const &b, list<Building> *const s, Real start_x)
+// Remove redundant empty buildings from the skyline.
+// If there are two adjacent empty buildings, they can be
+// turned into one.
+void
+Skyline::normalize ()
 {
-  while (!s->empty () && start_x < b.end_)
-    {
-      Building c = s->front ();
-      if (c.conceals (b, start_x))
-        return start_x;
+  bool last_empty = false;
+  list<Building>::iterator i;
 
-      Real i = b.intersection_x (c);
-      if (i > start_x && i <= b.end_ && i <= c.end_)
-        return i;
-
-      start_x = c.end_;
-      if (b.end_ > c.end_)
-        s->pop_front ();
+  for (i = buildings_.begin (); i != buildings_.end (); i++)
+    {
+      if (last_empty && i->y_intercept_ == -infinity_f)
+        {
+          list<Building>::iterator last = i;
+          last--;
+          last->end_ = i->end_;
+          buildings_.erase (i);
+          i = last;
+        }
+      last_empty = (i->y_intercept_ == -infinity_f);
     }
-  return b.end_;
-}
-
-bool
-Building::conceals (Building const &other, Real x) const
-{
-  if (slope_ == other.slope_)
-    return y_intercept_ > other.y_intercept_;
 
-  /* 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_);
+  assert (buildings_.front ().start_ == -infinity_f);
+  assert (buildings_.back ().end_ == infinity_f);
 }
 
 void
-Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
-                                 list<Building> *const result)
+Skyline::internal_merge_skyline (list<Building> *sb, list<Building> *sc,
+                                 list<Building> *const result) const
 {
-  if (s1->empty () || s2->empty ())
+  if (sb->empty () || sc->empty ())
     {
       programming_error ("tried to merge an empty skyline");
       return;
     }
 
-  Real x = -infinity_f;
-  while (!s1->empty ())
+  Building b = sb->front ();
+  for (; !sc->empty (); sc->pop_front ())
     {
-      if (s2->front ().conceals (s1->front (), x))
-        swap (s1, s2);
-
-      Building b = s1->front ();
-      Real end = first_intersection (b, s2, x);
-
-      if (s2->empty ())
+      /* Building b is continuing from the previous pass through the loop.
+         Building c is newly-considered, and starts no earlier than b started.
+         The comments draw b as if its roof had zero slope ----.
+         with dashes where b lies above c.
+         The roof of c could rise / or fall \ through the roof of b,
+         or the vertical sides | of c could intersect the roof of b.  */
+      Building c = sc->front ();
+      if (b.end_ < c.end_) /* finish with b */
         {
-          result->push_front (b);
-          break;
+          if (b.end_ <= b.start_) /* we are already finished with b */
+            ;
+          else if (c.above (b, c.start_)) /* -|   . | */
+            {
+              Building m (b);
+              m.end_ = c.start_;
+              if (m.end_ > m.start_)
+                result->push_back (m);
+              if (b.above (c, b.end_))    /* -|\--.   */
+                {
+                  Building n (c);
+                  n.end_ = b.start_ = b.intersection_x (c);
+                  result->push_back (n);
+                  result->push_back (b);
+                }
+            }
+          else
+            {
+              if (c.above (b, b.end_))    /* ---/ . | */
+                b.end_ = b.intersection_x (c);
+              else                        /* -----.   */
+                c.start_ = b.end_;
+              result->push_back (b);
+            }
+          /* 'c' continues further, so move it into 'b' for the next pass. */
+          b = c;
+          swap (sb, sc);
         }
-
-      /* only include buildings wider than epsilon */
-      if (end > x + EPS)
+      else /* b.end_ > c.end_ so finish with c */
         {
-          b.leading_part (end);
-          result->push_front (b);
+          if (c.above (b, c.start_))    /* -| |---. */
+            {
+              Building m (b);
+              m.end_ = c.start_;
+              if (m.end_ > m.start_)
+                result->push_back (m);
+              if (b.above (c, c.end_))  /* -| \---. */
+                c.end_ = b.intersection_x (c);
+            }
+          else if (c.above (b, c.end_)) /* ---/|--. */
+            {
+              Building m (b);
+              c.start_ = m.end_ = b.intersection_x (c);
+              result->push_back (m);
+            }
+          else  /* c is completely hidden by b */
+            continue;
+          result->push_back (c);
+          b.start_ = c.end_;
         }
-
-      if (end >= s1->front ().end_)
-        s1->pop_front ();
-
-      x = end;
     }
-  result->reverse ();
+  if (b.end_ > b.start_)
+    result->push_back (b);
 }
 
 static void
@@ -240,89 +283,89 @@ empty_skyline (list<Building> *const ret)
 }
 
 /*
-  Given Building 'b' with starting wall location 'start', extend each side
-  with a sloped roofline of width 'horizon_padding'; put the skyline in 'ret'
+  Given Building 'b', build a skyline containing only that building.
 */
 static void
-single_skyline (Building b, Real start, Real horizon_padding, list<Building> *const ret)
+single_skyline (Building b, 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));
+  assert (b.end_ >= b.start_);
 
-  if (b.end_ > start + EPS)
-    ret->push_front (b);
-
-  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, start - horizon_padding));
+  if (b.start_ != -infinity_f)
+    ret->push_back (Building (-infinity_f, -infinity_f,
+                              -infinity_f, b.start_));
+  ret->push_back (b);
+  if (b.end_ != infinity_f)
+    ret->push_back (Building (b.end_, -infinity_f,
+                              -infinity_f, infinity_f));
 }
 
 /* remove a non-overlapping set of boxes from BOXES and build a skyline
    out of them */
 static list<Building>
-non_overlapping_skyline (list<Box> *const boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
+non_overlapping_skyline (list<Building> *const buildings)
 {
   list<Building> result;
   Real last_end = -infinity_f;
-  list<Box>::iterator i = boxes->begin ();
-  while (i != boxes->end ())
+  Building last_building (-infinity_f, -infinity_f, -infinity_f, infinity_f);
+  list<Building>::iterator i = buildings->begin ();
+  while (i != buildings->end ())
     {
-      Interval iv = (*i)[horizon_axis];
+      Real x1 = i->start_;
+      Real y1 = i->height (i->start_);
+      Real x2 = i->end_;
+      Real y2 = i->height (i->end_);
+
+      // Drop buildings that will obviously have no effect.
+      if (last_building.height (x1) >= y1
+          && last_building.end_ >= x2
+          && last_building.height (x2) >= y2)
+        {
+          list<Building>::iterator j = i++;
+          buildings->erase (j);
+          continue;
+        }
 
-      if (iv[LEFT] - horizon_padding < last_end)
+      if (x1 < 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));
+      // Insert empty Buildings into any gaps. (TODO: is this needed? -KOH)
+      if (x1 > last_end)
+        result.push_back (Building (last_end, -infinity_f, -infinity_f, x1));
 
-      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));
+      result.push_back (*i);
+      last_building = *i;
+      last_end = i->end_;
 
-      list<Box>::iterator j = i++;
-      boxes->erase (j);
-      last_end = result.front ().end_;
+      list<Building>::iterator j = i++;
+      buildings->erase (j);
     }
+
   if (last_end < infinity_f)
-    result.push_front (Building (last_end, -infinity_f, -infinity_f, infinity_f));
-  result.reverse ();
+    result.push_back (Building (last_end, -infinity_f, -infinity_f, infinity_f));
   return result;
 }
 
-class LessThanBox
+class LessThanBuilding
 {
-  Axis a_;
-
 public:
-  LessThanBox (Axis a)
-  {
-    a_ = a;
-  }
-
-  bool operator () (Box const &b1, Box const &b2)
+  bool operator () (Building const &b1, Building const &b2)
   {
-    return b1[a_][LEFT] < b2[a_][LEFT];
+    return b1.start_ < b2.start_
+           || (b1.start_ == b2.start_ && b1.height (b1.start_) > b2.height (b1.start_));
   }
 };
 
+/**
+   BUILDINGS is a list of buildings, but they could be overlapping
+   and in any order.  The returned list of buildings is ordered and non-overlapping.
+*/
 list<Building>
-Skyline::internal_build_skyline (list<Box> *boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
+Skyline::internal_build_skyline (list<Building> *buildings) const
 {
-  vsize size = boxes->size ();
+  vsize size = buildings->size ();
 
   if (size == 0)
     {
@@ -333,16 +376,14 @@ Skyline::internal_build_skyline (list<Box> *boxes, Real horizon_padding, Axis ho
   else if (size == 1)
     {
       list<Building> result;
-      single_skyline (Building (boxes->front (), horizon_padding, horizon_axis, sky),
-                      boxes->front ()[horizon_axis][LEFT] - horizon_padding,
-                      horizon_padding, &result);
+      single_skyline (buildings->front (), &result);
       return result;
     }
 
   deque<list<Building> > partials;
-  boxes->sort (LessThanBox (horizon_axis));
-  while (!boxes->empty ())
-    partials.push_back (non_overlapping_skyline (boxes, horizon_padding, horizon_axis, sky));
+  buildings->sort (LessThanBuilding ());
+  while (!buildings->empty ())
+    partials.push_back (non_overlapping_skyline (buildings));
 
   /* we'd like to say while (partials->size () > 1) but that's O (n).
      Instead, we exit in the middle of the loop */
@@ -369,18 +410,6 @@ Skyline::Skyline ()
   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++)
-    {
-      buildings_.push_back (Building ((*i)));
-    }
-}
-
 Skyline::Skyline (Direction sky)
 {
   sky_ = sky;
@@ -388,68 +417,90 @@ Skyline::Skyline (Direction sky)
 }
 
 /*
-  build padded skyline from an existing skyline with padding
-  added to it.
-*/
+  Build skyline from a set of boxes.
 
-Skyline::Skyline (Skyline const &src, Real horizon_padding, Axis /*a*/)
+  Boxes should be non-empty on both axes.  Otherwise, they will be ignored
+ */
+Skyline::Skyline (vector<Box> const &boxes, Axis horizon_axis, Direction sky)
 {
-  /*
-     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;
+  list<Building> buildings;
+  sky_ = sky;
 
-  // 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_;
+  for (vsize i = 0; i < boxes.size (); i++)
+    if (!boxes[i].is_empty (X_AXIS)
+        && !boxes[i].is_empty (Y_AXIS))
+      buildings.push_front (Building (boxes[i], horizon_axis, sky));
+
+  buildings_ = internal_build_skyline (&buildings);
+  normalize ();
 }
 
 /*
-  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
-  width horizon_padding). That is, the total amount of horizontal expansion is
-  horizon_padding*4, half of which is sloped and half of which is flat.
+  build skyline from a set of line segments.
 
-  Boxes should have fatness in the horizon_axis (after they are expanded by
-  horizon_padding), otherwise they are ignored.
+  Segments can be articulated from left to right or right to left.
+  In the case of the latter, they will be stored internally as left to right.
  */
-Skyline::Skyline (vector<Box> const &boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
+Skyline::Skyline (vector<Drul_array<Offset> > const &segments, Axis horizon_axis, Direction sky)
 {
-  list<Box> filtered_boxes;
+  list<Building> buildings;
   sky_ = sky;
 
-  Axis vert_axis = other_axis (horizon_axis);
-  for (vsize i = 0; i < boxes.size (); i++)
+  for (vsize i = 0; i < segments.size (); i++)
+    {
+      Drul_array<Offset> const &seg = segments[i];
+      Offset left = seg[LEFT];
+      Offset right = seg[RIGHT];
+      if (left[horizon_axis] > right[horizon_axis])
+        swap (left, right);
+
+      Real x1 = left[horizon_axis];
+      Real x2 = right[horizon_axis];
+      Real y1 = left[other_axis (horizon_axis)] * sky;
+      Real y2 = right[other_axis (horizon_axis)] * sky;
+
+      if (x1 <= x2)
+        buildings.push_back (Building (x1, y1, y2, x2));
+    }
+
+  buildings_ = internal_build_skyline (&buildings);
+  normalize ();
+}
+
+Skyline::Skyline (vector<Skyline_pair> const &skypairs, Direction sky)
+{
+  sky_ = sky;
+
+  deque<Skyline> partials;
+  for (vsize i = 0; i < skypairs.size (); i++)
+    partials.push_back (Skyline ((skypairs[i])[sky]));
+
+  while (partials.size () > 1)
     {
-      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]);
+      Skyline one = partials.front ();
+      partials.pop_front ();
+      Skyline two = partials.front ();
+      partials.pop_front ();
+
+      one.merge (two);
+      partials.push_back (one);
     }
 
-  buildings_ = internal_build_skyline (&filtered_boxes, horizon_padding, horizon_axis, sky);
+  if (partials.size ())
+    buildings_.swap (partials.front ().buildings_);
+  else
+    buildings_.clear ();
 }
 
-Skyline::Skyline (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
+Skyline::Skyline (Box const &b, Axis horizon_axis, Direction sky)
 {
   sky_ = sky;
-  Building front (b, horizon_padding, horizon_axis, sky);
-  single_skyline (front, b[horizon_axis][LEFT] - horizon_padding,
-                  horizon_padding, &buildings_);
+  if (!b.is_empty (X_AXIS) && !b.is_empty (Y_AXIS))
+    {
+      Building front (b, horizon_axis, sky);
+      single_skyline (front, &buildings_);
+      normalize ();
+    }
 }
 
 void
@@ -457,14 +508,24 @@ Skyline::merge (Skyline const &other)
 {
   assert (sky_ == other.sky_);
 
+  if (other.is_empty ())
+    return;
+
+  if (is_empty ())
+    {
+      buildings_ = other.buildings_;
+      return;
+    }
+
   list<Building> other_bld (other.buildings_);
   list<Building> my_bld;
   my_bld.splice (my_bld.begin (), buildings_);
   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
+  normalize ();
 }
 
 void
-Skyline::insert (Box const &b, Real horizon_padding, Axis a)
+Skyline::insert (Box const &b, Axis a)
 {
   list<Building> other_bld;
   list<Building> my_bld;
@@ -477,15 +538,13 @@ Skyline::insert (Box const &b, Real horizon_padding, Axis a)
     }
 
   /* 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 ())
+  if (b.is_empty (X_AXIS) || b.is_empty (Y_AXIS))
     return;
 
   my_bld.splice (my_bld.begin (), buildings_);
-  single_skyline (Building (b, horizon_padding, a, sky_), b[a][LEFT] - horizon_padding,
-                  horizon_padding, &other_bld);
+  single_skyline (Building (b, a, sky_), &other_bld);
   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
+  normalize ();
 }
 
 void
@@ -502,6 +561,7 @@ Skyline::shift (Real s)
   list<Building>::iterator end = buildings_.end ();
   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
     {
+      i->start_ += s;
       i->end_ += s;
       i->y_intercept_ -= s * i->slope_;
     }
@@ -510,36 +570,113 @@ Skyline::shift (Real s)
 Real
 Skyline::distance (Skyline const &other, Real horizon_padding) const
 {
-  assert (sky_ == -other.sky_);
+  Real dummy;
+  return internal_distance (other, horizon_padding, &dummy);
+}
 
-  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)
+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
+{
+  if (horizon_padding == 0.0)
+    return internal_distance (other, touch_point);
+
+  // Note that it is not necessary to build a padded version of other,
+  // because the same effect can be achieved just by doubling horizon_padding.
+  Skyline padded_this = padded (horizon_padding);
+  return padded_this.internal_distance (other, touch_point);
+}
+
+Skyline
+Skyline::padded (Real horizon_padding) const
+{
+  if (horizon_padding < 0.0)
+    warning ("Cannot have negative horizon padding.  Junking.");
+
+  if (horizon_padding <= 0.0)
+    return *this;
+
+  list<Building> pad_buildings;
+  for (list<Building>::const_iterator i = buildings_.begin (); i != buildings_.end (); ++i)
     {
-      padded_this = new Skyline (*padded_this, horizon_padding, X_AXIS);
-      padded_other = new Skyline (*padded_other, horizon_padding, X_AXIS);
-      created_tmp_skylines = true;
+      if (i->start_ > -infinity_f)
+        {
+          Real height = i->height (i->start_);
+          if (height > -infinity_f)
+            {
+              // Add the sloped building that pads the left side of the current building.
+              Real start = i->start_ - 2 * horizon_padding;
+              Real end = i->start_ - horizon_padding;
+              pad_buildings.push_back (Building (start, height - horizon_padding, height, end));
+
+              // Add the flat building that pads the left side of the current building.
+              start = i->start_ - horizon_padding;
+              end = i->start_;
+              pad_buildings.push_back (Building (start, height, height, end));
+            }
+        }
+
+      if (i->end_ < infinity_f)
+        {
+          Real height = i->height (i->end_);
+          if (height > -infinity_f)
+            {
+              // Add the flat building that pads the right side of the current building.
+              Real start = i->end_;
+              Real end = start + horizon_padding;
+              pad_buildings.push_back (Building (start, height, height, end));
+
+              // Add the sloped building that pads the right side of the current building.
+              start = end;
+              end += horizon_padding;
+              pad_buildings.push_back (Building (start, height, height - horizon_padding, end));
+            }
+        }
     }
 
-  list<Building>::const_iterator i = padded_this->buildings_.begin ();
-  list<Building>::const_iterator j = padded_other->buildings_.begin ();
+  // The buildings may be overlapping, so resolve that.
+  list<Building> pad_skyline = internal_build_skyline (&pad_buildings);
+
+  // Merge the padding with the original, to make a new skyline.
+  Skyline padded (sky_);
+  list<Building> my_buildings = buildings_;
+  padded.buildings_.clear ();
+  internal_merge_skyline (&pad_skyline, &my_buildings, &padded.buildings_);
+  padded.normalize ();
+
+  return padded;
+}
+
+Real
+Skyline::internal_distance (Skyline const &other, Real *touch_point) const
+{
+  assert (sky_ == -other.sky_);
+
+  list<Building>::const_iterator i = buildings_.begin ();
+  list<Building>::const_iterator j = other.buildings_.begin ();
 
   Real dist = -infinity_f;
   Real start = -infinity_f;
-  while (i != padded_this->buildings_.end () && j != padded_other->buildings_.end ())
+  Real touch = -infinity_f;
+  while (i != buildings_.end () && j != other.buildings_.end ())
     {
       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
@@ -547,12 +684,7 @@ Skyline::distance (Skyline const &other, Real horizon_padding) const
       start = end;
     }
 
-  if (created_tmp_skylines)
-    {
-      delete padded_this;
-      delete padded_other;
-    }
-
+  *touch_point = touch;
   return dist;
 }
 
@@ -574,10 +706,53 @@ Skyline::height (Real airplane) const
 
 Real
 Skyline::max_height () const
+{
+  Real ret = -infinity_f;
+
+  list<Building>::const_iterator i;
+  for (i = buildings_.begin (); i != buildings_.end (); ++i)
+    {
+      ret = max (ret, i->height (i->start_));
+      ret = max (ret, i->height (i->end_));
+    }
+
+  return sky_ * ret;
+}
+
+Direction
+Skyline::direction () const
+{
+  return sky_;
+}
+
+Real
+Skyline::left () const
+{
+  for (list<Building>::const_iterator i (buildings_.begin ());
+       i != buildings_.end (); i++)
+    if (i->y_intercept_ > -infinity_f)
+      return i->start_;
+
+  return infinity_f;
+}
+
+Real
+Skyline::right () const
+{
+  for (list<Building>::const_reverse_iterator i (buildings_.rbegin ());
+       i != buildings_.rend (); ++i)
+    if (i->y_intercept_ > -infinity_f)
+      return i->end_;
+
+  return -infinity_f;
+}
+
+Real
+Skyline::max_height_position () const
 {
   Skyline s (-sky_);
   s.set_minimum_height (0);
-  return sky_ * distance (s);
+  return touching_point (s);
 }
 
 void
@@ -612,6 +787,8 @@ Skyline::to_points (Axis horizon_axis) const
 bool
 Skyline::is_empty () const
 {
+  if (!buildings_.size ())
+    return true;
   Building b = buildings_.front ();
   return b.end_ == infinity_f && b.y_intercept_ == -infinity_f;
 }
@@ -625,24 +802,71 @@ Skyline::clear ()
 
 /****************************************************************/
 
-IMPLEMENT_SIMPLE_SMOBS (Skyline);
-IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
-IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
+const char * const Skyline::type_p_name_ = "ly:skyline?";
+
+MAKE_SCHEME_CALLBACK_WITH_OPTARGS (Skyline, get_touching_point, 3, 1, "")
+SCM
+Skyline::get_touching_point (SCM skyline_scm, SCM other_skyline_scm, SCM horizon_padding_scm)
+{
+  LY_ASSERT_SMOB (Skyline, other_skyline_scm, 1);
+
+  Real horizon_padding = 0;
+  if (!SCM_UNBNDP (horizon_padding_scm))
+    {
+      LY_ASSERT_TYPE (scm_is_number, horizon_padding_scm, 3);
+      horizon_padding = scm_to_double (horizon_padding_scm);
+    }
+
+  Skyline *skyline = unsmob<Skyline> (skyline_scm);
+  Skyline *other_skyline = unsmob<Skyline> (other_skyline_scm);
+  return scm_from_double (skyline->touching_point (*other_skyline, horizon_padding));
+}
 
+MAKE_SCHEME_CALLBACK_WITH_OPTARGS (Skyline, get_distance, 3, 1, "")
 SCM
-Skyline::mark_smob (SCM)
+Skyline::get_distance (SCM skyline_scm, SCM other_skyline_scm, SCM horizon_padding_scm)
 {
-  ASSERT_LIVE_IS_ALLOWED ();
-  return SCM_EOL;
+  LY_ASSERT_SMOB (Skyline, other_skyline_scm, 1);
+
+  Real horizon_padding = 0;
+  if (!SCM_UNBNDP (horizon_padding_scm))
+    {
+      LY_ASSERT_TYPE (scm_is_number, horizon_padding_scm, 3);
+      horizon_padding = scm_to_double (horizon_padding_scm);
+    }
+
+  Skyline *skyline = unsmob<Skyline> (skyline_scm);
+  Skyline *other_skyline = unsmob<Skyline> (other_skyline_scm);
+  return scm_from_double (skyline->distance (*other_skyline, horizon_padding));
 }
 
-int
-Skyline::print_smob (SCM s, SCM port, scm_print_state *)
+MAKE_SCHEME_CALLBACK (Skyline, get_max_height, 1)
+SCM
+Skyline::get_max_height (SCM skyline_scm)
 {
-  Skyline *r = (Skyline *) SCM_CELL_WORD_1 (s);
-  (void) r;
+  return scm_from_double (unsmob<Skyline> (skyline_scm)->max_height ());
+}
 
-  scm_puts ("#<Skyline>", port);
+MAKE_SCHEME_CALLBACK (Skyline, get_max_height_position, 1)
+SCM
+Skyline::get_max_height_position (SCM skyline_scm)
+{
+  return scm_from_double (unsmob<Skyline> (skyline_scm)->max_height_position ());
+}
 
-  return 1;
+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 (unsmob<Skyline> (skyline_scm)->height (x));
+}
+
+LY_DEFINE (ly_skyline_empty_p, "ly:skyline-empty?",
+           1, 0, 0, (SCM sky),
+           "Return whether @var{sky} is empty.")
+{
+  Skyline *s = unsmob<Skyline> (sky);
+  LY_ASSERT_SMOB (Skyline, sky, 1);
+  return scm_from_bool (s->is_empty ());
 }