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