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 #define EPS 1e-10
  45
  46 static inline bool
  47 equal (Real x, Real y)
  48 {
  49   return abs (x - y) < EPS || (isinf (x) && isinf (y) && ((x > 0) == (y > 0)));
  50 }
  51
  52 void
  53 Skyline::print () const
  54 {
  55   for (list<Building>::const_iterator i = buildings_.begin ();
  56        i != buildings_.end (); i++)
  57     {
  58       (*i).print ();
  59     }
  60 }
  61
  62 bool
  63 Skyline::is_legal_skyline () const
  64 {
  65   list<Building>::const_iterator i;
  66   Real last_x = -infinity_f;
  67   Real last_h = -infinity_f;
  68   for (i = buildings_.begin (); i != buildings_.end (); i++)
  69     {
  70       if (i->iv_[LEFT] != last_x)
  71         return false;
  72       if (i != buildings_.begin () && !equal (i->height_[LEFT], last_h))
  73         return false;
  74       last_x = i->iv_[RIGHT];
  75       last_h = i->height_[RIGHT];
  76     }
  77   return last_x == infinity_f;
  78 }
  79
  80 Building::Building (Real start, Real start_height, Real end_height,
  81                     Real end, Real max_slope)
  82   : iv_ (start, end)
  83 {
  84   height_[LEFT] = start_height;
  85   height_[RIGHT] = end_height;
  86
  87   if (isinf (start))
  88     assert (isinf (start_height) || start_height == end_height);
  89   if (isinf (end))
  90     assert (isinf (end_height) || start_height == end_height);
  91
  92   precompute (max_slope);
  93 }
  94
  95 void
  96 Building::precompute (Real max_slope)
  97 {
  98   slope_ = (height_[RIGHT] - height_[LEFT]) / (iv_.length());
  99   if (height_[LEFT] == height_[RIGHT])
 100     slope_ = 0;
 101   if (isinf (slope_) || isnan (slope_))
 102     slope_ = max_slope * (height_[LEFT] < height_[RIGHT] ? 1 : -1);
 103
 104 #if 0
 105   /*
 106     this check is sensitive to roundoff errors when converting to/from
 107     sequences of points.
 108    */
 109   assert (abs (slope_) <= max_slope + EPS);
 110 #endif
 111
 112   if (isinf (iv_[START]))
 113     {
 114       if (isinf (iv_[STOP]))
 115         zero_height_ = height_[LEFT];
 116       else
 117         zero_height_ = height_[RIGHT] - slope_ * iv_[STOP];
 118     }
 119   else
 120     zero_height_ = height_[LEFT] - slope_ * iv_[START];
 121 }
 122
 123 Real
 124 Building::height (Real x) const
 125 {
 126   if (isinf (x))
 127     return (x > 0) ? height_[RIGHT] : height_[LEFT];
 128   return slope_*x + zero_height_;
 129 }
 130
 131 void
 132 Building::print () const
 133 {
 134   printf ("X[%f,%f] -> Y[%f,%f]\n",
 135           iv_[LEFT], iv_[RIGHT],
 136           height_[LEFT], height_[RIGHT]);
 137 }
 138
 139 Real
 140 Building::intersection (Building const &other) const
 141 {
 142   return (zero_height_ - other.zero_height_) / (other.slope_ - slope_);
 143 }
 144
 145 void
 146 Building::leading_part (Real chop, Real h)
 147 {
 148   assert (chop > iv_[LEFT] && chop <= iv_[RIGHT] && !equal (chop, iv_[LEFT]));
 149   assert (equal (h, height (chop)));
 150   iv_[RIGHT] = chop;
 151   height_[RIGHT] = h;
 152 }
 153
 154 static void
 155 skyline_trailing_part (list<Building> *sky, Real x)
 156 {
 157   if (equal (x, sky->front ().iv_[RIGHT]))
 158     sky->pop_front ();
 159   else
 160     assert (x < sky->front ().iv_[RIGHT]);
 161
 162   if (!sky->empty ())
 163     {
 164       sky->front ().iv_[LEFT] = x;
 165       sky->front ().height_[LEFT] = sky->front ().height (x);
 166     }
 167 }
 168
 169 bool
 170 Building::obstructs (Building const &other) const
 171 {
 172   if (equal (intersection (other), iv_[LEFT]) || equal (height_[LEFT], other.height_[LEFT]))
 173     return slope_ > other.slope_ || (slope_ == other.slope_ && zero_height_ > other.zero_height_);
 174   return height_[LEFT] > other.height_[LEFT];
 175 }
 176
 177 void
 178 Skyline::internal_merge_skyline (list<Building> *s1, list<Building> *s2,
 179                                  list<Building> *const result)
 180 {
 181   while (!s1->empty ())
 182     {
 183       if (s2->front ().obstructs (s1->front ()))
 184         swap (s1, s2);
 185
 186       Building b = s1->front ();
 187       while (s2->front ().iv_[RIGHT] < b.iv_[RIGHT]
 188              && s2->front ().height_[RIGHT] <= b.height (s2->front ().iv_[RIGHT]) + EPS)
 189         s2->pop_front ();
 190
 191       /* the front of s2 either intersects with b or it ends after b */
 192       Real end = infinity_f;
 193       Real s2_end_height = s2->front ().height_[RIGHT];
 194       Real s1_end_height = b.height (s2->front ().iv_[RIGHT]);
 195       if (s2_end_height > s1_end_height + EPS)
 196         end = b.intersection (s2->front ());
 197       end = min (end, b.iv_[RIGHT]);
 198       Real height = b.height (end);
 199
 200       b.leading_part (end, height);
 201       result->push_front (b);
 202
 203       skyline_trailing_part (s1, end);
 204       if (!s1->empty ())
 205         s1->front ().height_[LEFT] = height;
 206       skyline_trailing_part (s2, end);
 207     }
 208   result->reverse ();
 209 }
 210
 211 static void
 212 empty_skyline (list<Building> *const ret)
 213 {
 214   ret->push_front (Building (-infinity_f, -infinity_f, -infinity_f, infinity_f, 0));
 215 }
 216
 217 static void
 218 single_skyline (Building const &b, list<Building> *const ret, Real max_slope)
 219 {
 220   if (!isinf (b.iv_[RIGHT]))
 221     ret->push_front (Building (b.iv_[RIGHT], b.height_[RIGHT],
 222                                -infinity_f, infinity_f, max_slope));
 223   if (b.iv_[RIGHT] > b.iv_[LEFT])
 224     ret->push_front (b);
 225   if (!isinf (b.iv_[LEFT]))
 226     ret->push_front (Building (-infinity_f, -infinity_f,
 227                                b.height_[LEFT], b.iv_[LEFT], max_slope));
 228 }
 229
 230 void
 231 Skyline::internal_build_skyline (list<Building> *buildings, list<Building> *const result)
 232 {
 233   vsize size = buildings->size ();
 234
 235   if (size == 0)
 236     {
 237       empty_skyline (result);
 238       return;
 239     }
 240   else if (size == 1)
 241     {
 242       single_skyline (buildings->front (), result, max_slope_);
 243       return;
 244     }
 245
 246   list<Building> right_half;
 247   list<Building>::iterator i = buildings->begin ();
 248
 249   for (vsize s = 0; s < size/2; s++)
 250     i++;
 251   right_half.splice (right_half.end (), *buildings, i, buildings->end ());
 252
 253   list<Building> right;
 254   list<Building> left;
 255   internal_build_skyline (&right_half, &right);
 256   internal_build_skyline (buildings, &left);
 257   internal_merge_skyline (&right, &left, result);
 258 }
 259
 260 Skyline::Skyline ()
 261 {
 262   max_slope_ = 2;
 263   sky_ = UP;
 264   empty_skyline (&buildings_);
 265
 266
 267 }
 268
 269 Skyline::Skyline (Skyline const &src)
 270 {
 271   max_slope_ = src.max_slope_;
 272   sky_ = src.sky_;
 273
 274   for (list<Building>::const_iterator i = src.buildings_.begin ();
 275        i != src.buildings_.end (); i++)
 276     {
 277       buildings_.push_back (Building ((*i)));
 278     }
 279 }
 280
 281 Skyline::Skyline (Direction sky)
 282 {
 283   max_slope_ = 2;
 284   sky_ = sky;
 285   empty_skyline (&buildings_);
 286 }
 287
 288 /*
 289   build skyline from a set of boxes.
 290
 291   Boxes should have fatness in the horizon_axis, otherwise they are ignored.
 292  */
 293 Skyline::Skyline (vector<Box> const &boxes, Axis horizon_axis, Direction sky)
 294 {
 295   list<Building> bldgs;
 296   sky_ = sky;
 297   max_slope_ = 2;
 298
 299   for (vsize i = 0; i < boxes.size (); i++)
 300     {
 301       Interval iv = boxes[i][horizon_axis];
 302       Real height = sky * boxes[i][other_axis (horizon_axis)][sky];
 303       if (!iv.is_empty () && !isinf (height) && !equal (iv[LEFT], iv[RIGHT]))
 304         bldgs.push_front (Building (iv[LEFT], height, height, iv[RIGHT],
 305                                     max_slope_));
 306     }
 307
 308   internal_build_skyline (&bldgs, &buildings_);
 309   assert (is_legal_skyline ());
 310 }
 311
 312 Skyline::Skyline (vector<Offset> const &points, Real max_slope, Direction sky)
 313 {
 314   sky_ = sky;
 315   max_slope_ = max_slope;
 316
 317   for (vsize i = 1; i < points.size (); i++)
 318     {
 319       buildings_.push_back (Building (points[i-1][X_AXIS], sky * points[i-1][Y_AXIS],
 320                                       sky * points[i][Y_AXIS],
 321
 322                                       points[i][X_AXIS],
 323                                       max_slope));
 324
 325     }
 326
 327   assert (is_legal_skyline ());
 328 }
 329
 330 void
 331 Skyline::merge (Skyline const &other)
 332 {
 333   assert (sky_ == other.sky_);
 334
 335   list<Building> other_bld (other.buildings_);
 336   list<Building> my_bld;
 337   my_bld.splice (my_bld.begin (), buildings_);
 338   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 339   assert (is_legal_skyline ());
 340 }
 341
 342 void
 343 Skyline::insert (Box const &b, Axis a)
 344 {
 345   list<Building> other_bld;
 346   list<Building> my_bld;
 347   Interval iv = b[a];
 348   Real height = sky_ * b[other_axis (a)][sky_];
 349
 350   assert (!iv.is_empty ());
 351
 352   my_bld.splice (my_bld.begin (), buildings_);
 353   single_skyline (Building (iv[LEFT], height, height, iv[RIGHT], max_slope_), &other_bld, max_slope_);
 354   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 355   assert (is_legal_skyline ());
 356 }
 357
 358 void
 359 Skyline::raise (Real r)
 360 {
 361   list<Building>::iterator end = buildings_.end ();
 362   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
 363     {
 364       i->height_[LEFT] += sky_ * r;
 365       i->height_[RIGHT] += sky_ * r;
 366       i->zero_height_ += sky_ * r;
 367     }
 368   assert (is_legal_skyline ());
 369 }
 370
 371 Real
 372 Skyline::distance (Skyline const &other) const
 373 {
 374   assert (sky_ == -other.sky_);
 375   list<Building>::const_iterator i = buildings_.begin ();
 376   list<Building>::const_iterator j = other.buildings_.begin ();
 377
 378   Real dist = -infinity_f;
 379   while (i != buildings_.end () && j != other.buildings_.end ())
 380     {
 381       Interval iv = intersection (i->iv_, j->iv_);
 382       dist = max (dist, max (i->height (iv[LEFT]) + j->height (iv[LEFT]),
 383                              i->height (iv[RIGHT]) + j->height (iv[RIGHT])));
 384       if (i->iv_[RIGHT] <= j->iv_[RIGHT])
 385         i++;
 386       else
 387         j++;
 388     }
 389   return dist;
 390 }
 391
 392 Real
 393 Skyline::height (Real airplane) const
 394 {
 395   assert (!isinf (airplane));
 396
 397   list<Building>::const_iterator i;
 398   for (i = buildings_.begin (); i != buildings_.end (); i++)
 399     {
 400       if (i->iv_[RIGHT] >= airplane)
 401         return sky_ * i->height (airplane);
 402     }
 403
 404   assert (0);
 405   return 0;
 406 }
 407
 408 Real
 409 Skyline::max_height () const
 410 {
 411   Skyline s (-sky_);
 412   s.set_minimum_height (0);
 413   return sky_ * distance (s);
 414 }
 415
 416 void
 417 Skyline::set_minimum_height (Real h)
 418 {
 419   Skyline s (sky_);
 420   s.buildings_.front ().height_[LEFT] = h * sky_;
 421   s.buildings_.front ().height_[RIGHT] = h * sky_;
 422   s.buildings_.front ().zero_height_ = h * sky_;
 423   merge (s);
 424 }
 425
 426
 427 vector<Offset>
 428 Skyline::to_points () const
 429 {
 430   vector<Offset> out;
 431
 432   bool first = true;
 433   for (list<Building>::const_iterator i (buildings_.begin ());
 434        i != buildings_.end (); i++)
 435     {
 436       if (first)
 437         out.push_back (Offset ((*i).iv_[LEFT], sky_ * (*i).height_[LEFT]));
 438
 439       first = false;
 440       out.push_back (Offset ((*i).iv_[RIGHT], sky_ * (*i).height_[RIGHT]));
 441     }
 442
 443   return out;
 444 }
 445
 446 /****************************************************************/
 447
 448
 449 IMPLEMENT_SIMPLE_SMOBS (Skyline);
 450 IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
 451 IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
 452
 453 SCM
 454 Skyline::mark_smob (SCM)
 455 {
 456   return SCM_EOL;
 457 }
 458
 459 int
 460 Skyline::print_smob (SCM s, SCM port, scm_print_state *)
 461 {
 462   Skyline *r = (Skyline *) SCM_CELL_WORD_1 (s);
 463   (void) r;
 464
 465   scm_puts ("#<Skyline>", port);
 466
 467   return 1;
 468 }