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
  74 best_quant_score_idx (Array<Quant_score>  const & qscores)
  75 {
  76   Real best = 1e6;
  77   int best_idx = -1;
  78   for (int i = qscores.size (); i--;)
  79     {
  80       if (qscores[i].demerits < best)
  81         {
  82           best = qscores [i].demerits ;
  83           best_idx = i;
  84         }
  85     }
  86
  87   if (best_idx < 0)
  88     {
  89       programming_error ("Huh? No best beam quant score?");
  90       best_idx = 0;
  91     }
  92
  93
  94   return best_idx;
  95 }
  96
  97 MAKE_SCHEME_CALLBACK (Beam, quanting, 1);
  98 SCM
  99 Beam::quanting (SCM smob)
 100 {
 101   Grob *me = unsmob_grob (smob);
 102
 103   SCM s = me->get_property ("positions");
 104   Real yl = scm_to_double (scm_car (s));
 105   Real yr = scm_to_double (scm_cdr (s));
 106
 107
 108   /*
 109     Calculations are relative to a unit-scaled staff, i.e. the quants are
 110     divided by the current staff_space.
 111
 112    */
 113   Real ss = Staff_symbol_referencer::staff_space (me);
 114   Real thickness = Beam::get_thickness (me) / ss ;
 115   Real slt = Staff_symbol_referencer::line_thickness (me) / ss;
 116
 117   Real dy_mus = robust_scm2double (me->get_property ("least-squares-dy"), 0);
 118   Real straddle = 0.0;
 119   Real sit = (thickness - slt) / 2;
 120   Real inter = 0.5;
 121   Real hang = 1.0 - (thickness - slt) / 2;
 122   Real quants [] = {straddle, sit, inter, hang };
 123
 124
 125
 126   int num_quants = int (sizeof (quants)/sizeof (Real));
 127   Array<Real> quantsl;
 128   Array<Real> quantsr;
 129
 130   /*
 131     going to REGION_SIZE == 2, yields another 0.6 second with
 132     wtk1-fugue2.
 133
 134
 135     (result indexes between 70 and 575)  ? --hwn.
 136
 137   */
 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     Pointer_group_interface__extract_grobs (me, (Grob*)0, "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
 198   int region_size = REGION_SIZE;
 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);
 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
 247
 248
 249   Drul_array<int> edge_beam_counts
 250     (Stem::beam_multiplicity (stems[0]).length  () + 1,
 251      Stem::beam_multiplicity (stems.top ()).length  () + 1);
 252
 253   Real beam_translation = get_beam_translation (me) / ss;
 254
 255   Real reasonable_score = (is_knee) ? 200000 : 100;
 256   for (int i = qscores.size (); i--;)
 257     if (qscores[i].demerits < reasonable_score)
 258       {
 259         Real d = score_forbidden_quants (qscores[i].yl, qscores[i].yr,
 260                                      rad, slt, thickness, beam_translation,
 261                                      edge_beam_counts, ldir, rdir);
 262         qscores[i].demerits += d;
 263
 264 #if DEBUG_QUANTING
 265         qscores[i].score_card_ += to_string (" F %.2f", d);
 266 #endif
 267       }
 268
 269   for (int i = qscores.size (); i--;)
 270     if (qscores[i].demerits < reasonable_score)
 271       {
 272         Real d=score_stem_lengths (stems, stem_infos,
 273                                  base_lengths, stem_xposns,
 274                                  xl, xr,
 275                                  is_knee,
 276                                  qscores[i].yl, qscores[i].yr);
 277         qscores[i].demerits +=  d;
 278
 279 #if DEBUG_QUANTING
 280         qscores[i].score_card_ += to_string (" L %.2f", d);
 281 #endif
 282       }
 283
 284   int best_idx = best_quant_score_idx (qscores);
 285
 286 #if DEBUG_QUANTING
 287   SCM inspect_quants = me->get_property ("inspect-quants");
 288   if (to_boolean (me->get_paper ()->lookup_variable (ly_symbol2scm ("debug-beam-quanting")))
 289       && scm_is_pair (inspect_quants))
 290     {
 291       Drul_array<Real> ins = ly_scm2interval (inspect_quants);
 292
 293       int i = 0;
 294
 295       Real mindist = 1e6;
 296       for (; i < qscores.size (); i ++)
 297         {
 298           Real d =fabs (qscores[i].yl- ins[LEFT]) + fabs (qscores[i].yr - ins[RIGHT]);
 299           if (d < mindist)
 300             {
 301               best_idx = i;
 302               mindist= d;
 303             }
 304         }
 305       if (mindist > 1e5)
 306         programming_error ("Could not find quant.");
 307     }
 308 #endif
 309
 310   me->set_property ("positions",
 311                     ly_interval2scm (Drul_array<Real> (qscores[best_idx].yl,
 312                                                        qscores[best_idx].yr)));
 313 #if DEBUG_QUANTING
 314   if (to_boolean (me->get_paper ()->lookup_variable (ly_symbol2scm ("debug-beam-quanting"))))
 315     {
 316       qscores[best_idx].score_card_ += to_string ("i%d", best_idx);
 317
 318       // debug quanting
 319       me->set_property ("quant-score",
 320                         scm_makfrom0str (qscores[best_idx].score_card_.to_str0 ()));
 321     }
 322 #endif
 323
 324   return SCM_UNSPECIFIED;
 325 }
 326
 327 Real
 328 Beam::score_stem_lengths (Link_array<Grob> const &stems,
 329                           Array<Stem_info> const &stem_infos,
 330                           Array<Real> const &base_stem_ys,
 331                           Array<Real> const &stem_xs,
 332                           Real xl, Real xr,
 333                           bool knee,
 334                           Real yl, Real yr)
 335 {
 336   Real limit_penalty = STEM_LENGTH_LIMIT_PENALTY;
 337   Drul_array<Real> score (0, 0);
 338   Drul_array<int> count (0, 0);
 339
 340   for (int i=0; i < stems.size (); i++)
 341     {
 342       Grob* s = stems[i];
 343       if (Stem::is_invisible (s))
 344         continue;
 345
 346       Real x = stem_xs[i];
 347       Real dx = xr-xl;
 348       Real beam_y = dx ? yr *(x - xl)/dx + yl * ( xr - x)/dx : (yr + yl)/2;
 349       Real current_y = beam_y + base_stem_ys[i];
 350       Real length_pen = STEM_LENGTH_DEMERIT_FACTOR;
 351
 352       Stem_info info = stem_infos[i];
 353       Direction d = info.dir_;
 354
 355       score[d] += limit_penalty * (0 >? (d * (info.shortest_y_ - current_y)));
 356
 357       Real ideal_diff = d * (current_y - info.ideal_y_);
 358       Real ideal_score = shrink_extra_weight (ideal_diff, 1.5);
 359
 360       /* We introduce a power, to make the scoring strictly
 361          convex. Otherwise a symmetric knee beam (up/down/up/down)
 362          does not have an optimum in the middle. */
 363       if (knee)
 364         ideal_score = pow (ideal_score, 1.1);
 365
 366       score[d] += length_pen * ideal_score;
 367
 368       count[d] ++;
 369     }
 370
 371   Direction d = DOWN;
 372   do
 373     {
 374       score[d] /= (count[d] >? 1);
 375     }
 376   while (flip (&d) != DOWN);
 377
 378   return score[LEFT]+score[RIGHT];
 379 }
 380
 381 Real
 382 Beam::score_slopes_dy (Real yl, Real yr,
 383                        Real dy_mus, Real dy_damp,
 384                        Real dx,
 385                        bool xstaff)
 386 {
 387   Real dy = yr - yl;
 388   Real dem = 0.0;
 389
 390   /*
 391     DAMPING_DIRECTION_PENALTY is a very harsh measure, while for
 392     complex beaming patterns, horizontal is often a good choice.
 393
 394     TODO: find a way to incorporate the complexity of the beam in this
 395     penalty.
 396    */
 397   if (fabs (dy/dx) > ROUND_TO_ZERO_SLOPE
 398       && sign (dy_damp) != sign (dy))
 399     {
 400       dem += DAMPING_DIRECTION_PENALTY;
 401     }
 402
 403    dem += MUSICAL_DIRECTION_FACTOR * (0 >? (fabs (dy) - fabs (dy_mus)));
 404
 405
 406    Real slope_penalty = IDEAL_SLOPE_FACTOR;
 407
 408    /* Xstaff beams tend to use extreme slopes to get short stems. We
 409       put in a penalty here. */
 410    if (xstaff)
 411      slope_penalty *= 10;
 412
 413    /* Huh, why would a too steep beam be better than a too flat one ? */
 414    dem += shrink_extra_weight (fabs (dy_damp) - fabs (dy), 1.5)
 415      * slope_penalty;
 416
 417    return dem;
 418 }
 419
 420
 421 static Real
 422 my_modf (Real x)
 423 {
 424   return x - floor (x);
 425 }
 426
 427
 428 /*
 429   TODO: The fixed value SECONDARY_BEAM_DEMERIT is probably flawed:
 430    because for 32nd and 64th beams the forbidden quants are relatively
 431    more important than stem lengths.
 432 */
 433 Real
 434 Beam::score_forbidden_quants (Real yl, Real yr,
 435                               Real radius,
 436                               Real slt,
 437                               Real thickness, Real beam_translation,
 438                               Drul_array<int> beam_counts,
 439                               Direction ldir, Direction rdir)
 440 {
 441   Real dy = yr - yl;
 442   Drul_array<Real> y (yl,yr);
 443   Drul_array<Direction> dirs (ldir,rdir);
 444
 445   Real extra_demerit = SECONDARY_BEAM_DEMERIT / (beam_counts[LEFT] >? beam_counts[RIGHT]);
 446
 447   Direction d = LEFT;
 448   Real dem = 0.0;
 449
 450
 451   do
 452     {
 453       for (int j = 1; j <= beam_counts[d]; j++)
 454         {
 455           Direction stem_dir = dirs[d];
 456
 457           /*
 458             The 2.2 factor is to provide a little leniency for
 459             borderline cases. If we do 2.0, then the upper outer line
 460             will be in the gap of the (2,sit) quant, leading to a
 461             false demerit.
 462           */
 463           Real gap1 = y[d] - stem_dir * ((j-1) * beam_translation + thickness / 2 - slt/2.2 );
 464           Real gap2 = y[d] - stem_dir * (j * beam_translation - thickness / 2 + slt/2.2);
 465
 466           Interval gap;
 467           gap.add_point (gap1);
 468           gap.add_point (gap2);
 469
 470           for (Real k = - radius ;
 471                k <= radius + BEAM_EPS; k += 1.0)
 472             if (gap.contains (k))
 473               {
 474                 Real dist = fabs (gap[UP]-k) <? fabs (gap[DOWN] - k);
 475
 476                 /*
 477                   this parameter is tuned to grace-stem-length.ly
 478                 */
 479                 Real fixed_demerit = 0.4;
 480
 481                 dem += extra_demerit
 482                   * (fixed_demerit +
 483                      (1-fixed_demerit) * (dist / gap.length())* 2);
 484               }
 485         }
 486     }
 487   while ((flip (&d))!= LEFT);
 488
 489
 490   if ((beam_counts[LEFT] >? beam_counts[RIGHT]) >= 2)
 491     {
 492       Real straddle = 0.0;
 493       Real sit = (thickness - slt) / 2;
 494       Real inter = 0.5;
 495       Real hang = 1.0 - (thickness - slt) / 2;
 496
 497
 498       Direction d = LEFT;
 499       do
 500         {
 501           if (beam_counts[d] >= 2
 502               && fabs (y[d] - dirs[d] * beam_translation) < radius + inter)
 503             {
 504               if (dirs[d] == UP && dy <= BEAM_EPS
 505                   && fabs (my_modf (y[d]) - sit) < BEAM_EPS)
 506                 dem += extra_demerit;
 507
 508               if (dirs[d] == DOWN && dy >= BEAM_EPS
 509                   && fabs (my_modf (y[d]) - hang) < BEAM_EPS)
 510                 dem += extra_demerit;
 511             }
 512
 513           if (beam_counts[d] >= 3
 514               && fabs (y[d] - 2 * dirs[d] * beam_translation) < radius + inter)
 515             {
 516               if (dirs[d] == UP && dy <= BEAM_EPS
 517                   && fabs (my_modf (y[d]) - straddle) < BEAM_EPS)
 518                 dem += extra_demerit;
 519
 520               if (dirs[d] == DOWN && dy >= BEAM_EPS
 521                   && fabs (my_modf (y[d]) - straddle) < BEAM_EPS)
 522                 dem += extra_demerit;
 523             }
 524         }
 525       while (flip (&d) != LEFT);
 526     }
 527
 528   return dem;
 529 }
 530
 531