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