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