]> git.donarmstrong.com Git - lilypond.git/blobdiff - lily/skyline.cc
fix Skyline::distance calculation
[lilypond.git] / lily / skyline.cc
index 6c2007ed5babf021f6b49a8d8413dd4ee72c1a09..e52a6277d50adfeeb44b6f9369944f8c4a900909 100644 (file)
-/*
-  skyline.cc -- implement Skyline_entry and funcs.
+/* skyline.cc -- implement the Skyline class
 
-  source file of the GNU LilyPond music typesetter
-
-  (c) 2002--2006 Han-Wen Nienhuys <hanwen@xs4all.nl>
+   source file of the GNU LilyPond music typesetter
+   (c) 2006 Joe Neeman <joeneeman@gmail.com>
 */
 
 #include "skyline.hh"
 
-/*
-  A skyline is a shape of the form:
-
-
-  *                  ----
-  *                  |  |
-  *         ---------|  |
-  *         |           |
-  *         |           |
-  *         |          |______
-  * --------|                |___
-  *
-
-  This file deals with building such skyline structure, and computing
-  the minimum distance between two opposing skylines.
-
-  Invariants for a skyline:
-
-  skyline[...].width_ forms a partition of the real interval, where
-  the segments are adjacent, and ascending. Hence we have
-
-  skyline.back ().width_[RIGHT] = inf
-  skyline[0].width_[LEFT] = -inf
+#include "line-interface.hh"
+
+/* 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.
 */
 
-const Real EPS = 1e-12;
-
-/*
-  TODO: avoid unnecessary fragmentation.
+#define EPS 1e-10
 
-  This is O (n^2), searching and insertion.  Could be O (n log n) with
-  binsearch.
-*/
-void
-insert_extent_into_skyline (std::vector<Skyline_entry> *line, Box b, Axis line_axis,
-                           Direction d)
+static inline bool
+equal (Real x, Real y)
 {
-  Interval extent = b[line_axis];
-  if (extent.is_empty ())
-    return;
+  return abs (x - y) < EPS || (isinf (x) && isinf (y) && ((x > 0) == (y > 0)));
+}
 
-  Real stick_out = b[other_axis (line_axis)][d];
+bool
+Skyline::is_legal_skyline () const
+{
+  list<Building>::const_iterator i;
+  Real last_x = -infinity_f;
+  Real last_h = -infinity_f;
+  for (i = buildings_.begin (); i != buildings_.end (); i++)
+    {
+      if (i->iv_[LEFT] != last_x)
+       return false;
+      if (i != buildings_.begin () && !equal (i->start_height_, last_h))
+       return false;
+      last_x = i->iv_[RIGHT];
+      last_h = i->end_height_;
+    }
+  return last_x == infinity_f;
+}
 
-  /*
-    Intersect each segment of LINE with EXTENT, and if non-empty, insert relevant segments.
-  */
-  for (vsize i = line->size (); i--;)
+Building::Building (Real start, Real start_height, Real end_height, Real end, Real max_slope)
+  : iv_ (start, end)
+{
+  start_height_ = start_height;
+  end_height_ = end_height;
+
+  if (isinf (start))
+    assert (isinf (start_height) || start_height == end_height);
+  if (isinf (end))
+    assert (isinf (end_height) || start_height == end_height);
+
+  m_ = (end_height - start_height) / (end - start);
+  if (start_height == end_height)
+    m_ = 0;
+  if (isinf (m_) || isnan (m_))
+    m_ = max_slope * (start_height < end_height ? 1 : -1);
+  assert (abs (m_) <= max_slope);
+
+  if (isinf (start))
     {
-      Interval w = line->at (i).width_;
-      w.intersect (extent);
+      if (isinf (end))
+       b_ = start_height;
+      else
+       b_ = end_height - m_*end;
+    }
+  else
+    b_ = start_height - m_*start;
+}
 
-      if (extent[LEFT] >= w[RIGHT])
-       break;
+Real 
+Building::height (Real x) const
+{
+  if (isinf (x))
+    return (x > 0) ? end_height_ : start_height_;
+  return m_*x + b_;
+}
 
-      Real my_height = line->at (i).height_;
+Real
+Building::intersection (Building const &other) const
+{
+  return (b_ - other.b_) / (other.m_ - m_);
+}
 
-      if (!w.is_empty ()
-         && w.length () > EPS
-         && d * (my_height - stick_out) < 0)
-       {
-         Interval e1 (line->at (i).width_[LEFT], extent[LEFT]);
-         Interval e3 (extent[RIGHT], line->at (i).width_[RIGHT]);
+void
+Building::leading_part (Real chop, Real h)
+{
+  assert (chop > iv_[LEFT] && chop <= iv_[RIGHT] && !equal (chop, iv_[LEFT]));
+  assert (equal (h, height (chop)));
+  iv_[RIGHT] = chop;
+  end_height_ = h;
+}
 
-         if (!e3.is_empty () && e3.length () > EPS)
-           line->insert (line->begin () + i + 1, Skyline_entry (e3, my_height));
+static void
+skyline_trailing_part (list<Building> *sky, Real x)
+{
+  if (equal (x, sky->front ().iv_[RIGHT]))
+    sky->pop_front ();
+  else
+    assert (x < sky->front ().iv_[RIGHT]);
 
-         line->at (i).height_ = stick_out;
-         line->at (i).width_ = w;
-         if (!e1.is_empty () && e1.length () > EPS)
-           line->insert (line->begin () + i, Skyline_entry (e1, my_height));
-       }
+  if (!sky->empty ())
+    {
+      sky->front ().iv_[LEFT] = x;
+      sky->front ().start_height_ = sky->front ().height (x);
     }
 }
 
+bool
+Building::obstructs (Building const &other) const
+{
+  if (equal (intersection (other), iv_[LEFT]) || equal (start_height_, other.start_height_))
+    return m_ > other.m_ || (m_ == other.m_ && b_ > other.b_);
+  return start_height_ > other.start_height_;
+}
+
 void
-merge_skyline (std::vector<Skyline_entry> *a1,
-              std::vector<Skyline_entry> const &a2,
-              Direction dir)
+Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
+                                list<Building> *const result)
 {
-  for (vsize i = 0; i < a2.size (); i++)
+  while (!s1->empty ())
     {
-      Box b;
-      b[X_AXIS] = a2[i].width_;
-      b[Y_AXIS][dir] = a2[i].height_;
-      b[Y_AXIS][-dir] = dir * infinity_f;
-
-      insert_extent_into_skyline (a1, b, X_AXIS, dir);
+      if (s2->front ().obstructs (s1->front ()))
+       swap (s1, s2);
+
+      Building b = s1->front ();
+      while (s2->front ().iv_[RIGHT] < b.iv_[RIGHT]
+            && s2->front ().end_height_ <= b.height (s2->front ().iv_[RIGHT]) + EPS)
+       s2->pop_front ();
+
+      /* the front of s2 either intersects with b or it ends after b */
+      Real end = infinity_f;
+      Real s2_end_height = s2->front ().end_height_;
+      Real s1_end_height = b.height (s2->front ().iv_[RIGHT]);
+      if (s2_end_height > s1_end_height + EPS)
+       end = b.intersection (s2->front ());
+      end = min (end, b.iv_[RIGHT]);
+      Real height = b.height (end);
+
+      b.leading_part (end, height);
+      result->push_front (b);
+
+      skyline_trailing_part (s1, end);
+      if (!s1->empty ())
+       s1->front ().start_height_ = height;
+      skyline_trailing_part (s2, end);
     }
+  result->reverse ();
 }
 
-std::vector<Skyline_entry>
-empty_skyline (Direction d)
+static void
+empty_skyline (list<Building> *const ret)
 {
-  std::vector<Skyline_entry> skyline;
+  ret->push_front (Building (-infinity_f, -infinity_f, -infinity_f, infinity_f, 0));
+}
 
-  Interval i;
-  i.set_empty ();
-  i.swap ();
-  Skyline_entry e;
-  e.width_ = i;
-  e.height_ = -d * infinity_f;
-  skyline.push_back (e);
-  return skyline;
+static void
+single_skyline (Building const &b, list<Building> *const ret, Real max_slope)
+{
+  if (!isinf (b.iv_[RIGHT]))
+    ret->push_front (Building (b.iv_[RIGHT], b.end_height_, -infinity_f, infinity_f, max_slope));
+  ret->push_front (b);
+  if (!isinf (b.iv_[LEFT]))
+    ret->push_front (Building (-infinity_f, -infinity_f, b.start_height_, b.iv_[LEFT], max_slope));
 }
 
-std::vector<Skyline_entry>
-extents_to_skyline (std::vector<Box> const &extents, Axis a, Direction d)
+void
+Skyline::internal_build_skyline (list<Building> *buildings, list<Building> *const result)
 {
+  vsize size = buildings->size ();
 
-  std::vector<Skyline_entry> skyline = empty_skyline (d);
+  if (size == 0)
+    {
+      empty_skyline (result);
+      return;
+    }
+  else if (size == 1)
+    {
+      single_skyline (buildings->front (), result, max_slope_);
+      return;
+    }
 
-  /*
-    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).)
+  list<Building> right_half;
+  list<Building>::iterator i = buildings->begin ();
 
-    We could do a lot better (n log (n), using a balanced tree) but
-    that seems overkill for now.
-  */
-  for (vsize j = extents.size (); j--;)
-    insert_extent_into_skyline (&skyline, extents[j], a, d);
+  for (vsize s = 0; s < size/2; s++)
+    i++;
+  right_half.splice (right_half.end (), *buildings, i, buildings->end ());
 
-  return skyline;
+  list<Building> right;
+  list<Building> left;
+  internal_build_skyline (&right_half, &right);
+  internal_build_skyline (buildings, &left);
+  internal_merge_skyline (&right, &left, result);
 }
 
-/*
-  minimum distance that can be achieved between baselines. "Clouds" is
-  a skyline pointing down.
+Skyline::Skyline ()
+{
+  max_slope_ = 2;
+  sky_ = UP;
+  empty_skyline (&buildings_);
+}
 
-  This is an O (n) algorithm.
-*/
-Real
-skyline_meshing_distance (std::vector<Skyline_entry> const &buildings,
-                         std::vector<Skyline_entry> const &clouds)
+Skyline::Skyline (Direction sky)
 {
-  int i = buildings.size () -1;
-  int j = clouds.size () -1;
+  max_slope_ = 2;
+  sky_ = sky;
+  empty_skyline (&buildings_);
+}
 
-  Real distance = -infinity_f;
+Skyline::Skyline (vector<Box> const &boxes, Axis a, Direction sky)
+{
+  list<Building> bldgs;
+  sky_ = sky;
+  max_slope_ = 2;
 
-  while (i > 0 || j > 0)
+  for (vsize i = 0; i < boxes.size (); i++)
     {
-      Interval w = buildings[i].width_;
-      w.intersect (clouds[j].width_);
+      Interval iv = boxes[i][a];
+      Real height = sky * boxes[i][other_axis (a)][sky];
+      if (!iv.is_empty () && !isinf (height) && !equal (iv[LEFT], iv[RIGHT]))
+       bldgs.push_front (Building (iv[LEFT], height, height, iv[RIGHT], max_slope_));
+    }
+  internal_build_skyline (&bldgs, &buildings_);
+  assert (is_legal_skyline ());
+}
 
-      if (!w.is_empty ())
-       distance = max (distance, (buildings[i].height_ - clouds[j].height_));
+void
+Skyline::merge (Skyline const &other)
+{
+  assert (sky_ == other.sky_);
 
-      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--;
-    }
+  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_);
+  assert (is_legal_skyline ());
+}
 
-  return distance;
+void
+Skyline::insert (Box const &b, Axis a)
+{
+  list<Building> other_bld;
+  list<Building> my_bld;
+  Interval iv = b[a];
+  Real height = sky_ * b[other_axis (a)][sky_];
+
+  my_bld.splice (my_bld.begin (), buildings_);
+  single_skyline (Building (iv[LEFT], height, height, iv[RIGHT], max_slope_), &other_bld, max_slope_);
+  internal_merge_skyline (&other_bld, &my_bld, &buildings_);
+  assert (is_legal_skyline ());
 }
 
-Skyline_entry::Skyline_entry ()
+void
+Skyline::raise (Real r)
 {
-  height_ = 0.0;
+  list<Building>::iterator end = buildings_.end ();
+  for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
+    {
+      i->start_height_ += sky_ * r;
+      i->end_height_ += sky_ * r;
+      i->b_ += sky_ * r;
+    }
+  assert (is_legal_skyline ());
 }
 
-Skyline_entry::Skyline_entry (Interval i, Real r)
+Real
+Skyline::distance (Skyline const &other) const
 {
-  width_ = i;
-  height_ = r;
+  assert (sky_ == -other.sky_);
+  list<Building>::const_iterator i = buildings_.begin ();
+  list<Building>::const_iterator j = other.buildings_.begin ();
+
+  Real dist = -infinity_f;
+  while (i != buildings_.end () && j != 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++;
+      else
+       j++;
+    }
+  return dist;
 }
 
-void
-heighten_skyline (std::vector<Skyline_entry> *buildings, Real ground)
+Real
+Skyline::height (Real airplane) const
 {
-  for (vsize i = 0; i < buildings->size (); i++)
-    buildings->at (i).height_ += ground;
+  assert (!isinf (airplane));
+
+  list<Building>::const_iterator i;
+  for (i = buildings_.begin (); i != buildings_.end (); i++)
+    {
+      if (i->iv_[RIGHT] >= airplane)
+       return sky_ * i->height (airplane);
+    }
+  assert (0);
+  return 0;
 }
 
 Real
-skyline_height (std::vector<Skyline_entry> const &buildings,
-               Real airplane,
-               Direction sky_dir)
+Skyline::max_height () const
 {
-  Real h = - sky_dir * infinity_f;
+  Skyline s (-sky_);
+  s.set_minimum_height (0);
+  return sky_ * distance (s);
+}
 
-  /*
-    Ugh! linear, should be O(log n).
-   */
-  for (vsize i = 0; i < buildings.size (); i++)
-    if (buildings[i].width_.contains (airplane))
-      h = sky_dir * max (sky_dir * h,
-                        sky_dir * buildings[i].height_);
-  
-  return h;
+void
+Skyline::set_minimum_height (Real h)
+{
+  Skyline s (sky_);
+  s.buildings_.front ().start_height_ = h*sky_;
+  s.buildings_.front ().end_height_ = h*sky_;
+  s.buildings_.front ().b_ = h*sky_;
+  merge (s);
 }
 
+Stencil
+Skyline::stencil ()
+{
+  Stencil ret;
+  for (list<Building>::iterator i = buildings_.begin (); i != buildings_.end (); i++)
+    {
+      if (!isinf (i->iv_.length ()))
+       {
+         Stencil line = Line_interface::make_line (0.1,
+                                                   Offset (i->iv_[LEFT], sky_*i->start_height_),
+                                                   Offset (i->iv_[RIGHT], sky_*i->end_height_));
+         ret.add_stencil (line);
+       }
+    }
+  return ret;
+}