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, Real max_slope, 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, max_slope, &right);
 256   internal_build_skyline (buildings, max_slope, &left);
 257   internal_merge_skyline (&right, &left, result);
 258 }
 259
 260 Skyline::Skyline ()
 261 {
 262   sky_ = UP;
 263   empty_skyline (&buildings_);
 264 }
 265
 266 Skyline::Skyline (Skyline const &src)
 267 {
 268   sky_ = src.sky_;
 269
 270   for (list<Building>::const_iterator i = src.buildings_.begin ();
 271        i != src.buildings_.end (); i++)
 272     {
 273       buildings_.push_back (Building ((*i)));
 274     }
 275 }
 276
 277 Skyline::Skyline (Direction sky)
 278 {
 279   sky_ = sky;
 280   empty_skyline (&buildings_);
 281 }
 282
 283 /*
 284   build skyline from a set of boxes.
 285
 286   Boxes should have fatness in the horizon_axis, otherwise they are ignored.
 287  */
 288 Skyline::Skyline (vector<Box> const &boxes, Real max_slope, Axis horizon_axis, Direction sky)
 289 {
 290   list<Building> bldgs;
 291   sky_ = sky;
 292
 293   for (vsize i = 0; i < boxes.size (); i++)
 294     {
 295       Interval iv = boxes[i][horizon_axis];
 296       Real height = sky * boxes[i][other_axis (horizon_axis)][sky];
 297       if (!iv.is_empty () && !isinf (height) && !equal (iv[LEFT], iv[RIGHT]))
 298         bldgs.push_front (Building (iv[LEFT], height, height, iv[RIGHT],
 299                                     max_slope));
 300     }
 301
 302   internal_build_skyline (&bldgs, max_slope, &buildings_);
 303   assert (is_legal_skyline ());
 304 }
 305
 306 Skyline::Skyline (vector<Offset> const &points, Real max_slope, Direction sky)
 307 {
 308   sky_ = sky;
 309
 310   for (vsize i = 1; i < points.size (); i++)
 311     {
 312       buildings_.push_back (Building (points[i-1][X_AXIS],
 313                                       sky * points[i-1][Y_AXIS],
 314                                       sky * points[i][Y_AXIS],
 315                                       points[i][X_AXIS],
 316                                       max_slope));
 317
 318     }
 319
 320   assert (is_legal_skyline ());
 321 }
 322
 323 void
 324 Skyline::merge (Skyline const &other)
 325 {
 326   assert (sky_ == other.sky_);
 327
 328   list<Building> other_bld (other.buildings_);
 329   list<Building> my_bld;
 330   my_bld.splice (my_bld.begin (), buildings_);
 331   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 332   assert (is_legal_skyline ());
 333 }
 334
 335 void
 336 Skyline::insert (Box const &b, Real max_slope, Axis a)
 337 {
 338   list<Building> other_bld;
 339   list<Building> my_bld;
 340   Interval iv = b[a];
 341   Real height = sky_ * b[other_axis (a)][sky_];
 342
 343   assert (!iv.is_empty ());
 344
 345   my_bld.splice (my_bld.begin (), buildings_);
 346   single_skyline (Building (iv[LEFT], height, height, iv[RIGHT], max_slope), &other_bld, max_slope);
 347   internal_merge_skyline (&other_bld, &my_bld, &buildings_);
 348   assert (is_legal_skyline ());
 349 }
 350
 351 void
 352 Skyline::raise (Real r)
 353 {
 354   list<Building>::iterator end = buildings_.end ();
 355   for (list<Building>::iterator i = buildings_.begin (); i != end; i++)
 356     {
 357       i->height_[LEFT] += sky_ * r;
 358       i->height_[RIGHT] += sky_ * r;
 359       i->zero_height_ += sky_ * r;
 360     }
 361   assert (is_legal_skyline ());
 362 }
 363
 364 Real
 365 Skyline::distance (Skyline const &other) const
 366 {
 367   assert (sky_ == -other.sky_);
 368   list<Building>::const_iterator i = buildings_.begin ();
 369   list<Building>::const_iterator j = other.buildings_.begin ();
 370
 371   Real dist = -infinity_f;
 372   while (i != buildings_.end () && j != other.buildings_.end ())
 373     {
 374       Interval iv = intersection (i->iv_, j->iv_);
 375       dist = max (dist, max (i->height (iv[LEFT]) + j->height (iv[LEFT]),
 376                              i->height (iv[RIGHT]) + j->height (iv[RIGHT])));
 377       if (i->iv_[RIGHT] <= j->iv_[RIGHT])
 378         i++;
 379       else
 380         j++;
 381     }
 382   return dist;
 383 }
 384
 385 Real
 386 Skyline::height (Real airplane) const
 387 {
 388   assert (!isinf (airplane));
 389
 390   list<Building>::const_iterator i;
 391   for (i = buildings_.begin (); i != buildings_.end (); i++)
 392     {
 393       if (i->iv_[RIGHT] >= airplane)
 394         return sky_ * i->height (airplane);
 395     }
 396
 397   assert (0);
 398   return 0;
 399 }
 400
 401 Real
 402 Skyline::max_height () const
 403 {
 404   Skyline s (-sky_);
 405   s.set_minimum_height (0);
 406   return sky_ * distance (s);
 407 }
 408
 409 void
 410 Skyline::set_minimum_height (Real h)
 411 {
 412   Skyline s (sky_);
 413   s.buildings_.front ().height_[LEFT] = h * sky_;
 414   s.buildings_.front ().height_[RIGHT] = h * sky_;
 415   s.buildings_.front ().zero_height_ = h * sky_;
 416   merge (s);
 417 }
 418
 419
 420 vector<Offset>
 421 Skyline::to_points () const
 422 {
 423   vector<Offset> out;
 424
 425   bool first = true;
 426   for (list<Building>::const_iterator i (buildings_.begin ());
 427        i != buildings_.end (); i++)
 428     {
 429       if (first)
 430         out.push_back (Offset ((*i).iv_[LEFT], sky_ * (*i).height_[LEFT]));
 431
 432       first = false;
 433       out.push_back (Offset ((*i).iv_[RIGHT], sky_ * (*i).height_[RIGHT]));
 434     }
 435
 436   return out;
 437 }
 438
 439 /****************************************************************/
 440
 441
 442 IMPLEMENT_SIMPLE_SMOBS (Skyline);
 443 IMPLEMENT_TYPE_P (Skyline, "ly:skyline?");
 444 IMPLEMENT_DEFAULT_EQUAL_P (Skyline);
 445
 446 SCM
 447 Skyline::mark_smob (SCM)
 448 {
 449   return SCM_EOL;
 450 }
 451
 452 int
 453 Skyline::print_smob (SCM s, SCM port, scm_print_state *)
 454 {
 455   Skyline *r = (Skyline *) SCM_CELL_WORD_1 (s);
 456   (void) r;
 457
 458   scm_puts ("#<Skyline>", port);
 459
 460   return 1;
 461 }