]> git.donarmstrong.com Git - lilypond.git/blobdiff - lily/skyline.cc
Merge branch 'master' of git+ssh://jneem@git.sv.gnu.org/srv/git/lilypond
[lilypond.git] / lily / skyline.cc
index 41aaa19039fa9a97eac3f76562e5e83bfca313a4..e61abbaea6bed28b05f7b3a3f69aa532d8952be5 100644 (file)
-/*
-  skyline.cc -- implement Skyline_entry and funcs.
-
-  source file of the GNU LilyPond music typesetter
+/* skyline.cc -- implement the Skyline class
 
-  (c) 2002--2005 Han-Wen Nienhuys <hanwen@cs.uu.nl>
+   source file of the GNU LilyPond music typesetter
+   (c) 2006 Joe Neeman <joeneeman@gmail.com>
 */
 
 #include "skyline.hh"
 
+#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.
+*/
+
 /*
-  A skyline is a shape of the form:
+  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.
 
 
-  This file deals with building such skyline structure, and computing
-  the minimum distance between two opposing skylines.
+  * All the messing around with EPS is very fishy.  There are no
+  complicated numerical algorithms involved, so EPS should not be
+  necessary.
 
+  --hwn
+  
+  
+ */
 
-  Invariants for a skyline:
+#define EPS 1e-10
 
-  skyline[...].width_ forms a partition of the real interval, where
-  the segments are adjacent, and ascending. Hence we have
+static inline bool
+approx_equal (Real x, Real y)
+{
+  return abs (x - y) < EPS || (isinf (x) && isinf (y) && ((x > 0) == (y > 0)));
+}
 
-  skyline.top ().width_[RIGHT] = inf
-  skyline[0].width_[LEFT] = -inf
-*/
+static inline bool
+approx_greater_than (Real x, Real y)
+{
+  return x > y + EPS;
+}
 
-const Real EPS = 1e-12;
+static inline bool
+approx_less_than (Real x, Real y)
+{
+  return x < y - EPS;
+}
 
-/*
-  TODO: avoid unnecessary fragmentation.
+static inline bool
+approx_less_equal (Real x, Real y)
+{
+  return x <= y + EPS;
+}
+
+static inline bool
+approx_greater_equal (Real x, Real y)
+{
+  return x >= y - EPS;
+}
+
+void
+Skyline::print () const
+{
+  for (list<Building>::const_iterator i = buildings_.begin ();
+       i != buildings_.end (); i++)
+    {
+      (*i).print ();
+    }
+}
+
+bool
+Skyline::is_legal_skyline () 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;
+}
+
+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 ();
+}
+
+Building::Building (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
+{
+  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 ();
+}
 
-  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::precompute ()
 {
-  Interval extent = b[line_axis];
-  if (extent.is_empty ())
-    return;
+  slope_ = (height_[RIGHT] - height_[LEFT]) / (iv_.length());
+  if (height_[LEFT] == height_[RIGHT]) /* in case they're both infinity */
+    slope_ = 0;
 
-  Real stick_out = b[other_axis (line_axis)][d];
+  assert (!isinf (slope_) && !isnan (slope_));
 
-  /*
-    Intersect each segment of LINE with EXTENT, and if non-empty, insert relevant segments.
-  */
-  for (int i = line->size (); i--;)
+  if (isinf (iv_[START]))
     {
-      Interval w = line->elem (i).width_;
-      w.intersect (extent);
+      assert (slope_ == 0);
+      y_intercept_ = height_[LEFT];
+    }
+  else
+    y_intercept_ = height_[LEFT] - slope_ * iv_[START];
+}
+
+Real 
+Building::height (Real x) const
+{
+  if (isinf (x))
+    return (x > 0) ? height_[RIGHT] : height_[LEFT];
+  return slope_*x + y_intercept_;
+}
+
+void
+Building::print () const
+{
+  printf ("X[%f,%f] -> Y[%f,%f]\n",
+         iv_[LEFT], iv_[RIGHT],
+         height_[LEFT], height_[RIGHT]);
+}
 
-      if (extent[LEFT] >= w[RIGHT])
-       break;
+Real
+Building::intersection_x (Building const &other) const
+{
+  return (y_intercept_ - other.y_intercept_) / (other.slope_ - slope_);
+}
 
-      Real my_height = line->elem (i).height_;
+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);
+}
 
-      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]);
+Building
+Building::sloped_neighbour (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;
+  if (d == LEFT)
+    {
+      swap (left, right);
+      swap (left_height, right_height);
+    }
+  return Building (left, left_height, right_height, right);
+}
+
+bool
+Building::sane () const
+{
+  return approx_less_than (iv_[LEFT], iv_[RIGHT])
+    && !isinf (height_[RIGHT])
+    && !isinf (height_[LEFT]);
+}
 
-         if (!e3.is_empty () && e3.length () > EPS)
-           line->insert (Skyline_entry (e3, my_height), i + 1);
+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]);
 
-         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);
-       }
+  if (!sky->empty ())
+    {
+      sky->front ().iv_[LEFT] = x;
+      sky->front ().height_[LEFT] = sky->front ().height (x);
     }
 }
 
+bool
+Building::conceals_beginning (Building const &other) 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;
+}
+
+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]);
+}
+
 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)
 {
-  for (int 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;
+      if (s2->front ().conceals_beginning (s1->front ()))
+       swap (s1, s2);
+
+      Building b = s1->front ();
+      while (!s2->empty () && b.conceals (s2->front ()))
+       s2->pop_front ();
+      if (s2->empty ())
+       {
+         result->push_front (b);
+         break;
+       }
 
-      insert_extent_into_skyline (a1, b, X_AXIS, dir);
+      /* 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->reverse ();
 }
 
-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;
+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])
+    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]))
+    ret->push_front (Building (-infinity_f, -infinity_f,
+                              -infinity_f, b.iv_[LEFT]));
 }
 
-Array<Skyline_entry>
-extents_to_skyline (Array<Box> const &extents, Axis a, Direction d)
+void
+Skyline::internal_build_skyline (list<Building> *buildings, list<Building> *const result)
 {
+  vsize size = buildings->size ();
 
-  Array<Skyline_entry> skyline = empty_skyline (d);
+  if (size == 0)
+    {
+      empty_skyline (result);
+      return;
+    }
+  else if (size == 1)
+    {
+      single_skyline (buildings->front (), 0, result);
+      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 (int 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);
+}
+
+Skyline::Skyline ()
+{
+  sky_ = UP;
+  empty_skyline (&buildings_);  
+}
+
+Skyline::Skyline (Skyline const &src)
+{
+  sky_ = src.sky_;
+  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;
+  empty_skyline (&buildings_);
 }
 
 /*
-  minimum distance that can be achieved between baselines. "Clouds" is
-  a skyline pointing down.
+  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.
+
+  Boxes should have fatness in the horizon_axis (after they are expanded by
+  horizon_padding), otherwise they are ignored.
+ */
+Skyline::Skyline (vector<Box> const &boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
+{
+  list<Building> bldgs;
+  sky_ = sky;
+
+  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));
+           }
+       }
+    }
+  
+  internal_build_skyline (&bldgs, &buildings_);
+  assert (is_legal_skyline ());
+}
+
+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_);
+}
+
+void
+Skyline::merge (Skyline const &other)
+{
+  assert (sky_ == other.sky_);
+
+  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 ());
+}
+
+void
+Skyline::insert (Box const &b, Real horizon_padding, Axis a)
+{
+  list<Building> other_bld;
+  list<Building> my_bld;
+
+  my_bld.splice (my_bld.begin (), buildings_);
+  single_skyline (Building (b, 0, a, sky_), horizon_padding, &other_bld);
+  internal_merge_skyline (&other_bld, &my_bld, &buildings_);
+  assert (is_legal_skyline ());
+}
+
+void
+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 ());
+}
+
+void
+Skyline::shift (Real r)
+{
+  list<Building>::iterator end = buildings_.end ();
+  for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
+    {
+      i->iv_[LEFT] += r;
+      i->iv_[RIGHT] += r;
+    }
+}
 
-  This is an O (n) algorithm.
-*/
 Real
-skyline_meshing_distance (Array<Skyline_entry> const &buildings,
-                         Array<Skyline_entry> const &clouds)
+Skyline::distance (Skyline const &other) 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 distance = -infinity_f;
+  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;
+}
+
+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->iv_[RIGHT] >= airplane)
+       return sky_ * i->height (airplane);
+    }
+
+  assert (0);
+  return 0;
+}
+
+Real
+Skyline::max_height () const
+{
+  Skyline s (-sky_);
+  s.set_minimum_height (0);
+  return sky_ * distance (s);
+}
 
-      if (!w.is_empty ())
-       distance = max (distance, (buildings[i].height_ - clouds[j].height_));
+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
+{
+  vector<Offset> out;
 
-      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--;
+  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]));
     }
+  return out;
+}
 
-  return distance;
+Skyline_pair::Skyline_pair ()
+  : skylines_ (Skyline (DOWN), Skyline (UP))
+{
 }
 
-Skyline_entry::Skyline_entry ()
+Skyline_pair::Skyline_pair (vector<Box> const &boxes, Real padding, Axis a)
+  : skylines_ (Skyline (boxes, padding, a, DOWN), Skyline (boxes, padding, a, UP))
 {
-  height_ = 0.0;
 }
 
-Skyline_entry::Skyline_entry (Interval i, Real r)
+Skyline_pair::Skyline_pair (Box const &b, Real padding, Axis a)
+  : skylines_ (Skyline (b, padding, a, DOWN), Skyline (b, padding, a, UP))
 {
-  width_ = i;
-  height_ = r;
 }
 
 void
-heighten_skyline (Array<Skyline_entry> *buildings, Real ground)
+Skyline_pair::raise (Real r)
+{
+  skylines_[UP].raise (r);
+  skylines_[DOWN].raise (r);
+}
+
+void
+Skyline_pair::shift (Real r)
+{
+  skylines_[UP].shift (r);
+  skylines_[DOWN].shift (r);
+}
+
+void
+Skyline_pair::insert (Box const &b, Real padding, Axis a)
+{
+  skylines_[UP].insert (b, padding, a);
+  skylines_[DOWN].insert (b, padding, a);
+}
+
+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
+{
+  return skylines_[d];
+}
+
+/****************************************************************/
+
+
+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)
 {
-  for (int i = 0; i < buildings->size (); i++)
-    buildings->elem_ref (i).height_ += ground;
+  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;
+}
+
+SCM
+Skyline_pair::mark_smob (SCM)
+{
+  return SCM_EOL;
+}
+
+int
+Skyline_pair::print_smob (SCM s, SCM port, scm_print_state *)
+{
+  Skyline_pair *r = (Skyline_pair *) SCM_CELL_WORD_1 (s);
+  (void) r;
+
+  scm_puts ("#<Skyline-pair>", port);
+  return 1;
 }