lily/skyline.cc

   1 /*
   2   This file is part of LilyPond, the GNU music typesetter.
   3
   4   Copyright (C) 2006--2012 Joe Neeman <joeneeman@gmail.com>
   5
   6   LilyPond is free software: you can redistribute it and/or modify
   7   it under the terms of the GNU General Public License as published by
   8   the Free Software Foundation, either version 3 of the License, or
   9   (at your option) any later version.
  10
  11   LilyPond is distributed in the hope that it will be useful,
  12   but WITHOUT ANY WARRANTY; without even the implied warranty of
  13   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14   GNU General Public License for more details.
  15
  16   You should have received a copy of the GNU General Public License
  17   along with LilyPond.  If not, see <http://www.gnu.org/licenses/>.
  18 */
  19
  20 #include "skyline.hh"
  21 #include "skyline-pair.hh"
  22 #include <deque>
  23 #include <cstdio>
  24
  25 #include "ly-smobs.icc"
  26
  27 /* A skyline is a sequence of non-overlapping buildings: something like
  28    this:
  29                    _______
  30                   |       \                                 ________
  31                   |        \                       ________/        \
  32         /\        |          \                    /                  \
  33        /  --------             \                 /                    \
  34       /                          \              /                      \
  35      /                             ------------/                        ----
  36    --
  37    Each building has a starting position, and ending position, a starting
  38    height and an ending height.
  39
  40    The following invariants are observed:
  41     - the start of the first building is at -infinity
  42     - the end of the last building is at infinity
  43     - if a building has infinite length (ie. the first and last buildings),
  44       then its starting height and ending height are equal
  45     - the end of one building is the same as the beginning of the next
  46       building
  47
  48    We also allow skylines to point down (the structure is exactly the same,
  49    but we think of the part above the line as being filled with mass and the
  50    part below as being empty). ::distance finds the minimum distance between
  51    an UP skyline and a DOWN skyline.
  52
  53    Note that we store DOWN skylines upside-down. That is, in order to compare
  54    a DOWN skyline with an UP skyline, we need to flip the DOWN skyline first.
  55    This means that the merging routine doesn't need to be aware of direction,
  56    but the distance routine does.
  57 */
  58
  59 /* If we start including very thin buildings, numerical accuracy errors can
  60    arise. Therefore, we ignore all buildings that are less than epsilon wide. */
  61 #define EPS 1e-5
  62
  63 static void
  64 print_buildings (list<Building> const &b)
  65 {
  66   for (list<Building>::const_iterator i = b.begin (); i != b.end (); i++)
  67     i->print ();
  68 }
  69
  70 void
  71 Skyline::print () const
  72 {
  73   print_buildings (buildings_);
  74 }
  75
  76 void
  77 Skyline::print_points () const
  78 {
  79   vector<Offset> ps (to_points (X_AXIS));
  80
  81   for (vsize i = 0; i < ps.size (); i++)
  82     printf ("(%f,%f)%s", ps[i][X_AXIS], ps[i][Y_AXIS],
  83             (i % 2) == 1 ? "\n" : " ");
  84 }
  85
  86 Building::Building (Real start, Real start_height, Real end_height, Real end)
  87 {
  88   if (isinf (start) || isinf (end))
  89     assert (start_height == end_height);
  90
  91   start_ = start;
  92   end_ = end;
  93   precompute (start, start_height, end_height, end);
  94 }
  95
  96 Building::Building (Box const &b, Axis horizon_axis, Direction sky)
  97 {
  98   Real start = b[horizon_axis][LEFT];
  99   Real end = b[horizon_axis][RIGHT];
 100   Real height = sky * b[other_axis (horizon_axis)][sky];
 101
 102   start_ = start;
 103   end_ = end;
 104   precompute (start, height, height, end);
 105 }
 106
 107 void
 108 Building::precompute (Real start, Real start_height, Real end_height, Real end)
 109 {
 110   slope_ = 0.0; /* if they were both infinite, we would get nan, not 0, from the prev line */
 111   if (start_height != end_height)
 112     slope_ = (end_height - start_height) / (end - start);
 113
 114   assert (!isinf (slope_) && !isnan (slope_));
 115
 116   if (isinf (start))
 117     {
 118       assert (start_height == end_height);
 119       y_intercept_ = start_height;
 120     }
 121   else
 122     y_intercept_ = start_height - slope_ * start;
 123 }
 124
 125 inline Real
 126 Building::height (Real x) const
 127 {
 128   return isinf (x) ? y_intercept_ : slope_ * x + y_intercept_;
 129 }
 130
 131 void
 132 Building::print () const
 133 {
 134   printf ("%f x + %f from %f to %f\n", slope_, y_intercept_, start_, end_);
 135 }
 136
 137 inline Real
 138 Building::intersection_x (Building const &other) const
 139 {
 140   Real ret = (y_intercept_ - other.y_intercept_) / (other.slope_ - slope_);
 141   return isnan (ret) ? -infinity_f : ret;
 142 }
 143
 144 void
 145 Building::leading_part (Real chop)
 146 {
 147   assert (chop <= end_);
 148   end_ = chop;
 149 }
 150
 151 // Returns a shift s such that (x + s, y) intersects the roof of
 152 // this building.  If no such shift exists, returns infinity_f.
 153 Real
 154 Building::shift_to_intersect (Real x, Real y) const
 155 {
 156   // Solve for s: y = (x + s)*m + b
 157   Real ret = (y - y_intercept_ - slope_ * x) / slope_;
 158
 159   if (ret >= start_ && ret <= end_ && !isinf (ret))
 160     return ret;
 161   return infinity_f;
 162 }
 163
 164 // Returns the interval of horizontal shifts for which this
 165 // building (pointing up) overlaps the other building (pointing down).
 166 Interval
 167 Building::overlapping_shift_interval (Building const &other) const
 168 {
 169   Interval iv;
 170
 171   // If one building is empty, there will never be an overlap.
 172   if (y_intercept_ == -infinity_f || other.y_intercept_ == -infinity_f)
 173     return iv;
 174
 175   // There are two kinds of interesting positions:
 176   // - when the horizontal extents of the buildings just touch
 177   // - when an endpoint of one building intersects the roof of the other.
 178   // The interval we are looking for is the smallest one that
 179   // contains all of the interesting points.
 180
 181
 182   Real my_y1 = height (start_);
 183   Real my_y2 = height (end_);
 184   Real his_y1 = -other.height (other.start_); // "-" because OTHER points down
 185   Real his_y2 = -other.height (other.end_);
 186
 187   // If both buildings are infinite in the same direction,
 188   // the return value is either empty or full.
 189   if ((isinf (start_) && isinf (other.start_))
 190       || (isinf (end_) && isinf (other.end_)))
 191     return (y_intercept_ > other.y_intercept_)
 192            ? Interval (-infinity_f, infinity_f) : Interval ();
 193
 194   // ...when the horizontal extents of the buildings just touch...
 195   if (my_y1 >= his_y2)
 196     iv.add_point (other.end_ - start_);
 197   if (my_y2 >= his_y1)
 198     iv.add_point (other.start_ - end_);
 199
 200   // ...when an endpoint of one building intersects the roof of the other.
 201   Real p1 = shift_to_intersect (other.start_, his_y1);
 202   Real p2 = shift_to_intersect (other.end_, his_y2);
 203   // "-my_y1" because OTHER points down:
 204   Real p3 = other.shift_to_intersect (start_, -my_y1);
 205   Real p4 = other.shift_to_intersect (end_, -my_y2);
 206   if (!isinf (p1))
 207     iv.add_point (p1);
 208   if (!isinf (p2))
 209     iv.add_point (p2);
 210   if (!isinf (p3))
 211     iv.add_point (p3);
 212   if (!isinf (p4))
 213     iv.add_point (p4);
 214
 215   return iv;
 216 }
 217
 218 static Real
 219 first_intersection (Building const &b, list<Building> *const s, Real start_x)
 220 {
 221   while (!s->empty () && start_x < b.end_)
 222     {
 223       Building c = s->front ();
 224
 225       // conceals and intersection_x involve multiplication and
 226       // division. Avoid that, if we can.
 227       if (c.y_intercept_ == -infinity_f)
 228         {
 229           if (c.end_ > b.end_)
 230             return b.end_;
 231           start_x = c.end_;
 232           s->pop_front ();
 233           continue;
 234         }
 235
 236       if (c.conceals (b, start_x))
 237         return start_x;
 238
 239       Real i = b.intersection_x (c);
 240       if (i > start_x && i <= b.end_ && i <= c.end_)
 241         return i;
 242
 243       start_x = c.end_;
 244       if (b.end_ > c.end_)
 245         s->pop_front ();
 246     }
 247   return b.end_;
 248 }
 249
 250 bool
 251 Building::conceals (Building const &other, Real x) const
 252 {
 253   if (slope_ == other.slope_)
 254     return y_intercept_ > other.y_intercept_;
 255
 256   /* their slopes were not equal, so there is an intersection point */
 257   Real i = intersection_x (other);
 258   return (i <= x && slope_ > other.slope_)
 259          || (i > x && slope_ < other.slope_);
 260 }
 261
 262 // Remove redundant empty buildings from the skyline.
 263 // If there are two adjacent empty buildings, they can be
 264 // turned into one.
 265 void
 266 Skyline::normalize ()
 267 {
 268   bool last_empty = false;
 269   list<Building>::iterator i;
 270   for (i = buildings_.begin (); i != buildings_.end (); i++)
 271     {
 272       if (last_empty && i->y_intercept_ == -infinity_f)
 273         {
 274           list<Building>::iterator last = i;
 275           last--;
 276           last->end_ = i->end_;
 277           buildings_.erase (i);
 278           i = last;
 279         }
 280       last_empty = (i->y_intercept_ == -infinity_f);
 281     }
 282
 283   assert (buildings_.front ().start_ == -infinity_f);
 284   assert (buildings_.back ().end_ == infinity_f);
 285 }
 286
 287 void
 288 Skyline::deholify ()
 289 {
 290   // Since a skyline should always be normalized, we can
 291   // assume that there are never two adjacent empty buildings.
 292   // That is, if center is empty then left and right are not.
 293   list<Building>::iterator left = buildings_.begin ();
 294   list<Building>::iterator center = buildings_.begin ();
 295   list<Building>::iterator right;
 296
 297   for (right = buildings_.begin (); right != buildings_.end (); right++)
 298     {
 299       if (center != buildings_.begin () && center->y_intercept_ == -infinity_f)
 300         {
 301           Real p1 = left->height (left->end_);
 302           Real p2 = right->height (right->start_);
 303           *center = Building (center->start_, p1, p2, center->end_);
 304
 305           left = center;
 306           center = right;
 307         }
 308     }
 309 }
 310
 311 void
 312 Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
 313                                  list<Building> *const result) const
 314 {
 315   if (s1->empty () || s2->empty ())
 316     {
 317       programming_error ("tried to merge an empty skyline");
 318       return;
 319     }
 320
 321   Real x = -infinity_f;
 322   Real last_end = -infinity_f;
 323   while (!s1->empty ())
 324     {
 325       if (s2->front ().conceals (s1->front (), x))
 326         swap (s1, s2);
 327
 328       Building b = s1->front ();
 329       Building c = s2->front ();
 330
 331       // Optimization: if the other skyline is empty at this point,
 332       // we can avoid testing some intersections. Just grab as many
 333       // buildings from s1 as we can, and shove them onto the output.
 334       if (c.y_intercept_ == -infinity_f
 335           && c.end_ >= b.end_)
 336         {
 337           list<Building>::iterator i = s1->begin ();
 338           i++;
 339           while (i != s1->end () && i->end_ <= c.end_)
 340             i++;
 341
 342           s1->front ().start_ = x;
 343           result->splice (result->end (), *s1, s1->begin (), i);
 344           x = result->back ().end_;
 345           last_end = x;
 346           continue;
 347         }
 348
 349       Real end = first_intersection (b, s2, x);
 350       if (s2->empty ())
 351         {
 352           b.start_ = last_end;
 353           result->push_back (b);
 354           break;
 355         }
 356
 357       /* only include buildings wider than epsilon */
 358       if (end > x + EPS)
 359         {
 360           b.leading_part (end);
 361           b.start_ = last_end;
 362           last_end = b.end_;
 363           result->push_back (b);
 364         }
 365
 366       if (end >= s1->front ().end_)
 367         s1->pop_front ();
 368
 369       x = end;
 370     }
 371 }
 372
 373 static void
 374 empty_skyline (list<Building> *const ret)
 375 {
 376   ret->push_front (Building (-infinity_f, -infinity_f, -infinity_f, infinity_f));
 377 }
 378
 379 /*
 380   Given Building 'b', build a skyline containing only that building.
 381 */
 382 static void
 383 single_skyline (Building b, list<Building> *const ret)
 384 {
 385   if (b.end_ > b.start_ + EPS)
 386     {
 387       ret->push_back (Building (-infinity_f, -infinity_f,
 388                                 -infinity_f, b.start_));
 389       ret->push_back (b);
 390       ret->push_back (Building (b.end_, -infinity_f,
 391                                 -infinity_f, infinity_f));
 392     }
 393   else
 394     {
 395       empty_skyline (ret);
 396     }
 397 }
 398
 399 /* remove a non-overlapping set of boxes from BOXES and build a skyline
 400    out of them */
 401 static list<Building>
 402 non_overlapping_skyline (list<Building> *const buildings)
 403 {
 404   list<Building> result;
 405   Real last_end = -infinity_f;
 406   Building last_building (-infinity_f, -infinity_f, -infinity_f, infinity_f);
 407   list<Building>::iterator i = buildings->begin ();
 408   while (i != buildings->end ())
 409     {
 410       Real x1 = i->start_;
 411       Real y1 = i->height (i->start_);
 412       Real x2 = i->end_;
 413       Real y2 = i->height (i->end_);
 414
 415       // Drop buildings that will obviously have no effect.
 416       if (last_building.height (x1) >= y1
 417           && last_building.end_ >= x2
 418           && last_building.height (x2) >= y2)
 419         {
 420           list<Building>::iterator j = i++;
 421           buildings->erase (j);
 422           continue;
 423         }
 424
 425       if (x1 < last_end)
 426         {
 427           i++;
 428           continue;
 429         }
 430
 431       if (x1 > last_end + EPS)
 432         result.push_back (Building (last_end, -infinity_f, -infinity_f, x1));
 433
 434       result.push_back (*i);
 435       last_building = *i;
 436       last_end = i->end_;
 437
 438       list<Building>::iterator j = i++;
 439       buildings->erase (j);
 440     }
 441
 442   if (last_end < infinity_f)
 443     result.push_back (Building (last_end, -infinity_f, -infinity_f, infinity_f));
 444   return result;
 445 }
 446
 447 class LessThanBuilding
 448 {
 449 public:
 450   bool operator () (Building const &b1, Building const &b2)
 451   {
 452     return b1.start_ < b2.start_
 453            || (b1.start_ == b2.start_ && b1.height (b1.start_) > b2.height (b1.start_));
 454   }
 455 };
 456
 457 /**
 458    BUILDINGS is a list of buildings, but they could be overlapping
 459    and in any order.  The returned list of buildings is ordered and non-overlapping.
 460 */
 461 list<Building>
 462 Skyline::internal_build_skyline (list<Building> *buildings) const
 463 {
 464   vsize size = buildings->size ();
 465
 466   if (size == 0)
 467     {
 468       list<Building> result;
 469       empty_skyline (&result);
 470       return result;
 471     }
 472   else if (size == 1)
 473     {
 474       list<Building> result;
 475       single_skyline (buildings->front (), &result);
 476       return result;
 477     }
 478
 479   deque<list<Building> > partials;
 480   buildings->sort (LessThanBuilding ());
 481   while (!buildings->empty ())
 482     partials.push_back (non_overlapping_skyline (buildings));
 483
 484   /* we'd like to say while (partials->size () > 1) but that's O (n).
 485      Instead, we exit in the middle of the loop */
 486   while (!partials.empty ())
 487     {
 488       list<Building> merged;
 489       list<Building> one = partials.front ();
 490       partials.pop_front ();
 491       if (partials.empty ())
 492         return one;
 493
 494       list<Building> two = partials.front ();
 495       partials.pop_front ();
 496       internal_merge_skyline (&one, &two, &merged);
 497       partials.push_back (merged);
 498     }
 499   assert (0);
 500   return list<Building> ();
 501 }
 502
 503 Skyline::Skyline ()
 504 {
 505   sky_ = UP;
 506   empty_skyline (&buildings_);
 507 }
 508
 509 Skyline::Skyline (Skyline const &src)
 510 {
 511   sky_ = src.sky_;
 512
 513   /* doesn't a list's copy constructor do this? -- jneem */
 514   for (list<Building>::const_iterator i = src.buildings_.begin ();
 515        i != src.buildings_.end (); i++)
 516     {
 517       buildings_.push_back (Building ((*i)));
 518     }
 519 }
 520
 521 Skyline::Skyline (Direction sky)
 522 {
 523   sky_ = sky;
 524   empty_skyline (&buildings_);
 525 }
 526
 527 /*
 528   Build skyline from a set of boxes.
 529
 530   Boxes should have fatness in the horizon_axis, otherwise they are ignored.
 531  */
 532 Skyline::Skyline (vector<Box> const &boxes, Axis horizon_axis, Direction sky)
 533 {
 534   list<Building> buildings;
 535   sky_ = sky;
 536
 537   Axis vert_axis = other_axis (horizon_axis);
 538   for (vsize i = 0; i < boxes.size (); i++)
 539     {
 540       Interval iv = boxes[i][horizon_axis];
 541       if (iv.length () > EPS && !boxes[i][vert_axis].is_empty ())
 542         buildings.push_front (Building (boxes[i], horizon_axis, sky));
 543     }
 544
 545   buildings_ = internal_build_skyline (&buildings);
 546   normalize ();
 547 }
 548
 549 /*
 550   build skyline from a set of line segments.
 551
 552   Buildings should have fatness in the horizon_axis, otherwise they are ignored.
 553  */
 554 Skyline::Skyline (vector<Drul_array<Offset> > const &segments, Axis horizon_axis, Direction sky)
 555 {
 556   list<Building> buildings;
 557   sky_ = sky;
 558
 559   for (vsize i = 0; i < segments.size (); i++)
 560     {
 561       Drul_array<Offset> const &seg = segments[i];
 562       Offset left = seg[LEFT];
 563       Offset right = seg[RIGHT];
 564       if (left[horizon_axis] > right[horizon_axis])
 565         swap (left, right);
 566
 567       Real x1 = left[horizon_axis];
 568       Real x2 = right[horizon_axis];
 569       Real y1 = left[other_axis (horizon_axis)] * sky;
 570       Real y2 = right[other_axis (horizon_axis)] * sky;
 571
 572       if (x1 + EPS < x2)
 573         buildings.push_back (Building (x1, y1, y2, x2));
 574     }
 575
 576   buildings_ = internal_build_skyline (&buildings);
 577   normalize ();
 578 }
 579
 580 Skyline::Skyline (vector<Skyline_pair> const &skypairs, Direction sky)
 581 {
 582   sky_ = sky;
 583
 584   deque<Skyline> partials;
 585   for (vsize i = 0; i < skypairs.size (); i++)
 586     partials.push_back (Skyline ((skypairs[i])[sky]));
 587
 588   while (partials.size () > 1)
 589     {
 590       Skyline one = partials.front ();
 591       partials.pop_front ();
 592       Skyline two = partials.front ();
 593       partials.pop_front ();
 594
 595       one.merge (two);
 596       partials.push_back (one);
 597     }
 598
 599   if (partials.size ())
 600     buildings_.swap (partials.front ().buildings_);
 601   else
 602     buildings_.clear ();
 603 }
 604
 605 Skyline::Skyline (Box const &b, Axis horizon_axis, Direction sky)
 606 {
 607   sky_ = sky;
 608   Building front (b, horizon_axis, sky);
 609   single_skyline (front, &buildings_);
 610 }
 611
 612 void
 613 Skyline::merge (Skyline const &other)
 614 {
 615   assert (sky_ == other.sky_);
 616
 617   if (other.is_empty ())
 618     return;
 619
 620   if (is_empty ())
 621     {
 622       buildings_ = other.buildings_;
 623       return;
 624     }
 625
 626   list<Building> other_bld (other.buildings_);
 627   list<Building> my_bld;
 628   my_bld.splice (my_bld.begin (), buildings_);
 629   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 630   normalize ();
 631 }
 632
 633 void
 634 Skyline::insert (Box const &b, Axis a)
 635 {
 636   list<Building> other_bld;
 637   list<Building> my_bld;
 638
 639   if (isnan (b[other_axis (a)][LEFT])
 640       || isnan (b[other_axis (a)][RIGHT]))
 641     {
 642       programming_error ("insane box for skyline");
 643       return;
 644     }
 645
 646   /* do the same filtering as in Skyline (vector<Box> const&, etc.) */
 647   Interval iv = b[a];
 648   if (iv.length () <= EPS || b[other_axis (a)].is_empty ())
 649     return;
 650
 651   my_bld.splice (my_bld.begin (), buildings_);
 652   single_skyline (Building (b, a, sky_), &other_bld);
 653   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 654   normalize ();
 655 }
 656
 657 void
 658 Skyline::raise (Real r)
 659 {
 660   list<Building>::iterator end = buildings_.end ();
 661   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
 662     i->y_intercept_ += sky_ * r;
 663 }
 664
 665 void
 666 Skyline::shift (Real s)
 667 {
 668   list<Building>::iterator end = buildings_.end ();
 669   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
 670     {
 671       i->start_ += s;
 672       i->end_ += s;
 673       i->y_intercept_ -= s * i->slope_;
 674     }
 675 }
 676
 677 Real
 678 Skyline::distance (Skyline const &other, Real horizon_padding) const
 679 {
 680   Real dummy;
 681   return internal_distance (other, horizon_padding, &dummy);
 682 }
 683
 684 Real
 685 Skyline::touching_point (Skyline const &other, Real horizon_padding) const
 686 {
 687   Real touch;
 688   internal_distance (other, horizon_padding, &touch);
 689   return touch;
 690 }
 691
 692 Real
 693 Skyline::internal_distance (Skyline const &other, Real horizon_padding, Real *touch_point) const
 694 {
 695   if (horizon_padding == 0.0)
 696     return internal_distance (other, touch_point);
 697
 698   // Note that it is not necessary to build a padded version of other,
 699   // because the same effect can be achieved just by doubling horizon_padding.
 700   Skyline padded_this = padded (horizon_padding);
 701   return padded_this.internal_distance (other, touch_point);
 702 }
 703
 704 Skyline
 705 Skyline::padded (Real horizon_padding) const
 706 {
 707   list<Building> pad_buildings;
 708   for (list<Building>::const_iterator i = buildings_.begin (); i != buildings_.end (); ++i)
 709     {
 710       if (i->start_ > -infinity_f)
 711         {
 712           Real height = i->height (i->start_);
 713           if (height > -infinity_f)
 714             {
 715               // Add the sloped building that pads the left side of the current building.
 716               Real start = i->start_ - 2 * horizon_padding;
 717               Real end = i->start_ - horizon_padding;
 718               pad_buildings.push_back (Building (start, height - horizon_padding, height, end));
 719
 720               // Add the flat building that pads the left side of the current building.
 721               start = i->start_ - horizon_padding;
 722               end = i->start_;
 723               pad_buildings.push_back (Building (start, height, height, end));
 724             }
 725         }
 726
 727       if (i->end_ < infinity_f)
 728         {
 729           Real height = i->height (i->end_);
 730           if (height > -infinity_f)
 731             {
 732               // Add the flat building that pads the right side of the current building.
 733               Real start = i->end_;
 734               Real end = start + horizon_padding;
 735               pad_buildings.push_back (Building (start, height, height, end));
 736
 737               // Add the sloped building that pads the right side of the current building.
 738               start = end;
 739               end += horizon_padding;
 740               pad_buildings.push_back (Building (start, height, height - horizon_padding, end));
 741             }
 742         }
 743     }
 744
 745   // The buildings may be overlapping, so resolve that.
 746   list<Building> pad_skyline = internal_build_skyline (&pad_buildings);
 747
 748   // Merge the padding with the original, to make a new skyline.
 749   Skyline padded (sky_);
 750   list<Building> my_buildings = buildings_;
 751   padded.buildings_.clear ();
 752   internal_merge_skyline (&pad_skyline, &my_buildings, &padded.buildings_);
 753   padded.normalize ();
 754
 755   return padded;
 756 }
 757
 758 Real
 759 Skyline::internal_distance (Skyline const &other, Real *touch_point) const
 760 {
 761   assert (sky_ == -other.sky_);
 762
 763   list<Building>::const_iterator i = buildings_.begin ();
 764   list<Building>::const_iterator j = other.buildings_.begin ();
 765
 766   Real dist = -infinity_f;
 767   Real start = -infinity_f;
 768   Real touch = -infinity_f;
 769   while (i != buildings_.end () && j != other.buildings_.end ())
 770     {
 771       Real end = min (i->end_, j->end_);
 772       Real start_dist = i->height (start) + j->height (start);
 773       Real end_dist = i->height (end) + j->height (end);
 774       dist = max (dist, max (start_dist, end_dist));
 775
 776       if (end_dist == dist)
 777         touch = end;
 778       else if (start_dist == dist)
 779         touch = start;
 780
 781       if (i->end_ <= j->end_)
 782         i++;
 783       else
 784         j++;
 785       start = end;
 786     }
 787
 788   *touch_point = touch;
 789   return dist;
 790 }
 791
 792 Real
 793 Skyline::height (Real airplane) const
 794 {
 795   assert (!isinf (airplane));
 796
 797   list<Building>::const_iterator i;
 798   for (i = buildings_.begin (); i != buildings_.end (); i++)
 799     {
 800       if (i->end_ >= airplane)
 801         return sky_ * i->height (airplane);
 802     }
 803
 804   assert (0);
 805   return 0;
 806 }
 807
 808 Real
 809 Skyline::max_height () const
 810 {
 811   Real ret = -infinity_f;
 812
 813   list<Building>::const_iterator i;
 814   for (i = buildings_.begin (); i != buildings_.end (); ++i)
 815     {
 816       ret = max (ret, i->height (i->start_));
 817       ret = max (ret, i->height (i->end_));
 818     }
 819
 820   return sky_ * ret;
 821 }
 822
 823 Real
 824 Skyline::max_height_position () const
 825 {
 826   Skyline s (-sky_);
 827   s.set_minimum_height (0);
 828   return touching_point (s);
 829 }
 830
 831 void
 832 Skyline::set_minimum_height (Real h)
 833 {
 834   Skyline s (sky_);
 835   s.buildings_.front ().y_intercept_ = h * sky_;
 836   merge (s);
 837 }
 838
 839 vector<Offset>
 840 Skyline::to_points (Axis horizon_axis) const
 841 {
 842   vector<Offset> out;
 843
 844   Real start = -infinity_f;
 845   for (list<Building>::const_iterator i (buildings_.begin ());
 846        i != buildings_.end (); i++)
 847     {
 848       out.push_back (Offset (start, sky_ * i->height (start)));
 849       out.push_back (Offset (i->end_, sky_ * i->height (i->end_)));
 850       start = i->end_;
 851     }
 852
 853   if (horizon_axis == Y_AXIS)
 854     for (vsize i = 0; i < out.size (); i++)
 855       out[i] = out[i].swapped ();
 856
 857   return out;
 858 }
 859
 860 bool
 861 Skyline::is_empty () const
 862 {
 863   if (!buildings_.size ())
 864     return true;
 865   Building b = buildings_.front ();
 866   return b.end_ == infinity_f && b.y_intercept_ == -infinity_f;
 867 }
 868
 869 void
 870 Skyline::clear ()
 871 {
 872   buildings_.clear ();
 873   empty_skyline (&buildings_);
 874 }
 875
 876 /****************************************************************/
 877
 878 IMPLEMENT_SIMPLE_SMOBS (Skyline);
 879 IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
 880 IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
 881
 882 SCM
 883 Skyline::mark_smob (SCM s)
 884 {
 885   ASSERT_LIVE_IS_ALLOWED (s);
 886   return SCM_EOL;
 887 }
 888
 889 int
 890 Skyline::print_smob (SCM s, SCM port, scm_print_state *)
 891 {
 892   Skyline *r = (Skyline *) SCM_CELL_WORD_1 (s);
 893   (void) r;
 894
 895   scm_puts ("#<Skyline>", port);
 896
 897   return 1;
 898 }
 899
 900 MAKE_SCHEME_CALLBACK_WITH_OPTARGS (Skyline, get_touching_point, 3, 1, "")
 901 SCM
 902 Skyline::get_touching_point (SCM skyline_scm, SCM other_skyline_scm, SCM horizon_padding_scm)
 903 {
 904   LY_ASSERT_SMOB (Skyline, other_skyline_scm, 1);
 905
 906   Real horizon_padding = 0;
 907   if (horizon_padding_scm != SCM_UNDEFINED)
 908     {
 909       LY_ASSERT_TYPE (scm_is_number, horizon_padding_scm, 3);
 910       horizon_padding = scm_to_double (horizon_padding_scm);
 911     }
 912
 913   Skyline *skyline = Skyline::unsmob (skyline_scm);
 914   Skyline *other_skyline = Skyline::unsmob (other_skyline_scm);
 915   return scm_from_double (skyline->touching_point (*other_skyline, horizon_padding));
 916 }
 917
 918 MAKE_SCHEME_CALLBACK_WITH_OPTARGS (Skyline, get_distance, 3, 1, "")
 919 SCM
 920 Skyline::get_distance (SCM skyline_scm, SCM other_skyline_scm, SCM horizon_padding_scm)
 921 {
 922   LY_ASSERT_SMOB (Skyline, other_skyline_scm, 1);
 923
 924   Real horizon_padding = 0;
 925   if (horizon_padding_scm != SCM_UNDEFINED)
 926     {
 927       LY_ASSERT_TYPE (scm_is_number, horizon_padding_scm, 3);
 928       horizon_padding = scm_to_double (horizon_padding_scm);
 929     }
 930
 931   Skyline *skyline = Skyline::unsmob (skyline_scm);
 932   Skyline *other_skyline = Skyline::unsmob (other_skyline_scm);
 933   return scm_from_double (skyline->distance (*other_skyline, horizon_padding));
 934 }
 935
 936 MAKE_SCHEME_CALLBACK (Skyline, get_max_height, 1)
 937 SCM
 938 Skyline::get_max_height (SCM skyline_scm)
 939 {
 940   return scm_from_double (Skyline::unsmob (skyline_scm)->max_height ());
 941 }
 942
 943 MAKE_SCHEME_CALLBACK (Skyline, get_max_height_position, 1)
 944 SCM
 945 Skyline::get_max_height_position (SCM skyline_scm)
 946 {
 947   return scm_from_double (Skyline::unsmob (skyline_scm)->max_height_position ());
 948 }
 949
 950 MAKE_SCHEME_CALLBACK (Skyline, get_height, 2)
 951 SCM
 952 Skyline::get_height (SCM skyline_scm, SCM x_scm)
 953 {
 954   Real x = robust_scm2double (x_scm, 0.0);
 955   return scm_from_double (Skyline::unsmob (skyline_scm)->height (x));
 956 }