lily/beam-quanting.cc

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