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