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--2008 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: The default values should be copied to 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   Real ss = Staff_symbol_referencer::staff_space (me);
 118   Real thickness = Beam::get_thickness (me) / ss;
 119   Real slt = Staff_symbol_referencer::line_thickness (me) / ss;
 120
 121   Real dy_mus = robust_scm2double (me->get_property ("least-squares-dy"), 0);
 122   Real straddle = 0.0;
 123   Real sit = (thickness - slt) / 2;
 124   Real inter = 0.5;
 125   Real hang = 1.0 - (thickness - slt) / 2;
 126   Real quants [] = {straddle, sit, inter, hang };
 127
 128   int num_quants = int (sizeof (quants) / sizeof (Real));
 129   vector<Real> quantsl;
 130   vector<Real> quantsr;
 131
 132   /*
 133     going to REGION_SIZE == 2, yields another 0.6 second with
 134     wtk1-fugue2.
 135
 136     (result indexes between 70 and 575)  ? --hwn.
 137
 138   */
 139
 140   /*
 141     Do stem computations.  These depend on YL and YR linearly, so we can
 142     precompute for every stem 2 factors.
 143   */
 144   vector<Grob*> stems
 145     = extract_grob_array (me, "stems");
 146   vector<Stem_info> stem_infos;
 147   vector<Real> base_lengths;
 148   vector<Real> stem_xposns;
 149
 150   Drul_array<bool> dirs_found (0, 0);
 151   Grob *common[2];
 152   for (int a = 2; a--;)
 153     common[a] = common_refpoint_of_array (stems, me, Axis (a));
 154
 155   Grob *fvs = first_normal_stem (me);
 156   Grob *lvs = last_normal_stem (me);
 157   Real xl = fvs ? fvs->relative_coordinate (common[X_AXIS], X_AXIS) : 0.0;
 158   Real xr = fvs ? lvs->relative_coordinate (common[X_AXIS], X_AXIS) : 0.0;
 159
 160   /*
 161     We store some info to quickly interpolate.  The stemlength are
 162     affine linear in YL and YR. If YL == YR == 0, then we might have
 163     stem_y != 0.0, when we're cross staff.
 164
 165   */
 166   for (vsize i = 0; i < stems.size (); i++)
 167     {
 168       Grob *s = stems[i];
 169
 170       Stem_info si (Stem::get_stem_info (s));
 171       si.scale (1 / ss);
 172       stem_infos.push_back (si);
 173       dirs_found[stem_infos.back ().dir_] = true;
 174
 175       bool f = to_boolean (s->get_property ("french-beaming"))
 176         && s != lvs && s != fvs;
 177
 178       if (Stem::is_normal_stem (s))
 179         {
 180           base_lengths.push_back (calc_stem_y (me, s, common, xl, xr, CENTER,
 181                                                Interval (0, 0), f) / ss);
 182         }
 183       else
 184         {
 185           base_lengths.push_back (0);
 186         }
 187
 188       stem_xposns.push_back (s->relative_coordinate (common[X_AXIS], X_AXIS));
 189     }
 190   bool xstaff = Align_interface::has_interface (common[Y_AXIS]);
 191
 192   Direction ldir = Direction (stem_infos[0].dir_);
 193   Direction rdir = Direction (stem_infos.back ().dir_);
 194   bool is_knee = dirs_found[LEFT] && dirs_found[RIGHT];
 195
 196   int region_size = (int) parameters.REGION_SIZE;
 197
 198   /*
 199     Knees are harder, lets try some more possibilities for knees.
 200   */
 201   if (is_knee)
 202     region_size += 2;
 203
 204   /*
 205     Asymetry ? should run to <= region_size ?
 206   */
 207   for (int i = -region_size; i < region_size; i++)
 208     for (int j = 0; j < num_quants; j++)
 209       {
 210         quantsl.push_back (i + quants[j] + int (yl));
 211         quantsr.push_back (i + quants[j] + int (yr));
 212       }
 213
 214   vector<Quant_score> qscores;
 215
 216   for (vsize l = 0; l < quantsl.size (); l++)
 217     for (vsize r = 0; r < quantsr.size (); r++)
 218       {
 219         Quant_score qs;
 220         qs.yl = quantsl[l];
 221         qs.yr = quantsr[r];
 222         qs.demerits = 0.0;
 223
 224         qscores.push_back (qs);
 225       }
 226
 227   /* This is a longish function, but we don't separate this out into
 228      neat modular separate subfunctions, as the subfunctions would be
 229      called for many values of YL, YR. By precomputing various
 230      parameters outside of the loop, we can save a lot of time. */
 231   for (vsize i = qscores.size (); i--;)
 232     {
 233       Real d = score_slopes_dy (qscores[i].yl, qscores[i].yr,
 234                                 dy_mus, yr- yl,
 235                                 xr - xl,
 236                                 xstaff, &parameters);
 237       qscores[i].demerits += d;
 238
 239 #if DEBUG_BEAM_SCORING
 240       qscores[i].score_card_ += to_string ("S%.2f", d);
 241 #endif
 242     }
 243
 244   Real rad = Staff_symbol_referencer::staff_radius (me);
 245   Drul_array<int> edge_beam_counts
 246     (Stem::beam_multiplicity (stems[0]).length () + 1,
 247      Stem::beam_multiplicity (stems.back ()).length () + 1);
 248
 249   Real beam_translation = get_beam_translation (me) / ss;
 250
 251   Real reasonable_score = (is_knee) ? 200000 : 100;
 252   for (vsize i = qscores.size (); i--;)
 253     if (qscores[i].demerits < reasonable_score)
 254       {
 255         Real d = score_forbidden_quants (qscores[i].yl, qscores[i].yr,
 256                                          rad, slt, thickness, beam_translation,
 257                                          edge_beam_counts, ldir, rdir, &parameters);
 258         qscores[i].demerits += d;
 259
 260 #if DEBUG_BEAM_SCORING
 261         qscores[i].score_card_ += to_string (" F %.2f", d);
 262 #endif
 263       }
 264
 265   for (vsize i = qscores.size (); i--;)
 266     if (qscores[i].demerits < reasonable_score)
 267       {
 268         Real d = score_stem_lengths (stems, stem_infos,
 269                                      base_lengths, stem_xposns,
 270                                      xl, xr,
 271                                      is_knee,
 272                                      qscores[i].yl, qscores[i].yr, &parameters);
 273         qscores[i].demerits += d;
 274
 275 #if DEBUG_BEAM_SCORING
 276         qscores[i].score_card_ += to_string (" L %.2f", d);
 277 #endif
 278       }
 279
 280   int best_idx = best_quant_score_idx (qscores);
 281
 282 #if DEBUG_BEAM_SCORING
 283   SCM inspect_quants = me->get_property ("inspect-quants");
 284   if (to_boolean (me->layout ()->lookup_variable (ly_symbol2scm ("debug-beam-scoring")))
 285       && scm_is_pair (inspect_quants))
 286     {
 287       Drul_array<Real> ins = ly_scm2interval (inspect_quants);
 288
 289       Real mindist = 1e6;
 290       for (vsize i = 0; i < qscores.size (); i++)
 291         {
 292           Real d = fabs (qscores[i].yl- ins[LEFT]) + fabs (qscores[i].yr - ins[RIGHT]);
 293           if (d < mindist)
 294             {
 295               best_idx = i;
 296               mindist = d;
 297             }
 298         }
 299       if (mindist > 1e5)
 300         programming_error ("cannot find quant");
 301     }
 302 #endif
 303
 304   Interval final_positions;
 305   if (best_idx < 0)
 306     {
 307       warning (_ ("no feasible beam position"));
 308       final_positions = Interval (0, 0);
 309     }
 310   else
 311     {
 312       final_positions = Drul_array<Real> (qscores[best_idx].yl,
 313                                           qscores[best_idx].yr);
 314     }
 315
 316 #if DEBUG_BEAM_SCORING
 317   if (best_idx >= 0
 318       && to_boolean (me->layout ()->lookup_variable (ly_symbol2scm ("debug-beam-scoring"))))
 319     {
 320       qscores[best_idx].score_card_ += to_string ("i%d", best_idx);
 321
 322       // debug quanting
 323       me->set_property ("quant-score",
 324                         ly_string2scm (qscores[best_idx].score_card_));
 325     }
 326 #endif
 327
 328   return ly_interval2scm (final_positions);
 329 }
 330
 331 Real
 332 Beam::score_stem_lengths (vector<Grob*> const &stems,
 333                           vector<Stem_info> const &stem_infos,
 334                           vector<Real> const &base_stem_ys,
 335                           vector<Real> const &stem_xs,
 336                           Real xl, Real xr,
 337                           bool knee,
 338                           Real yl, Real yr,
 339
 340                           Beam_quant_parameters const *parameters)
 341 {
 342   Real limit_penalty = parameters->STEM_LENGTH_LIMIT_PENALTY;
 343   Drul_array<Real> score (0, 0);
 344   Drul_array<int> count (0, 0);
 345
 346   for (vsize i = 0; i < stems.size (); i++)
 347     {
 348       Grob *s = stems[i];
 349       if (!Stem::is_normal_stem (s))
 350         continue;
 351
 352       Real x = stem_xs[i];
 353       Real dx = xr - xl;
 354       Real beam_y = dx ? yr * (x - xl) / dx + yl * (xr - x) / dx : (yr + yl) / 2;
 355       Real current_y = beam_y + base_stem_ys[i];
 356       Real length_pen = parameters->STEM_LENGTH_DEMERIT_FACTOR;
 357
 358       Stem_info info = stem_infos[i];
 359       Direction d = info.dir_;
 360
 361       score[d] += limit_penalty * max (0.0, (d * (info.shortest_y_ - current_y)));
 362
 363       Real ideal_diff = d * (current_y - info.ideal_y_);
 364       Real ideal_score = shrink_extra_weight (ideal_diff, 1.5);
 365
 366       /* We introduce a power, to make the scoring strictly
 367          convex. Otherwise a symmetric knee beam (up/down/up/down)
 368          does not have an optimum in the middle. */
 369       if (knee)
 370         ideal_score = pow (ideal_score, 1.1);
 371
 372       score[d] += length_pen * ideal_score;
 373
 374       count[d]++;
 375     }
 376
 377   Direction d = DOWN;
 378   do
 379     score[d] /= max (count[d], 1);
 380   while (flip (&d) != DOWN);
 381
 382   return score[LEFT] + score[RIGHT];
 383 }
 384
 385 Real
 386 Beam::score_slopes_dy (Real yl, Real yr,
 387                        Real dy_mus, Real dy_damp,
 388                        Real dx,
 389                        bool xstaff,
 390
 391                        Beam_quant_parameters const *parameters)
 392 {
 393   Real dy = yr - yl;
 394   Real dem = 0.0;
 395
 396   /*
 397     DAMPING_DIRECTION_PENALTY is a very harsh measure, while for
 398     complex beaming patterns, horizontal is often a good choice.
 399
 400     TODO: find a way to incorporate the complexity of the beam in this
 401     penalty.
 402   */
 403   if (sign (dy_damp) != sign (dy))
 404     {
 405       if (!dy)
 406         {
 407           if (fabs (dy_damp / dx) > parameters->ROUND_TO_ZERO_SLOPE)
 408             dem += parameters->DAMPING_DIRECTION_PENALTY;
 409           else
 410             dem += parameters->HINT_DIRECTION_PENALTY;
 411         }
 412       else
 413         dem += parameters->DAMPING_DIRECTION_PENALTY;
 414     }
 415
 416   dem += parameters->MUSICAL_DIRECTION_FACTOR
 417     * max (0.0, (fabs (dy) - fabs (dy_mus)));
 418
 419   Real slope_penalty = parameters->IDEAL_SLOPE_FACTOR;
 420
 421   /* Xstaff beams tend to use extreme slopes to get short stems. We
 422      put in a penalty here. */
 423   if (xstaff)
 424     slope_penalty *= 10;
 425
 426   /* Huh, why would a too steep beam be better than a too flat one ? */
 427   dem += shrink_extra_weight (fabs (dy_damp) - fabs (dy), 1.5)
 428     * slope_penalty;
 429
 430   return dem;
 431 }
 432
 433 static Real
 434 my_modf (Real x)
 435 {
 436   return x - floor (x);
 437 }
 438
 439 /*
 440   TODO: The fixed value SECONDARY_BEAM_DEMERIT is probably flawed:
 441   because for 32nd and 64th beams the forbidden quants are relatively
 442   more important than stem lengths.
 443 */
 444 Real
 445 Beam::score_forbidden_quants (Real yl, Real yr,
 446                               Real radius,
 447                               Real slt,
 448                               Real thickness, Real beam_translation,
 449                               Drul_array<int> beam_counts,
 450                               Direction ldir, Direction rdir,
 451
 452                               Beam_quant_parameters const *parameters)
 453 {
 454   Real dy = yr - yl;
 455   Drul_array<Real> y (yl, yr);
 456   Drul_array<Direction> dirs (ldir, rdir);
 457
 458   Real extra_demerit = parameters->SECONDARY_BEAM_DEMERIT / (max (beam_counts[LEFT], beam_counts[RIGHT]));
 459
 460   Direction d = LEFT;
 461   Real dem = 0.0;
 462   Real eps = parameters->BEAM_EPS;
 463
 464   do
 465     {
 466       for (int j = 1; j <= beam_counts[d]; j++)
 467         {
 468           Direction stem_dir = dirs[d];
 469
 470           /*
 471             The 2.2 factor is to provide a little leniency for
 472             borderline cases. If we do 2.0, then the upper outer line
 473             will be in the gap of the (2, sit) quant, leading to a
 474             false demerit.
 475           */
 476           Real gap1 = y[d] - stem_dir * ((j - 1) * beam_translation + thickness / 2 - slt / 2.2);
 477           Real gap2 = y[d] - stem_dir * (j * beam_translation - thickness / 2 + slt / 2.2);
 478
 479           Interval gap;
 480           gap.add_point (gap1);
 481           gap.add_point (gap2);
 482
 483           for (Real k = -radius;
 484                k <= radius + eps; k += 1.0)
 485             if (gap.contains (k))
 486               {
 487                 Real dist = min (fabs (gap[UP] - k), fabs (gap[DOWN] - k));
 488
 489                 /*
 490                   this parameter is tuned to grace-stem-length.ly
 491                 */
 492                 Real fixed_demerit = 0.4;
 493
 494                 dem += extra_demerit
 495                   * (fixed_demerit
 496                      + (1 - fixed_demerit) * (dist / gap.length ()) * 2);
 497               }
 498         }
 499     }
 500   while ((flip (&d)) != LEFT);
 501
 502   if (max (beam_counts[LEFT], beam_counts[RIGHT]) >= 2)
 503     {
 504       Real straddle = 0.0;
 505       Real sit = (thickness - slt) / 2;
 506       Real inter = 0.5;
 507       Real hang = 1.0 - (thickness - slt) / 2;
 508
 509       Direction d = LEFT;
 510       do
 511         {
 512           if (beam_counts[d] >= 2
 513               && fabs (y[d] - dirs[d] * beam_translation) < radius + inter)
 514             {
 515               if (dirs[d] == UP && dy <= eps
 516                   && fabs (my_modf (y[d]) - sit) < eps)
 517                 dem += extra_demerit;
 518
 519               if (dirs[d] == DOWN && dy >= eps
 520                   && fabs (my_modf (y[d]) - hang) < eps)
 521                 dem += extra_demerit;
 522             }
 523
 524           if (beam_counts[d] >= 3
 525               && fabs (y[d] - 2 * dirs[d] * beam_translation) < radius + inter)
 526             {
 527               if (dirs[d] == UP && dy <= eps
 528                   && fabs (my_modf (y[d]) - straddle) < eps)
 529                 dem += extra_demerit;
 530
 531               if (dirs[d] == DOWN && dy >= eps
 532                   && fabs (my_modf (y[d]) - straddle) < eps)
 533                 dem += extra_demerit;
 534             }
 535         }
 536       while (flip (&d) != LEFT);
 537     }
 538
 539   return dem;
 540 }
 541