lily/beam-quanting.cc

   1 /*
   2   beam-quanting.cc -- implement Beam quanting functions
   3
   4   source file of the GNU LilyPond music typesetter
   5
   6   (c) 1997--2007 Han-Wen Nienhuys <hanwen@xs4all.nl>
   7   Jan Nieuwenhuizen <janneke@gnu.org>
   8 */
   9
  10 #include "beam.hh"
  11
  12 #include <algorithm>
  13 using namespace std;
  14
  15 #include "grob.hh"
  16 #include "align-interface.hh"
  17 #include "international.hh"
  18 #include "output-def.hh"
  19 #include "pointer-group-interface.hh"
  20 #include "staff-symbol-referencer.hh"
  21 #include "stem.hh"
  22 #include "warn.hh"
  23 #include "main.hh"
  24
  25 Real
  26 get_detail (SCM alist, SCM sym, Real def)
  27 {
  28   SCM entry = scm_assq (sym, alist);
  29
  30   if (scm_is_pair (entry))
  31     return robust_scm2double (scm_cdr (entry), def);
  32   return def;
  33 }
  34
  35 void
  36 Beam_quant_parameters::fill (Grob *him)
  37 {
  38   SCM details = him->get_property ("details");
  39
  40   /*
  41     TODO: put in define-grobs.scm
  42    */
  43   INTER_QUANT_PENALTY = get_detail (details, ly_symbol2scm ("inter-quant-penalty"), 1000.0);
  44   SECONDARY_BEAM_DEMERIT = get_detail (details, ly_symbol2scm ("secondary-beam-demerit"), 10.0);
  45   STEM_LENGTH_DEMERIT_FACTOR = get_detail (details, ly_symbol2scm ("stem-length-demerit-factor"), 5);
  46   REGION_SIZE = get_detail (details, ly_symbol2scm ("region-size"), 2);
  47   BEAM_EPS = get_detail (details, ly_symbol2scm ("beam-eps"), 1e-3);
  48   STEM_LENGTH_LIMIT_PENALTY = get_detail (details, ly_symbol2scm ("stem-length-limit-penalty"), 5000);
  49   DAMPING_DIRECTION_PENALTY = get_detail (details, ly_symbol2scm ("damping-direction-penalty"), 800);
  50   HINT_DIRECTION_PENALTY = get_detail (details, ly_symbol2scm ("hint-direction-penalty"), 20);
  51   MUSICAL_DIRECTION_FACTOR = get_detail (details, ly_symbol2scm ("musical-direction-factor"), 400);
  52   IDEAL_SLOPE_FACTOR = get_detail (details, ly_symbol2scm ("ideal-slope-factor"), 10);
  53   ROUND_TO_ZERO_SLOPE = get_detail (details, ly_symbol2scm ("round-to-zero-slope"), 0.02);
  54 }
  55
  56 static Real
  57 shrink_extra_weight (Real x, Real fac)
  58 {
  59   return fabs (x) * ((x < 0) ? fac : 1.0);
  60 }
  61
  62 struct Quant_score
  63 {
  64   Real yl;
  65   Real yr;
  66   Real demerits;
  67
  68 #if DEBUG_BEAM_SCORING
  69   string score_card_;
  70 #endif
  71 };
  72
  73 /*
  74   TODO:
  75
  76   - Make all demerits customisable
  77
  78   - One sensible check per demerit (what's this --hwn)
  79
  80   - Add demerits for quants per se, as to forbid a specific quant
  81   entirely
  82 */
  83
  84 int
  85 best_quant_score_idx (vector<Quant_score> const &qscores)
  86 {
  87   Real best = 1e6;
  88   int best_idx = -1;
  89   for (vsize i = qscores.size (); i--;)
  90     {
  91       if (qscores[i].demerits < best)
  92         {
  93           best = qscores [i].demerits;
  94           best_idx = i;
  95         }
  96     }
  97
  98   return best_idx;
  99 }
 100
 101 MAKE_SCHEME_CALLBACK (Beam, quanting, 2);
 102 SCM
 103 Beam::quanting (SCM smob, SCM posns)
 104 {
 105   Grob *me = unsmob_grob (smob);
 106
 107   Beam_quant_parameters parameters;
 108   parameters.fill (me);
 109
 110   Real yl = scm_to_double (scm_car (posns));
 111   Real yr = scm_to_double (scm_cdr (posns));
 112
 113   /*
 114     Calculations are relative to a unit-scaled staff, i.e. the quants are
 115     divided by the current staff_space.
 116
 117   */
 118   Real ss = Staff_symbol_referencer::staff_space (me);
 119   Real thickness = Beam::get_thickness (me) / ss;
 120   Real slt = Staff_symbol_referencer::line_thickness (me) / ss;
 121
 122   Real dy_mus = robust_scm2double (me->get_property ("least-squares-dy"), 0);
 123   Real straddle = 0.0;
 124   Real sit = (thickness - slt) / 2;
 125   Real inter = 0.5;
 126   Real hang = 1.0 - (thickness - slt) / 2;
 127   Real quants [] = {straddle, sit, inter, hang };
 128
 129   int num_quants = int (sizeof (quants) / sizeof (Real));
 130   vector<Real> quantsl;
 131   vector<Real> quantsr;
 132
 133   /*
 134     going to REGION_SIZE == 2, yields another 0.6 second with
 135     wtk1-fugue2.
 136
 137
 138     (result indexes between 70 and 575)  ? --hwn.
 139
 140   */
 141
 142   /*
 143     Do stem computations.  These depend on YL and YR linearly, so we can
 144     precompute for every stem 2 factors.
 145   */
 146   vector<Grob*> stems
 147     = extract_grob_array (me, "stems");
 148   vector<Stem_info> stem_infos;
 149   vector<Real> base_lengths;
 150   vector<Real> stem_xposns;
 151
 152   Drul_array<bool> dirs_found (0, 0);
 153   Grob *common[2];
 154   for (int a = 2; a--;)
 155     common[a] = common_refpoint_of_array (stems, me, Axis (a));
 156
 157   Grob *fvs = first_normal_stem (me);
 158   Grob *lvs = last_normal_stem (me);
 159   Real xl = fvs ? fvs->relative_coordinate (common[X_AXIS], X_AXIS) : 0.0;
 160   Real xr = fvs ? lvs->relative_coordinate (common[X_AXIS], X_AXIS) : 0.0;
 161
 162   /*
 163     We store some info to quickly interpolate.  The stemlength are
 164     affine linear in YL and YR. If YL == YR == 0, then we might have
 165     stem_y != 0.0, when we're cross staff.
 166
 167   */
 168   for (vsize i = 0; i < stems.size (); i++)
 169     {
 170       Grob *s = stems[i];
 171
 172       Stem_info si (Stem::get_stem_info (s));
 173       si.scale (1 / ss);
 174       stem_infos.push_back (si);
 175       dirs_found[stem_infos.back ().dir_] = true;
 176
 177       bool f = to_boolean (s->get_property ("french-beaming"))
 178         && s != lvs && s != fvs;
 179
 180       if (Stem::is_normal_stem (s))
 181         {
 182           base_lengths.push_back (calc_stem_y (me, s, common, xl, xr, CENTER,
 183                                                Interval (0, 0), f) / ss);
 184         }
 185       else
 186         {
 187           base_lengths.push_back (0);
 188         }
 189
 190       stem_xposns.push_back (s->relative_coordinate (common[X_AXIS], X_AXIS));
 191     }
 192   bool xstaff = Align_interface::has_interface (common[Y_AXIS]);
 193
 194   Direction ldir = Direction (stem_infos[0].dir_);
 195   Direction rdir = Direction (stem_infos.back ().dir_);
 196   bool is_knee = dirs_found[LEFT] && dirs_found[RIGHT];
 197
 198   int region_size = (int) parameters.REGION_SIZE;
 199
 200   /*
 201     Knees are harder, lets try some more possibilities for knees.
 202   */
 203   if (is_knee)
 204     region_size += 2;
 205
 206   /*
 207     Asymetry ? should run to <= region_size ?
 208   */
 209   for (int i = -region_size; i < region_size; i++)
 210     for (int j = 0; j < num_quants; j++)
 211       {
 212         quantsl.push_back (i + quants[j] + int (yl));
 213         quantsr.push_back (i + quants[j] + int (yr));
 214       }
 215
 216   vector<Quant_score> qscores;
 217
 218   for (vsize l = 0; l < quantsl.size (); l++)
 219     for (vsize r = 0; r < quantsr.size (); r++)
 220       {
 221         Quant_score qs;
 222         qs.yl = quantsl[l];
 223         qs.yr = quantsr[r];
 224         qs.demerits = 0.0;
 225
 226         qscores.push_back (qs);
 227       }
 228
 229   /* This is a longish function, but we don't separate this out into
 230      neat modular separate subfunctions, as the subfunctions would be
 231      called for many values of YL, YR. By precomputing various
 232      parameters outside of the loop, we can save a lot of time. */
 233   for (vsize i = qscores.size (); i--;)
 234     {
 235       Real d = score_slopes_dy (qscores[i].yl, qscores[i].yr,
 236                                 dy_mus, yr- yl,
 237                                 xr - xl,
 238                                 xstaff, &parameters);
 239       qscores[i].demerits += d;
 240
 241 #if DEBUG_BEAM_SCORING
 242       qscores[i].score_card_ += to_string ("S%.2f", d);
 243 #endif
 244     }
 245
 246   Real rad = Staff_symbol_referencer::staff_radius (me);
 247   Drul_array<int> edge_beam_counts
 248     (Stem::beam_multiplicity (stems[0]).length () + 1,
 249      Stem::beam_multiplicity (stems.back ()).length () + 1);
 250
 251   Real beam_translation = get_beam_translation (me) / ss;
 252
 253   Real reasonable_score = (is_knee) ? 200000 : 100;
 254   for (vsize i = qscores.size (); i--;)
 255     if (qscores[i].demerits < reasonable_score)
 256       {
 257         Real d = score_forbidden_quants (qscores[i].yl, qscores[i].yr,
 258                                          rad, slt, thickness, beam_translation,
 259                                          edge_beam_counts, ldir, rdir, &parameters);
 260         qscores[i].demerits += d;
 261
 262 #if DEBUG_BEAM_SCORING
 263         qscores[i].score_card_ += to_string (" F %.2f", d);
 264 #endif
 265       }
 266
 267   for (vsize i = qscores.size (); i--;)
 268     if (qscores[i].demerits < reasonable_score)
 269       {
 270         Real d = score_stem_lengths (stems, stem_infos,
 271                                      base_lengths, stem_xposns,
 272                                      xl, xr,
 273                                      is_knee,
 274                                      qscores[i].yl, qscores[i].yr, &parameters);
 275         qscores[i].demerits += d;
 276
 277 #if DEBUG_BEAM_SCORING
 278         qscores[i].score_card_ += to_string (" L %.2f", d);
 279 #endif
 280       }
 281
 282   int best_idx = best_quant_score_idx (qscores);
 283
 284 #if DEBUG_BEAM_SCORING
 285   SCM inspect_quants = me->get_property ("inspect-quants");
 286   if (to_boolean (me->layout ()->lookup_variable (ly_symbol2scm ("debug-beam-scoring")))
 287       && scm_is_pair (inspect_quants))
 288     {
 289       Drul_array<Real> ins = ly_scm2interval (inspect_quants);
 290
 291       Real mindist = 1e6;
 292       for (vsize i = 0; i < qscores.size (); i++)
 293         {
 294           Real d = fabs (qscores[i].yl- ins[LEFT]) + fabs (qscores[i].yr - ins[RIGHT]);
 295           if (d < mindist)
 296             {
 297               best_idx = i;
 298               mindist = d;
 299             }
 300         }
 301       if (mindist > 1e5)
 302         programming_error ("cannot find quant");
 303     }
 304 #endif
 305
 306   Interval final_positions;
 307   if (best_idx < 0)
 308     {
 309       warning (_ ("no feasible beam position"));
 310       final_positions = Interval (0, 0);
 311     }
 312   else
 313     {
 314       final_positions = Drul_array<Real> (qscores[best_idx].yl,
 315                                           qscores[best_idx].yr);
 316     }
 317
 318 #if DEBUG_BEAM_SCORING
 319   if (best_idx >= 0
 320       && to_boolean (me->layout ()->lookup_variable (ly_symbol2scm ("debug-beam-scoring"))))
 321     {
 322       qscores[best_idx].score_card_ += to_string ("i%d", best_idx);
 323
 324       // debug quanting
 325       me->set_property ("quant-score",
 326                         ly_string2scm (qscores[best_idx].score_card_));
 327     }
 328 #endif
 329
 330   return ly_interval2scm (final_positions);
 331 }
 332
 333 Real
 334 Beam::score_stem_lengths (vector<Grob*> const &stems,
 335                           vector<Stem_info> const &stem_infos,
 336                           vector<Real> const &base_stem_ys,
 337                           vector<Real> const &stem_xs,
 338                           Real xl, Real xr,
 339                           bool knee,
 340                           Real yl, Real yr,
 341
 342                           Beam_quant_parameters const *parameters)
 343 {
 344   Real limit_penalty = parameters->STEM_LENGTH_LIMIT_PENALTY;
 345   Drul_array<Real> score (0, 0);
 346   Drul_array<int> count (0, 0);
 347
 348   for (vsize i = 0; i < stems.size (); i++)
 349     {
 350       Grob *s = stems[i];
 351       if (!Stem::is_normal_stem (s))
 352         continue;
 353
 354       Real x = stem_xs[i];
 355       Real dx = xr - xl;
 356       Real beam_y = dx ? yr * (x - xl) / dx + yl * (xr - x) / dx : (yr + yl) / 2;
 357       Real current_y = beam_y + base_stem_ys[i];
 358       Real length_pen = parameters->STEM_LENGTH_DEMERIT_FACTOR;
 359
 360       Stem_info info = stem_infos[i];
 361       Direction d = info.dir_;
 362
 363       score[d] += limit_penalty * max (0.0, (d * (info.shortest_y_ - current_y)));
 364
 365       Real ideal_diff = d * (current_y - info.ideal_y_);
 366       Real ideal_score = shrink_extra_weight (ideal_diff, 1.5);
 367
 368       /* We introduce a power, to make the scoring strictly
 369          convex. Otherwise a symmetric knee beam (up/down/up/down)
 370          does not have an optimum in the middle. */
 371       if (knee)
 372         ideal_score = pow (ideal_score, 1.1);
 373
 374       score[d] += length_pen * ideal_score;
 375
 376       count[d]++;
 377     }
 378
 379   Direction d = DOWN;
 380   do
 381     score[d] /= max (count[d], 1);
 382   while (flip (&d) != DOWN);
 383
 384   return score[LEFT] + score[RIGHT];
 385 }
 386
 387 Real
 388 Beam::score_slopes_dy (Real yl, Real yr,
 389                        Real dy_mus, Real dy_damp,
 390                        Real dx,
 391                        bool xstaff,
 392
 393                        Beam_quant_parameters const *parameters)
 394 {
 395   Real dy = yr - yl;
 396   Real dem = 0.0;
 397
 398   /*
 399     DAMPING_DIRECTION_PENALTY is a very harsh measure, while for
 400     complex beaming patterns, horizontal is often a good choice.
 401
 402     TODO: find a way to incorporate the complexity of the beam in this
 403     penalty.
 404   */
 405   if (sign (dy_damp) != sign (dy))
 406     {
 407       if (!dy)
 408         {
 409           if (fabs (dy_damp / dx) > parameters->ROUND_TO_ZERO_SLOPE)
 410             dem += parameters->DAMPING_DIRECTION_PENALTY;
 411           else
 412             dem += parameters->HINT_DIRECTION_PENALTY;
 413         }
 414       else
 415         dem += parameters->DAMPING_DIRECTION_PENALTY;
 416     }
 417
 418   dem += parameters->MUSICAL_DIRECTION_FACTOR
 419     * max (0.0, (fabs (dy) - fabs (dy_mus)));
 420
 421   Real slope_penalty = parameters->IDEAL_SLOPE_FACTOR;
 422
 423   /* Xstaff beams tend to use extreme slopes to get short stems. We
 424      put in a penalty here. */
 425   if (xstaff)
 426     slope_penalty *= 10;
 427
 428   /* Huh, why would a too steep beam be better than a too flat one ? */
 429   dem += shrink_extra_weight (fabs (dy_damp) - fabs (dy), 1.5)
 430     * slope_penalty;
 431
 432   return dem;
 433 }
 434
 435 static Real
 436 my_modf (Real x)
 437 {
 438   return x - floor (x);
 439 }
 440
 441 /*
 442   TODO: The fixed value SECONDARY_BEAM_DEMERIT is probably flawed:
 443   because for 32nd and 64th beams the forbidden quants are relatively
 444   more important than stem lengths.
 445 */
 446 Real
 447 Beam::score_forbidden_quants (Real yl, Real yr,
 448                               Real radius,
 449                               Real slt,
 450                               Real thickness, Real beam_translation,
 451                               Drul_array<int> beam_counts,
 452                               Direction ldir, Direction rdir,
 453
 454                               Beam_quant_parameters const *parameters)
 455 {
 456   Real dy = yr - yl;
 457   Drul_array<Real> y (yl, yr);
 458   Drul_array<Direction> dirs (ldir, rdir);
 459
 460   Real extra_demerit = parameters->SECONDARY_BEAM_DEMERIT / (max (beam_counts[LEFT], beam_counts[RIGHT]));
 461
 462   Direction d = LEFT;
 463   Real dem = 0.0;
 464   Real eps = parameters->BEAM_EPS;
 465
 466   do
 467     {
 468       for (int j = 1; j <= beam_counts[d]; j++)
 469         {
 470           Direction stem_dir = dirs[d];
 471
 472           /*
 473             The 2.2 factor is to provide a little leniency for
 474             borderline cases. If we do 2.0, then the upper outer line
 475             will be in the gap of the (2, sit) quant, leading to a
 476             false demerit.
 477           */
 478           Real gap1 = y[d] - stem_dir * ((j - 1) * beam_translation + thickness / 2 - slt / 2.2);
 479           Real gap2 = y[d] - stem_dir * (j * beam_translation - thickness / 2 + slt / 2.2);
 480
 481           Interval gap;
 482           gap.add_point (gap1);
 483           gap.add_point (gap2);
 484
 485           for (Real k = -radius;
 486                k <= radius + eps; k += 1.0)
 487             if (gap.contains (k))
 488               {
 489                 Real dist = min (fabs (gap[UP] - k), fabs (gap[DOWN] - k));
 490
 491                 /*
 492                   this parameter is tuned to grace-stem-length.ly
 493                 */
 494                 Real fixed_demerit = 0.4;
 495
 496                 dem += extra_demerit
 497                   * (fixed_demerit
 498                      + (1 - fixed_demerit) * (dist / gap.length ()) * 2);
 499               }
 500         }
 501     }
 502   while ((flip (&d)) != LEFT);
 503
 504   if (max (beam_counts[LEFT], beam_counts[RIGHT]) >= 2)
 505     {
 506       Real straddle = 0.0;
 507       Real sit = (thickness - slt) / 2;
 508       Real inter = 0.5;
 509       Real hang = 1.0 - (thickness - slt) / 2;
 510
 511       Direction d = LEFT;
 512       do
 513         {
 514           if (beam_counts[d] >= 2
 515               && fabs (y[d] - dirs[d] * beam_translation) < radius + inter)
 516             {
 517               if (dirs[d] == UP && dy <= eps
 518                   && fabs (my_modf (y[d]) - sit) < eps)
 519                 dem += extra_demerit;
 520
 521               if (dirs[d] == DOWN && dy >= eps
 522                   && fabs (my_modf (y[d]) - hang) < eps)
 523                 dem += extra_demerit;
 524             }
 525
 526           if (beam_counts[d] >= 3
 527               && fabs (y[d] - 2 * dirs[d] * beam_translation) < radius + inter)
 528             {
 529               if (dirs[d] == UP && dy <= eps
 530                   && fabs (my_modf (y[d]) - straddle) < eps)
 531                 dem += extra_demerit;
 532
 533               if (dirs[d] == DOWN && dy >= eps
 534                   && fabs (my_modf (y[d]) - straddle) < eps)
 535                 dem += extra_demerit;
 536             }
 537         }
 538       while (flip (&d) != LEFT);
 539     }
 540
 541   return dem;
 542 }
 543