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