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--2006 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 "align-interface.hh"
  16 #include "international.hh"
  17 #include "output-def.hh"
  18 #include "pointer-group-interface.hh"
  19 #include "staff-symbol-referencer.hh"
  20 #include "stem.hh"
  21 #include "warn.hh"
  22
  23 Real
  24 get_detail (SCM alist, SCM sym, Real def)
  25 {
  26   SCM entry = scm_assq (sym, alist);
  27
  28   if (scm_is_pair (entry))
  29     return robust_scm2double (scm_cdr (entry), def);
  30   return def;
  31 }
  32
  33 void
  34 Beam_quant_parameters::fill (Grob *him)
  35 {
  36   SCM details = him->get_property ("details");
  37
  38   INTER_QUANT_PENALTY = get_detail (details, ly_symbol2scm ("inter-quant-penalty"), 1000.0);
  39   SECONDARY_BEAM_DEMERIT = get_detail (details, ly_symbol2scm ("secondary-beam-demerit"), 10.0);
  40   STEM_LENGTH_DEMERIT_FACTOR = get_detail (details, ly_symbol2scm ("stem-length-demerit-factor"), 5);
  41   REGION_SIZE = get_detail (details, ly_symbol2scm ("region-size"), 2);
  42   BEAM_EPS = get_detail (details, ly_symbol2scm ("beam-eps"), 1e-3);
  43
  44   // possibly ridiculous, but too short stems just won't do
  45   STEM_LENGTH_LIMIT_PENALTY = get_detail (details, ly_symbol2scm ("stem-length-limit-penalty"), 5000);
  46   DAMPING_DIRECTION_PENALTY = get_detail (details, ly_symbol2scm ("damping-direction-penalty"), 800);
  47   MUSICAL_DIRECTION_FACTOR = get_detail (details, ly_symbol2scm ("musical-direction-factor"), 400);
  48   IDEAL_SLOPE_FACTOR = get_detail (details, ly_symbol2scm ("ideal-slope-factor"), 10);
  49   ROUND_TO_ZERO_SLOPE = get_detail (details, ly_symbol2scm ("round-to-zero-slope"), 0.02);
  50 }
  51
  52 static Real
  53 shrink_extra_weight (Real x, Real fac)
  54 {
  55   return fabs (x) * ((x < 0) ? fac : 1.0);
  56 }
  57
  58 struct Quant_score
  59 {
  60   Real yl;
  61   Real yr;
  62   Real demerits;
  63
  64 #if DEBUG_QUANTING
  65   std::string score_card_;
  66 #endif
  67 };
  68
  69 /*
  70   TODO:
  71
  72   - Make all demerits customisable
  73
  74   - One sensible check per demerit (what's this --hwn)
  75
  76   - Add demerits for quants per se, as to forbid a specific quant
  77   entirely
  78 */
  79
  80 int
  81 best_quant_score_idx (std::vector<Quant_score> const &qscores)
  82 {
  83   Real best = 1e6;
  84   int best_idx = -1;
  85   for (vsize i = qscores.size (); i--;)
  86     {
  87       if (qscores[i].demerits < best)
  88         {
  89           best = qscores [i].demerits;
  90           best_idx = i;
  91         }
  92     }
  93
  94   return best_idx;
  95 }
  96
  97 MAKE_SCHEME_CALLBACK (Beam, quanting, 2);
  98 SCM
  99 Beam::quanting (SCM smob, SCM posns)
 100 {
 101   Grob *me = unsmob_grob (smob);
 102
 103   Beam_quant_parameters parameters;
 104   parameters.fill (me);
 105
 106   Real yl = scm_to_double (scm_car (posns));
 107   Real yr = scm_to_double (scm_cdr (posns));
 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   std::vector<Real> quantsl;
 127   std::vector<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   std::vector<Stem_info> stem_infos;
 145   std::vector<Real> base_lengths;
 146   std::vector<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 (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       base_lengths.push_back (calc_stem_y (me, s, common, xl, xr,
 179                                       Interval (0, 0), f) / ss);
 180       stem_xposns.push_back (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.back ().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_back (i + quants[j] + int (yl));
 209         quantsr.push_back (i + quants[j] + int (yr));
 210       }
 211
 212   std::vector<Quant_score> qscores;
 213
 214   for (vsize l = 0; l < quantsl.size (); l++)
 215     for (vsize 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_back (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 (vsize 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.back ()).length () + 1);
 246
 247   Real beam_translation = get_beam_translation (me) / ss;
 248
 249   Real reasonable_score = (is_knee) ? 200000 : 100;
 250   for (vsize 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 (vsize 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->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       Real mindist = 1e6;
 288       for (vsize i = 0; i < qscores.size (); i++)
 289         {
 290           Real d = fabs (qscores[i].yl- ins[LEFT]) + fabs (qscores[i].yr - ins[RIGHT]);
 291           if (d < mindist)
 292             {
 293               best_idx = i;
 294               mindist = d;
 295             }
 296         }
 297       if (mindist > 1e5)
 298         programming_error ("can't find quant");
 299     }
 300 #endif
 301
 302   Interval final_positions;
 303   if (best_idx < 0)
 304     {
 305       warning (_ ("no feasible beam position"));
 306       final_positions = Interval (0, 0);
 307     }
 308   else
 309     {
 310       final_positions = Drul_array<Real> (qscores[best_idx].yl,
 311                                           qscores[best_idx].yr);
 312     }
 313
 314 #if DEBUG_QUANTING
 315   if (best_idx >= 0
 316       && to_boolean (me->layout ()->lookup_variable (ly_symbol2scm ("debug-beam-quanting"))))
 317     {
 318       qscores[best_idx].score_card_ += to_string ("i%d", best_idx);
 319
 320       // debug quanting
 321       me->set_property ("quant-score",
 322                         scm_makfrom0str (qscores[best_idx].score_card_.c_str ()));
 323     }
 324 #endif
 325
 326   return ly_interval2scm (final_positions);
 327 }
 328
 329 Real
 330 Beam::score_stem_lengths (Link_array__Grob_ const &stems,
 331                           std::vector<Stem_info> const &stem_infos,
 332                           std::vector<Real> const &base_stem_ys,
 333                           std::vector<Real> const &stem_xs,
 334                           Real xl, Real xr,
 335                           bool knee,
 336                           Real yl, Real yr,
 337
 338                           Beam_quant_parameters const *parameters)
 339 {
 340   Real limit_penalty = parameters->STEM_LENGTH_LIMIT_PENALTY;
 341   Drul_array<Real> score (0, 0);
 342   Drul_array<int> count (0, 0);
 343
 344   for (vsize i = 0; i < stems.size (); i++)
 345     {
 346       Grob *s = stems[i];
 347       if (Stem::is_invisible (s))
 348         continue;
 349
 350       Real x = stem_xs[i];
 351       Real dx = xr - xl;
 352       Real beam_y = dx ? yr * (x - xl) / dx + yl * (xr - x) / dx : (yr + yl) / 2;
 353       Real current_y = beam_y + base_stem_ys[i];
 354       Real length_pen = parameters->STEM_LENGTH_DEMERIT_FACTOR;
 355
 356       Stem_info info = stem_infos[i];
 357       Direction d = info.dir_;
 358
 359       score[d] += limit_penalty * max (0.0, (d * (info.shortest_y_ - current_y)));
 360
 361       Real ideal_diff = d * (current_y - info.ideal_y_);
 362       Real ideal_score = shrink_extra_weight (ideal_diff, 1.5);
 363
 364       /* We introduce a power, to make the scoring strictly
 365          convex. Otherwise a symmetric knee beam (up/down/up/down)
 366          does not have an optimum in the middle. */
 367       if (knee)
 368         ideal_score = pow (ideal_score, 1.1);
 369
 370       score[d] += length_pen * ideal_score;
 371
 372       count[d]++;
 373     }
 374
 375   Direction d = DOWN;
 376   do
 377     score[d] /= max (count[d], 1);
 378   while (flip (&d) != DOWN)
 379     ;
 380
 381   return score[LEFT] + score[RIGHT];
 382 }
 383
 384 Real
 385 Beam::score_slopes_dy (Real yl, Real yr,
 386                        Real dy_mus, Real dy_damp,
 387                        Real dx,
 388                        bool xstaff,
 389
 390                        Beam_quant_parameters const *parameters)
 391 {
 392   Real dy = yr - yl;
 393   Real dem = 0.0;
 394
 395   /*
 396     DAMPING_DIRECTION_PENALTY is a very harsh measure, while for
 397     complex beaming patterns, horizontal is often a good choice.
 398
 399     TODO: find a way to incorporate the complexity of the beam in this
 400     penalty.
 401   */
 402   if (fabs (dy / dx) > parameters->ROUND_TO_ZERO_SLOPE
 403       && sign (dy_damp) != sign (dy))
 404     dem += parameters->DAMPING_DIRECTION_PENALTY;
 405
 406   dem += parameters->MUSICAL_DIRECTION_FACTOR
 407     * max (0.0, (fabs (dy) - fabs (dy_mus)));
 408
 409   Real slope_penalty = parameters->IDEAL_SLOPE_FACTOR;
 410
 411   /* Xstaff beams tend to use extreme slopes to get short stems. We
 412      put in a penalty here. */
 413   if (xstaff)
 414     slope_penalty *= 10;
 415
 416   /* Huh, why would a too steep beam be better than a too flat one ? */
 417   dem += shrink_extra_weight (fabs (dy_damp) - fabs (dy), 1.5)
 418     * slope_penalty;
 419
 420   return dem;
 421 }
 422
 423 static Real
 424 my_modf (Real x)
 425 {
 426   return x - floor (x);
 427 }
 428
 429 /*
 430   TODO: The fixed value SECONDARY_BEAM_DEMERIT is probably flawed:
 431   because for 32nd and 64th beams the forbidden quants are relatively
 432   more important than stem lengths.
 433 */
 434 Real
 435 Beam::score_forbidden_quants (Real yl, Real yr,
 436                               Real radius,
 437                               Real slt,
 438                               Real thickness, Real beam_translation,
 439                               Drul_array<int> beam_counts,
 440                               Direction ldir, Direction rdir,
 441
 442                               Beam_quant_parameters const *parameters)
 443 {
 444   Real dy = yr - yl;
 445   Drul_array<Real> y (yl, yr);
 446   Drul_array<Direction> dirs (ldir, rdir);
 447
 448   Real extra_demerit = parameters->SECONDARY_BEAM_DEMERIT / (max (beam_counts[LEFT], beam_counts[RIGHT]));
 449
 450   Direction d = LEFT;
 451   Real dem = 0.0;
 452   Real eps = parameters->BEAM_EPS;
 453
 454   do
 455     {
 456       for (int j = 1; j <= beam_counts[d]; j++)
 457         {
 458           Direction stem_dir = dirs[d];
 459
 460           /*
 461             The 2.2 factor is to provide a little leniency for
 462             borderline cases. If we do 2.0, then the upper outer line
 463             will be in the gap of the (2, sit) quant, leading to a
 464             false demerit.
 465           */
 466           Real gap1 = y[d] - stem_dir * ((j - 1) * beam_translation + thickness / 2 - slt / 2.2);
 467           Real gap2 = y[d] - stem_dir * (j * beam_translation - thickness / 2 + slt / 2.2);
 468
 469           Interval gap;
 470           gap.add_point (gap1);
 471           gap.add_point (gap2);
 472
 473           for (Real k = -radius;
 474                k <= radius + eps; k += 1.0)
 475             if (gap.contains (k))
 476               {
 477                 Real dist = min (fabs (gap[UP] - k), fabs (gap[DOWN] - k));
 478
 479                 /*
 480                   this parameter is tuned to grace-stem-length.ly
 481                 */
 482                 Real fixed_demerit = 0.4;
 483
 484                 dem += extra_demerit
 485                   * (fixed_demerit
 486                      + (1 - fixed_demerit) * (dist / gap.length ()) * 2);
 487               }
 488         }
 489     }
 490   while ((flip (&d)) != LEFT);
 491
 492   if (max (beam_counts[LEFT], beam_counts[RIGHT]) >= 2)
 493     {
 494       Real straddle = 0.0;
 495       Real sit = (thickness - slt) / 2;
 496       Real inter = 0.5;
 497       Real hang = 1.0 - (thickness - slt) / 2;
 498
 499       Direction d = LEFT;
 500       do
 501         {
 502           if (beam_counts[d] >= 2
 503               && fabs (y[d] - dirs[d] * beam_translation) < radius + inter)
 504             {
 505               if (dirs[d] == UP && dy <= eps
 506                   && fabs (my_modf (y[d]) - sit) < eps)
 507                 dem += extra_demerit;
 508
 509               if (dirs[d] == DOWN && dy >= eps
 510                   && fabs (my_modf (y[d]) - hang) < eps)
 511                 dem += extra_demerit;
 512             }
 513
 514           if (beam_counts[d] >= 3
 515               && fabs (y[d] - 2 * dirs[d] * beam_translation) < radius + inter)
 516             {
 517               if (dirs[d] == UP && dy <= eps
 518                   && fabs (my_modf (y[d]) - straddle) < eps)
 519                 dem += extra_demerit;
 520
 521               if (dirs[d] == DOWN && dy >= eps
 522                   && fabs (my_modf (y[d]) - straddle) < eps)
 523                 dem += extra_demerit;
 524             }
 525         }
 526       while (flip (&d) != LEFT);
 527     }
 528
 529   return dem;
 530 }
 531