lily/skyline.cc

   1 /* skyline.cc -- implement the Skyline class
   2
   3    source file of the GNU LilyPond music typesetter
   4
   5    (c) 2006 Joe Neeman <joeneeman@gmail.com>
   6 */
   7
   8 #include "skyline.hh"
   9
  10 #include "ly-smobs.icc"
  11
  12 /* A skyline is a sequence of non-overlapping buildings: something like
  13    this:
  14                    _______
  15                   |       \                                 ________
  16                   |        \                       ________/        \
  17         /\        |          \                    /                  \
  18        /  --------             \                 /                    \
  19       /                          \              /                      \
  20      /                             ------------/                        ----
  21    --
  22    Each building has a starting position, and ending position, a starting
  23    height and an ending height.
  24
  25    The following invariants are observed:
  26     - the start of the first building is at -infinity
  27     - the end of the last building is at infinity
  28     - if a building has infinite length (ie. the first and last buildings),
  29       then its starting height and ending height are equal
  30     - the end of one building is the same as the beginning of the next
  31       building
  32
  33    We also allow skylines to point down (the structure is exactly the same,
  34    but we think of the part above the line as being filled with mass and the
  35    part below as being empty). ::distance finds the minimum distance between
  36    an UP skyline and a DOWN skyline.
  37
  38    Note that we store DOWN skylines upside-down. That is, in order to compare
  39    a DOWN skyline with an UP skyline, we need to flip the DOWN skyline first.
  40    This means that the merging routine doesn't need to be aware of direction,
  41    but the distance routine does.
  42 */
  43
  44 /*
  45   FIXME:
  46
  47   * Consider to use
  48
  49   typedef list<Skyline_point> Skyline;
  50   struct Skyline_point
  51   {
  52     Real x;
  53     Drul_array<Real> ys;
  54   };
  55
  56   this is a cleaner representation, as it doesn't duplicate the X, and
  57   doesn't need bogus buildings at infinity  --hwn.
  58
  59
  60   * All the messing around with EPS is very fishy.  There are no
  61   complicated numerical algorithms involved, so EPS should not be
  62   necessary.
  63
  64   --hwn
  65
  66
  67  */
  68
  69 #define EPS 1e-10
  70
  71 static inline bool
  72 approx_equal (Real x, Real y)
  73 {
  74   return abs (x - y) < EPS || (isinf (x) && isinf (y) && ((x > 0) == (y > 0)));
  75 }
  76
  77 static inline bool
  78 approx_greater_than (Real x, Real y)
  79 {
  80   return x > y + EPS;
  81 }
  82
  83 static inline bool
  84 approx_less_than (Real x, Real y)
  85 {
  86   return x < y - EPS;
  87 }
  88
  89 static inline bool
  90 approx_less_equal (Real x, Real y)
  91 {
  92   return x <= y + EPS;
  93 }
  94
  95 static inline bool
  96 approx_greater_equal (Real x, Real y)
  97 {
  98   return x >= y - EPS;
  99 }
 100
 101 void
 102 Skyline::print () const
 103 {
 104   for (list<Building>::const_iterator i = buildings_.begin ();
 105        i != buildings_.end (); i++)
 106     {
 107       (*i).print ();
 108     }
 109 }
 110
 111 bool
 112 Skyline::is_legal_skyline () const
 113 {
 114   list<Building>::const_iterator i;
 115   Real last_x = -infinity_f;
 116   for (i = buildings_.begin (); i != buildings_.end (); i++)
 117     {
 118       if (i->iv_[LEFT] != last_x)
 119         return false;
 120       last_x = i->iv_[RIGHT];
 121       if (isinf (i->iv_.length ()) && i->height_[LEFT] != i->height_[RIGHT])
 122         return false;
 123     }
 124   return last_x == infinity_f;
 125 }
 126
 127 Building::Building (Real start, Real start_height, Real end_height, Real end)
 128   : iv_ (start, end)
 129 {
 130   height_[LEFT] = start_height;
 131   height_[RIGHT] = end_height;
 132
 133   if (isinf (start) || isinf (end))
 134     assert (start_height == end_height);
 135
 136   precompute ();
 137 }
 138
 139 Building::Building (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
 140 {
 141   Real height = sky * b[other_axis (horizon_axis)][sky];
 142
 143   iv_ = b[horizon_axis];
 144   iv_.widen (horizon_padding + EPS);
 145   height_[LEFT] = height;
 146   height_[RIGHT] = height;
 147
 148   if (sane ())
 149     precompute ();
 150 }
 151
 152 void
 153 Building::precompute ()
 154 {
 155   slope_ = (height_[RIGHT] - height_[LEFT]) / (iv_.length());
 156   if (height_[LEFT] == height_[RIGHT]) /* in case they're both infinity */
 157     slope_ = 0;
 158
 159   assert (!isinf (slope_) && !isnan (slope_));
 160
 161   if (isinf (iv_[START]))
 162     {
 163       assert (slope_ == 0);
 164       y_intercept_ = height_[LEFT];
 165     }
 166   else
 167     y_intercept_ = height_[LEFT] - slope_ * iv_[START];
 168 }
 169
 170 Real
 171 Building::height (Real x) const
 172 {
 173   if (isinf (x))
 174     return (x > 0) ? height_[RIGHT] : height_[LEFT];
 175   return slope_*x + y_intercept_;
 176 }
 177
 178 void
 179 Building::print () const
 180 {
 181   printf ("X[%f,%f] -> Y[%f,%f]\n",
 182           iv_[LEFT], iv_[RIGHT],
 183           height_[LEFT], height_[RIGHT]);
 184 }
 185
 186 Real
 187 Building::intersection_x (Building const &other) const
 188 {
 189   return (y_intercept_ - other.y_intercept_) / (other.slope_ - slope_);
 190 }
 191
 192 void
 193 Building::leading_part (Real chop)
 194 {
 195   assert (chop > iv_[LEFT] && chop <= iv_[RIGHT] && !approx_equal (chop, iv_[LEFT]));
 196   iv_[RIGHT] = chop;
 197   height_[RIGHT] = height (chop);
 198 }
 199
 200 Building
 201 Building::sloped_neighbour (Real horizon_padding, Direction d) const
 202 {
 203   Real left = iv_[d];
 204   Real right = iv_[d] + d * horizon_padding;
 205   Real left_height = height_[d];
 206   Real right_height = height_[d] - horizon_padding;
 207   if (d == LEFT)
 208     {
 209       swap (left, right);
 210       swap (left_height, right_height);
 211     }
 212   return Building (left, left_height, right_height, right);
 213 }
 214
 215 bool
 216 Building::sane () const
 217 {
 218   return approx_less_than (iv_[LEFT], iv_[RIGHT])
 219     && !isinf (height_[RIGHT])
 220     && !isinf (height_[LEFT]);
 221 }
 222
 223 static void
 224 skyline_trailing_part (list<Building> *sky, Real x)
 225 {
 226   if (approx_equal (x, sky->front ().iv_[RIGHT]))
 227     sky->pop_front ();
 228   else
 229     assert (x < sky->front ().iv_[RIGHT]);
 230
 231   if (!sky->empty ())
 232     {
 233       sky->front ().iv_[LEFT] = x;
 234       sky->front ().height_[LEFT] = sky->front ().height (x);
 235     }
 236 }
 237
 238 bool
 239 Building::conceals_beginning (Building const &other) const
 240 {
 241   bool w = false;
 242   Real h = other.height (iv_[LEFT]);
 243   if (approx_equal (height_[LEFT], h))
 244     w = slope_ > other.slope_;
 245   else if (height_[LEFT] > h)
 246     w = true;
 247   else
 248     w = false;
 249
 250   return w;
 251 }
 252
 253 bool
 254 Building::conceals (Building const &other) const
 255 {
 256   assert (iv_[LEFT] <= other.iv_[LEFT]);
 257   return (iv_[RIGHT] >= other.iv_[RIGHT])
 258     && approx_greater_equal (height (other.iv_[LEFT]), other.height_[LEFT])
 259     && approx_greater_equal (height (other.iv_[RIGHT]), other.height_[RIGHT]);
 260 }
 261
 262 void
 263 Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
 264                                  list<Building> *const result)
 265 {
 266   while (!s1->empty ())
 267     {
 268       if (s2->front ().conceals_beginning (s1->front ()))
 269         swap (s1, s2);
 270
 271       Building b = s1->front ();
 272       while (!s2->empty () && b.conceals (s2->front ()))
 273         s2->pop_front ();
 274       if (s2->empty ())
 275         {
 276           result->push_front (b);
 277           break;
 278         }
 279
 280       /* s2 either intersects with b or it ends after b */
 281       Real end = infinity_f;
 282       Real s2_start_height = s2->front ().height_[LEFT];
 283       Real s2_end_height = s2->front ().height_[RIGHT];
 284       Real s1_start_height = b.height (s2->front ().iv_[LEFT]);
 285       Real s1_end_height = b.height (s2->front ().iv_[RIGHT]);
 286       if (approx_greater_than (s2_start_height, s1_start_height))
 287         end = s2->front ().iv_[LEFT];
 288       else if (approx_greater_than (s2_end_height, s1_end_height))
 289         end = b.intersection_x (s2->front ());
 290       end = min (end, b.iv_[RIGHT]);
 291
 292       b.leading_part (end);
 293       result->push_front (b);
 294
 295       skyline_trailing_part (s1, end);
 296       skyline_trailing_part (s2, end);
 297     }
 298   result->reverse ();
 299 }
 300
 301 static void
 302 empty_skyline (list<Building> *const ret)
 303 {
 304   ret->push_front (Building (-infinity_f, -infinity_f, -infinity_f, infinity_f));
 305 }
 306
 307 static void
 308 single_skyline (Building b, Real horizon_padding, list<Building> *const ret)
 309 {
 310   b.iv_.widen (horizon_padding);
 311
 312   if (!isinf (b.iv_[RIGHT]))
 313     ret->push_front (Building (b.iv_[RIGHT], -infinity_f,
 314                                -infinity_f, infinity_f));
 315   if (horizon_padding > 0 && !isinf (b.iv_.length ()))
 316     ret->push_front (b.sloped_neighbour (horizon_padding, RIGHT));
 317
 318   if (b.iv_[RIGHT] > b.iv_[LEFT])
 319     ret->push_front (b);
 320
 321   if (horizon_padding > 0 && !isinf (b.iv_.length ()))
 322     ret->push_front (b.sloped_neighbour (horizon_padding, LEFT));
 323   if (!isinf (b.iv_[LEFT]))
 324     ret->push_front (Building (-infinity_f, -infinity_f,
 325                                -infinity_f, b.iv_[LEFT]));
 326 }
 327
 328 void
 329 Skyline::internal_build_skyline (list<Building> *buildings, list<Building> *const result)
 330 {
 331   vsize size = buildings->size ();
 332
 333   if (size == 0)
 334     {
 335       empty_skyline (result);
 336       return;
 337     }
 338   else if (size == 1)
 339     {
 340       single_skyline (buildings->front (), 0, result);
 341       return;
 342     }
 343
 344   list<Building> right_half;
 345   list<Building>::iterator i = buildings->begin ();
 346
 347   for (vsize s = 0; s < size/2; s++)
 348     i++;
 349   right_half.splice (right_half.end (), *buildings, i, buildings->end ());
 350
 351   list<Building> right;
 352   list<Building> left;
 353   internal_build_skyline (&right_half, &right);
 354   internal_build_skyline (buildings, &left);
 355   internal_merge_skyline (&right, &left, result);
 356 }
 357
 358 Skyline::Skyline ()
 359 {
 360   sky_ = UP;
 361   empty_skyline (&buildings_);
 362 }
 363
 364 Skyline::Skyline (Skyline const &src)
 365 {
 366   sky_ = src.sky_;
 367
 368   for (list<Building>::const_iterator i = src.buildings_.begin ();
 369        i != src.buildings_.end (); i++)
 370     {
 371       buildings_.push_back (Building ((*i)));
 372     }
 373 }
 374
 375 Skyline::Skyline (Direction sky)
 376 {
 377   sky_ = sky;
 378   empty_skyline (&buildings_);
 379 }
 380
 381 /*
 382   build skyline from a set of boxes. If horizon_padding > 0, expand all the boxes
 383   by that amount and add 45-degree sloped boxes to the edges of each box (of
 384   width horizon_padding). That is, the total amount of horizontal expansion is
 385   horizon_padding*4, half of which is sloped and half of which is flat.
 386
 387   Boxes should have fatness in the horizon_axis (after they are expanded by
 388   horizon_padding), otherwise they are ignored.
 389  */
 390 Skyline::Skyline (vector<Box> const &boxes, Real horizon_padding, Axis horizon_axis, Direction sky)
 391 {
 392   list<Building> bldgs;
 393   sky_ = sky;
 394
 395   for (vsize i = 0; i < boxes.size (); i++)
 396     {
 397       Building front (boxes[i], horizon_padding, horizon_axis, sky);
 398       if (front.sane ())
 399         {
 400           bldgs.push_front (front);
 401           if (horizon_padding > 0 && !isinf (front.iv_.length ()))
 402             {
 403               bldgs.push_front (front.sloped_neighbour (horizon_padding, LEFT));
 404               bldgs.push_front (front.sloped_neighbour (horizon_padding, RIGHT));
 405             }
 406         }
 407     }
 408
 409   internal_build_skyline (&bldgs, &buildings_);
 410   assert (is_legal_skyline ());
 411 }
 412
 413 Skyline::Skyline (Box const &b, Real horizon_padding, Axis horizon_axis, Direction sky)
 414 {
 415   sky_ = sky;
 416   Building front (b, 0, horizon_axis, sky);
 417   single_skyline (front, horizon_padding, &buildings_);
 418 }
 419
 420 void
 421 Skyline::merge (Skyline const &other)
 422 {
 423   assert (sky_ == other.sky_);
 424
 425   list<Building> other_bld (other.buildings_);
 426   list<Building> my_bld;
 427   my_bld.splice (my_bld.begin (), buildings_);
 428   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 429   assert (is_legal_skyline ());
 430 }
 431
 432 void
 433 Skyline::insert (Box const &b, Real horizon_padding, Axis a)
 434 {
 435   list<Building> other_bld;
 436   list<Building> my_bld;
 437
 438   my_bld.splice (my_bld.begin (), buildings_);
 439   single_skyline (Building (b, 0, a, sky_), horizon_padding, &other_bld);
 440   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 441   assert (is_legal_skyline ());
 442 }
 443
 444 void
 445 Skyline::raise (Real r)
 446 {
 447   list<Building>::iterator end = buildings_.end ();
 448   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
 449     {
 450       i->height_[LEFT] += sky_ * r;
 451       i->height_[RIGHT] += sky_ * r;
 452       i->y_intercept_ += sky_ * r;
 453     }
 454   assert (is_legal_skyline ());
 455 }
 456
 457 void
 458 Skyline::shift (Real r)
 459 {
 460   list<Building>::iterator end = buildings_.end ();
 461   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
 462     {
 463       i->iv_[LEFT] += r;
 464       i->iv_[RIGHT] += r;
 465     }
 466 }
 467
 468 Real
 469 Skyline::distance (Skyline const &other) const
 470 {
 471   assert (sky_ == -other.sky_);
 472   list<Building>::const_iterator i = buildings_.begin ();
 473   list<Building>::const_iterator j = other.buildings_.begin ();
 474
 475   Real dist = -infinity_f;
 476   while (i != buildings_.end () && j != other.buildings_.end ())
 477     {
 478       Interval iv = intersection (i->iv_, j->iv_);
 479       dist = max (dist, max (i->height (iv[LEFT]) + j->height (iv[LEFT]),
 480                              i->height (iv[RIGHT]) + j->height (iv[RIGHT])));
 481       if (i->iv_[RIGHT] <= j->iv_[RIGHT])
 482         i++;
 483       else
 484         j++;
 485     }
 486   return dist;
 487 }
 488
 489 Real
 490 Skyline::height (Real airplane) const
 491 {
 492   assert (!isinf (airplane));
 493
 494   list<Building>::const_iterator i;
 495   for (i = buildings_.begin (); i != buildings_.end (); i++)
 496     {
 497       if (i->iv_[RIGHT] >= airplane)
 498         return sky_ * i->height (airplane);
 499     }
 500
 501   assert (0);
 502   return 0;
 503 }
 504
 505 Real
 506 Skyline::max_height () const
 507 {
 508   Skyline s (-sky_);
 509   s.set_minimum_height (0);
 510   return sky_ * distance (s);
 511 }
 512
 513 void
 514 Skyline::set_minimum_height (Real h)
 515 {
 516   Skyline s (sky_);
 517   s.buildings_.front ().height_[LEFT] = h * sky_;
 518   s.buildings_.front ().height_[RIGHT] = h * sky_;
 519   s.buildings_.front ().y_intercept_ = h * sky_;
 520   merge (s);
 521 }
 522
 523
 524 vector<Offset>
 525 Skyline::to_points () const
 526 {
 527   vector<Offset> out;
 528
 529   for (list<Building>::const_iterator i (buildings_.begin ());
 530        i != buildings_.end (); i++)
 531     {
 532       if (!isinf (i->iv_[LEFT]) && !isinf (i->height_[LEFT]))
 533         out.push_back (Offset (i->iv_[LEFT], sky_ * i->height_[LEFT]));
 534       if (!isinf (i->iv_[RIGHT]) && !isinf (i->height_[RIGHT]))
 535         out.push_back (Offset (i->iv_[RIGHT], sky_ * i->height_[RIGHT]));
 536     }
 537   return out;
 538 }
 539
 540 Skyline_pair::Skyline_pair ()
 541   : skylines_ (Skyline (DOWN), Skyline (UP))
 542 {
 543 }
 544
 545 Skyline_pair::Skyline_pair (vector<Box> const &boxes, Real padding, Axis a)
 546   : skylines_ (Skyline (boxes, padding, a, DOWN), Skyline (boxes, padding, a, UP))
 547 {
 548 }
 549
 550 Skyline_pair::Skyline_pair (Box const &b, Real padding, Axis a)
 551   : skylines_ (Skyline (b, padding, a, DOWN), Skyline (b, padding, a, UP))
 552 {
 553 }
 554
 555 void
 556 Skyline_pair::raise (Real r)
 557 {
 558   skylines_[UP].raise (r);
 559   skylines_[DOWN].raise (r);
 560 }
 561
 562 void
 563 Skyline_pair::shift (Real r)
 564 {
 565   skylines_[UP].shift (r);
 566   skylines_[DOWN].shift (r);
 567 }
 568
 569 void
 570 Skyline_pair::insert (Box const &b, Real padding, Axis a)
 571 {
 572   skylines_[UP].insert (b, padding, a);
 573   skylines_[DOWN].insert (b, padding, a);
 574 }
 575
 576 void
 577 Skyline_pair::merge (Skyline_pair const &other)
 578 {
 579   skylines_[UP].merge (other[UP]);
 580   skylines_[DOWN].merge (other[DOWN]);
 581 }
 582
 583 Skyline&
 584 Skyline_pair::operator [] (Direction d)
 585 {
 586   return skylines_[d];
 587 }
 588
 589 Skyline const&
 590 Skyline_pair::operator [] (Direction d) const
 591 {
 592   return skylines_[d];
 593 }
 594
 595 /****************************************************************/
 596
 597
 598 IMPLEMENT_SIMPLE_SMOBS (Skyline);
 599 IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
 600 IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
 601
 602 IMPLEMENT_SIMPLE_SMOBS (Skyline_pair);
 603 IMPLEMENT_TYPE_P (Skyline_pair, "ly:skyline-pair?");
 604 IMPLEMENT_DEFAULT_EQUAL_P (Skyline_pair);
 605
 606 SCM
 607 Skyline::mark_smob (SCM)
 608 {
 609   return SCM_EOL;
 610 }
 611
 612 int
 613 Skyline::print_smob (SCM s, SCM port, scm_print_state *)
 614 {
 615   Skyline *r = (Skyline *) SCM_CELL_WORD_1 (s);
 616   (void) r;
 617
 618   scm_puts ("#<Skyline>", port);
 619
 620   return 1;
 621 }
 622
 623 SCM
 624 Skyline_pair::mark_smob (SCM)
 625 {
 626   return SCM_EOL;
 627 }
 628
 629 int
 630 Skyline_pair::print_smob (SCM s, SCM port, scm_print_state *)
 631 {
 632   Skyline_pair *r = (Skyline_pair *) SCM_CELL_WORD_1 (s);
 633   (void) r;
 634
 635   scm_puts ("#<Skyline-pair>", port);
 636   return 1;
 637 }