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