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