]> git.donarmstrong.com Git - lilypond.git/blobdiff - lily/skyline.cc
Merge remote-tracking branch 'origin/translation'
[lilypond.git] / lily / skyline.cc
index b3b6349461a25369fffab2aaf0f954a05f359c71..33d2823d47e0d752934558f75ea74e90e85ead41 100644 (file)
 /*
-  skyline.cc -- implement Skyline_entry and funcs.
+  This file is part of LilyPond, the GNU music typesetter.
 
-  source file of the GNU LilyPond music typesetter
+  Copyright (C) 2006--2014 Joe Neeman <joeneeman@gmail.com>
 
-  (c) 2002--2005 Han-Wen Nienhuys <hanwen@cs.uu.nl>
+  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.
+
+  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 "skyline-pair.hh"
+#include <deque>
+#include <cstdio>
+
+#include "ly-smobs.icc"
+
+/* A skyline is a sequence of non-overlapping buildings: something like
+   this:
+                   _______
+                  |       \                                 ________
+                  |        \                       ________/        \
+        /\        |          \                    /                  \
+       /  --------             \                 /                    \
+      /                          \              /                      \
+     /                             ------------/                        ----
+   --
+   Each building has a starting position, and ending position, a starting
+   height and an ending height.
+
+   The following invariants are observed:
+    - the start of the first building is at -infinity
+    - the end of the last building is at infinity
+    - if a building has infinite length (ie. the first and last buildings),
+      then its starting height and ending height are equal
+    - the end of one building is the same as the beginning of the next
+      building
+
+   We also allow skylines to point down (the structure is exactly the same,
+   but we think of the part above the line as being filled with mass and the
+   part below as being empty). ::distance finds the minimum distance between
+   an UP skyline and a DOWN skyline.
+
+   Note that we store DOWN skylines upside-down. That is, in order to compare
+   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.
+
+   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.
+*/
 
-/*
-  A skyline is a shape of the form:
+static void
+print_buildings (list<Building> const &b)
+{
+  for (list<Building>::const_iterator i = b.begin (); i != b.end (); i++)
+    i->print ();
+}
 
+void
+Skyline::print () const
+{
+  print_buildings (buildings_);
+}
 
-  ----
-  |  |
-  ---------|  |
-  |          |
-  |          |
-  |          |______
-  --------|                 |___
+void
+Skyline::print_points () const
+{
+  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)
+{
+  if (isinf (start) || isinf (end))
+    assert (start_height == end_height);
 
-  This file deals with building such skyline structure, and computing
-  the minimum distance between two opposing skylines.
+  start_ = start;
+  end_ = end;
+  precompute (start, start_height, end_height, end);
+}
 
+Building::Building (Box const &b, Axis horizon_axis, Direction sky)
+{
+  Real start = b[horizon_axis][LEFT];
+  Real end = b[horizon_axis][RIGHT];
+  Real height = sky * b[other_axis (horizon_axis)][sky];
 
-  Invariants for a skyline:
+  start_ = start;
+  end_ = end;
+  precompute (start, height, height, end);
+}
 
-  skyline[...].width_ forms a partition of the real interval, where
-  the segments are adjacent, and ascending. Hence we have
+void
+Building::precompute (Real start, Real start_height, Real end_height, Real end)
+{
+  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);
 
-  skyline.top ().width_[RIGHT] = inf
-  skyline[0].width_[LEFT] = -inf
-*/
+  assert (!isinf (slope_) && !isnan (slope_));
 
-const Real EPS = 1e-12;
+  if (isinf (start))
+    {
+      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;
+}
 
-/*
-  TODO: avoid unnecessary fragmentation.
+inline Real
+Building::height (Real x) const
+{
+  return isinf (x) ? y_intercept_ : slope_ * x + y_intercept_;
+}
 
-  This is O (n^2), searching and insertion.  Could be O (n log n) with
-  binsearch.
-*/
 void
-insert_extent_into_skyline (Array<Skyline_entry> *line, Box b, Axis line_axis,
-                           Direction d)
+Building::print () const
 {
-  Interval extent = b[line_axis];
-  if (extent.is_empty ())
-    return;
+  printf ("%f x + %f from %f to %f\n", slope_, y_intercept_, start_, end_);
+}
 
-  Real stick_out = b[other_axis (line_axis)][d];
+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)
+{
+  assert (chop <= end_);
+  end_ = chop;
+}
 
-  /*
-    Intersect each segment of LINE with EXTENT, and if non-empty, insert relevant segments.
-  */
-  for (int i = line->size (); i--;)
+// 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
+{
+  // 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;
+}
+
+static Real
+first_intersection (Building const &b, list<Building> *s, Real start_x)
+/* Return the first x >= start_x where skyline s above Building b.
+ * Removes buildings from s that are concealed by b. */
+{
+  while (!s->empty () && start_x < b.end_)
     {
-      Interval w = line->elem (i).width_;
-      w.intersect (extent);
+      Building c = s->front ();
+
+      // conceals and intersection_x involve multiplication and
+      // division. Avoid that, if we can.
+      if (c.y_intercept_ == -infinity_f)
+        {
+          if (c.end_ > b.end_)
+            return b.end_;
+          start_x = c.end_;
+          s->pop_front ();
+          continue;
+        }
+
+      if (c.conceals (b, start_x))
+        return start_x;
+
+      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 ();
+    }
+  return b.end_;
+}
+
+bool
+Building::conceals (Building const &other, Real x) const
+{
+  if (slope_ == other.slope_)
+    return y_intercept_ > other.y_intercept_;
 
-      if (extent[LEFT] >= w[RIGHT])
-       break;
+  /* 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_);
+}
 
-      Real my_height = line->elem (i).height_;
+// Remove redundant empty buildings from the skyline.
+// If there are two adjacent empty buildings, they can be
+// turned into one.
+void
+Skyline::normalize ()
+{
+  bool last_empty = false;
+  list<Building>::iterator i;
 
-      if (!w.is_empty ()
-         && w.length () > EPS
-         && d * (my_height - stick_out) < 0)
-       {
-         Interval e1 (line->elem (i).width_[LEFT], extent[LEFT]);
-         Interval e3 (extent[RIGHT], line->elem (i).width_[RIGHT]);
+  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);
+    }
 
-         if (!e3.is_empty () && e3.length () > EPS)
-           line->insert (Skyline_entry (e3, my_height), i + 1);
+  assert (buildings_.front ().start_ == -infinity_f);
+  assert (buildings_.back ().end_ == infinity_f);
+}
 
-         line->elem_ref (i).height_ = stick_out;
-         line->elem_ref (i).width_ = w;
-         if (!e1.is_empty () && e1.length () > EPS)
-           line->insert (Skyline_entry (e1, my_height), i);
-       }
+void
+Skyline::deholify ()
+{
+  // Since a skyline should always be normalized, we can
+  // assume that there are never two adjacent empty buildings.
+  // That is, if center is empty then left and right are not.
+  list<Building>::iterator left = buildings_.begin ();
+  list<Building>::iterator center = buildings_.begin ();
+  list<Building>::iterator right;
+
+  for (right = buildings_.begin (); right != buildings_.end (); right++)
+    {
+      if (center != buildings_.begin () && center->y_intercept_ == -infinity_f)
+        {
+          Real p1 = left->height (left->end_);
+          Real p2 = right->height (right->start_);
+          *center = Building (center->start_, p1, p2, center->end_);
+
+          left = center;
+          center = right;
+        }
     }
 }
 
 void
-merge_skyline (Array<Skyline_entry> *a1,
-              Array<Skyline_entry> const &a2,
-              Direction dir)
+Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
+                                 list<Building> *const result) const
 {
-  for (int i = 0; i < a2.size (); i++)
+  if (s1->empty () || s2->empty ())
     {
-      Box b;
-      b[X_AXIS] = a2[i].width_;
-      b[Y_AXIS][dir] = a2[i].height_;
-      b[Y_AXIS][-dir] = dir * infinity_f;
+      programming_error ("tried to merge an empty skyline");
+      return;
+    }
 
-      insert_extent_into_skyline (a1, b, X_AXIS, dir);
+  Real x = -infinity_f;
+  Real last_end = -infinity_f;
+  while (!s1->empty ())
+    {
+      if (s2->front ().conceals (s1->front (), x))
+        swap (s1, s2);
+
+      Building b = s1->front ();
+      Building c = s2->front ();
+
+      // Optimization: if the other skyline is empty at this point,
+      // we can avoid testing some intersections. Just grab as many
+      // buildings from s1 as we can, and shove them onto the output.
+      if (c.y_intercept_ == -infinity_f
+          && c.end_ >= b.end_)
+        {
+          list<Building>::iterator i = s1->begin ();
+          i++;
+          while (i != s1->end () && i->end_ <= c.end_)
+            i++;
+
+          s1->front ().start_ = x;
+          result->splice (result->end (), *s1, s1->begin (), i);
+          x = result->back ().end_;
+          last_end = x;
+          continue;
+        }
+      // first_intersection() removes buildings from s2 if b hides them
+      Real end = first_intersection (b, s2, x);
+      if (s2->empty ())
+        {
+          b.start_ = last_end;
+          result->push_back (b);
+          break;
+        }
+
+      // Should be (end > x), during ver2.19.  end == x happens fairly often,
+      // and we do not need to keep vertical segments within a skyline.
+      if (end >= x)
+        {
+          b.leading_part (end);
+          b.start_ = last_end;
+          last_end = b.end_;
+          result->push_back (b);
+        }
+
+      if (end >= s1->front ().end_)
+        s1->pop_front ();
+      // Should add during ver2.19 (to avoid an endless loop
+      // when  merging identical skylines with a vertical segment)
+      // if (end >= s2->front().end_) s2->pop_front();
+
+      x = end;
     }
 }
 
-Array<Skyline_entry>
-empty_skyline (Direction d)
+static void
+empty_skyline (list<Building> *const ret)
 {
-  Array<Skyline_entry> skyline;
+  ret->push_front (Building (-infinity_f, -infinity_f, -infinity_f, infinity_f));
+}
 
-  Interval i;
-  i.set_empty ();
-  i.swap ();
-  Skyline_entry e;
-  e.width_ = i;
-  e.height_ = -d * infinity_f;
-  skyline.push (e);
-  return skyline;
+/*
+  Given Building 'b', build a skyline containing only that building.
+*/
+static void
+single_skyline (Building b, list<Building> *const ret)
+{
+  assert (b.end_ >= b.start_);
+
+  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<Building> *const buildings)
+{
+  list<Building> result;
+  Real last_end = -infinity_f;
+  Building last_building (-infinity_f, -infinity_f, -infinity_f, infinity_f);
+  list<Building>::iterator i = buildings->begin ();
+  while (i != buildings->end ())
+    {
+      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 (x1 < last_end)
+        {
+          i++;
+          continue;
+        }
+
+      // 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));
+
+      result.push_back (*i);
+      last_building = *i;
+      last_end = i->end_;
+
+      list<Building>::iterator j = i++;
+      buildings->erase (j);
+    }
+
+  if (last_end < infinity_f)
+    result.push_back (Building (last_end, -infinity_f, -infinity_f, infinity_f));
+  return result;
+}
+
+class LessThanBuilding
+{
+public:
+  bool operator () (Building const &b1, Building const &b2)
+  {
+    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<Building> *buildings) const
+{
+  vsize size = buildings->size ();
+
+  if (size == 0)
+    {
+      list<Building> result;
+      empty_skyline (&result);
+      return result;
+    }
+  else if (size == 1)
+    {
+      list<Building> result;
+      single_skyline (buildings->front (), &result);
+      return result;
+    }
+
+  deque<list<Building> > partials;
+  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 */
+  while (!partials.empty ())
+    {
+      list<Building> merged;
+      list<Building> one = partials.front ();
+      partials.pop_front ();
+      if (partials.empty ())
+        return one;
+
+      list<Building> two = partials.front ();
+      partials.pop_front ();
+      internal_merge_skyline (&one, &two, &merged);
+      partials.push_back (merged);
+    }
+  assert (0);
+  return list<Building> ();
 }
 
-Array<Skyline_entry>
-extents_to_skyline (Array<Box> const &extents, Axis a, Direction d)
+Skyline::Skyline ()
 {
+  sky_ = UP;
+  empty_skyline (&buildings_);
+}
 
-  Array<Skyline_entry> skyline = empty_skyline (d);
+Skyline::Skyline (Direction sky)
+{
+  sky_ = sky;
+  empty_skyline (&buildings_);
+}
+
+/*
+  Build skyline from a set of boxes.
 
-  /*
-    This makes a cubic algorithm (array  insertion is O (n),
-    searching the array dumbly is O (n), and for n items, we get O (n^3).)
+  Boxes should be non-empty on both axes.  Otherwise, they will be ignored
+ */
+Skyline::Skyline (vector<Box> const &boxes, Axis horizon_axis, Direction sky)
+{
+  list<Building> buildings;
+  sky_ = sky;
 
-    We could do a lot better (n log (n), using a balanced tree) but
-    that seems overkill for now.
-  */
-  for (int j = extents.size (); j--;)
-    insert_extent_into_skyline (&skyline, extents[j], a, d);
+  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));
 
-  return skyline;
+  buildings_ = internal_build_skyline (&buildings);
+  normalize ();
 }
 
 /*
-  minimum distance that can be achieved between baselines. "Clouds" is
-  a skyline pointing down.
+  build skyline from a set of line segments.
+
+  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<Drul_array<Offset> > const &segments, Axis horizon_axis, Direction sky)
+{
+  list<Building> buildings;
+  sky_ = sky;
+
+  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)
+    {
+      Skyline one = partials.front ();
+      partials.pop_front ();
+      Skyline two = partials.front ();
+      partials.pop_front ();
+
+      one.merge (two);
+      partials.push_back (one);
+    }
+
+  if (partials.size ())
+    buildings_.swap (partials.front ().buildings_);
+  else
+    buildings_.clear ();
+}
+
+Skyline::Skyline (Box const &b, Axis horizon_axis, Direction sky)
+{
+  sky_ = sky;
+  if (!b.is_empty (X_AXIS) && !b.is_empty (Y_AXIS))
+    {
+      Building front (b, horizon_axis, sky);
+      single_skyline (front, &buildings_);
+      normalize ();
+    }
+}
+
+void
+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, 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.) */
+  if (b.is_empty (X_AXIS) || b.is_empty (Y_AXIS))
+    return;
+
+  my_bld.splice (my_bld.begin (), buildings_);
+  single_skyline (Building (b, a, sky_), &other_bld);
+  internal_merge_skyline (&other_bld, &my_bld, &buildings_);
+  normalize ();
+}
+
+void
+Skyline::raise (Real r)
+{
+  list<Building>::iterator end = buildings_.end ();
+  for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
+    i->y_intercept_ += sky_ * r;
+}
+
+void
+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_;
+    }
+}
+
+Real
+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
+{
+  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)
+    {
+      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));
+            }
+        }
+    }
+
+  // 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;
+}
 
-  This is an O (n) algorithm.
-*/
 Real
-skyline_meshing_distance (Array<Skyline_entry> const &buildings,
-                         Array<Skyline_entry> const &clouds)
+Skyline::internal_distance (Skyline const &other, Real *touch_point) const
 {
-  int i = buildings.size () -1;
-  int j = clouds.size () -1;
+  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;
+  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
+        j++;
+      start = end;
+    }
+
+  *touch_point = touch;
+  return dist;
+}
 
-  Real distance = -infinity_f;
+Real
+Skyline::height (Real airplane) const
+{
+  assert (!isinf (airplane));
 
-  while (i > 0 || j > 0)
+  list<Building>::const_iterator i;
+  for (i = buildings_.begin (); i != buildings_.end (); i++)
     {
-      Interval w = buildings[i].width_;
-      w.intersect (clouds[j].width_);
+      if (i->end_ >= airplane)
+        return sky_ * i->height (airplane);
+    }
 
-      if (!w.is_empty ())
-       distance = distance >? (buildings[i].height_ - clouds[j].height_);
+  assert (0);
+  return 0;
+}
 
-      if (i > 0 && buildings[i].width_[LEFT] >= clouds[j].width_[LEFT])
-       {
-         i--;
-       }
-      else if (j > 0 && buildings[i].width_[LEFT] <= clouds[j].width_[LEFT])
-       {
-         j--;
-       }
+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 distance;
+  return sky_ * ret;
 }
 
-Skyline_entry::Skyline_entry ()
+Direction
+Skyline::direction () const
 {
-  height_ = 0.0;
+  return sky_;
 }
 
-Skyline_entry::Skyline_entry (Interval i, Real r)
+Real
+Skyline::left () const
 {
-  width_ = i;
-  height_ = r;
+  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 touching_point (s);
+}
+
+void
+Skyline::set_minimum_height (Real h)
+{
+  Skyline s (sky_);
+  s.buildings_.front ().y_intercept_ = h * sky_;
+  merge (s);
+}
+
+vector<Offset>
+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++)
+    {
+      out.push_back (Offset (start, sky_ * i->height (start)));
+      out.push_back (Offset (i->end_, sky_ * i->height (i->end_)));
+      start = i->end_;
+    }
+
+  if (horizon_axis == Y_AXIS)
+    for (vsize i = 0; i < out.size (); i++)
+      out[i] = out[i].swapped ();
+
+  return out;
+}
+
+bool
+Skyline::is_empty () const
+{
+  if (!buildings_.size ())
+    return true;
+  Building b = buildings_.front ();
+  return b.end_ == infinity_f && b.y_intercept_ == -infinity_f;
 }
 
 void
-heighten_skyline (Array<Skyline_entry> *buildings, Real ground)
+Skyline::clear ()
+{
+  buildings_.clear ();
+  empty_skyline (&buildings_);
+}
+
+/****************************************************************/
+
+IMPLEMENT_SIMPLE_SMOBS (Skyline);
+IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
+IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
+
+SCM
+Skyline::mark_smob (SCM s)
+{
+  ASSERT_LIVE_IS_ALLOWED (s);
+  return SCM_EOL;
+}
+
+int
+Skyline::print_smob (SCM s, SCM port, scm_print_state *)
+{
+  Skyline *r = (Skyline *) SCM_CELL_WORD_1 (s);
+  (void) r;
+
+  scm_puts ("#<Skyline>", port);
+
+  return 1;
+}
+
+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 (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));
+}
+
+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)
+{
+  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->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));
+}
+
+LY_DEFINE (ly_skyline_empty_p, "ly:skyline-empty?",
+           1, 0, 0, (SCM sky),
+           "Return whether @var{sky} is empty.")
 {
-  for (int i = 0; i < buildings->size (); i++)
-    buildings->elem_ref (i).height_ += ground;
+  Skyline *s = Skyline::unsmob (sky);
+  LY_ASSERT_SMOB (Skyline, sky, 1);
+  return scm_from_bool (s->is_empty ());
 }