2 bezier.cc -- implement Bezier and Bezier_bow
4 source file of the GNU LilyPond music typesetter
6 (c) 1998--2000 Jan Nieuwenhuizen <janneke@gnu.org>
12 #include "libc-extension.hh"
14 #include "polynomial.hh"
17 binomial_coefficient (Real over
, int under
)
23 x
*= over
/ Real (under
);
32 flip (Array
<Offset
>* arr_p
, Axis a
)
34 for (int i
= 0; i
< arr_p
->size (); i
++)
35 (*arr_p
)[i
][a
] = - (*arr_p
)[i
][a
];
39 rotate (Array
<Offset
>* arr_p
, Real phi
)
41 Offset
rot (complex_exp (Offset (0, phi
)));
42 for (int i
= 0; i
< arr_p
->size (); i
++)
43 (*arr_p
)[i
] = complex_multiply (rot
, (*arr_p
)[i
]);
47 translate (Array
<Offset
>* arr_p
, Offset o
)
49 for (int i
= 0; i
< arr_p
->size (); i
++)
55 Formula of the bezier 3-spline
57 sum_{j=0}^3 (3 over j) z_j (1-t)^(3-j) t^j
61 Bezier::get_other_coordinate (Axis a
, Real x
) const
63 Axis other
= Axis ((a
+1)%NO_AXES
);
64 Array
<Real
> ts
= solve_point (a
, x
);
66 Offset c
= curve_point (ts
[0]);
67 assert (fabs (c
[a
] - x
) < 1e-8);
74 Bezier::curve_point (Real t
)const
77 Real one_min_tj
= (1-t
)*(1-t
)*(1-t
);
80 for (int j
=0 ; j
< 4; j
++)
82 o
+= control_
[j
] * binomial_coefficient (3, j
)
83 * pow (t
,j
) * pow (1-t
, 3-j
);
91 assert (fabs (o
[X_AXIS
] - polynomial (X_AXIS
).eval (t
))< 1e-8);
92 assert (fabs (o
[Y_AXIS
] - polynomial (Y_AXIS
).eval (t
))< 1e-8);
100 Bezier::polynomial (Axis a
)const
103 for (int j
=0; j
<= 3; j
++)
106 * Polynomial::power (j
, Polynomial (0,1))*
107 Polynomial::power (3 - j
, Polynomial (1,-1))*
108 binomial_coefficient(3, j
);
115 Remove all numbers outside [0,1] from SOL
118 filter_solutions (Array
<Real
> sol
)
120 for (int i
= sol
.size (); i
--;)
121 if (sol
[i
] < 0 || sol
[i
] >1)
127 find t such that derivative is proportional to DERIV
130 Bezier::solve_derivative (Offset deriv
)const
132 Polynomial xp
=polynomial (X_AXIS
);
133 Polynomial yp
=polynomial (Y_AXIS
);
137 Polynomial combine
= xp
* deriv
[Y_AXIS
] - yp
* deriv
[X_AXIS
];
139 return filter_solutions (combine
.solve ());
144 Find t such that curve_point (t)[AX] == COORDINATE
147 Bezier::solve_point (Axis ax
, Real coordinate
) const
149 Polynomial
p(polynomial (ax
));
150 p
.coefs_
[0] -= coordinate
;
152 Array
<Real
> sol (p
.solve ());
153 return filter_solutions (sol
);
157 Compute the bounding box dimensions in direction of A.
160 Bezier::extent (Axis a
)const
162 int o
= (a
+1)%NO_AXES
;
166 Array
<Real
> sols (solve_derivative (d
));
169 for (int i
= sols
.size (); i
--;)
171 Offset
o (curve_point (sols
[i
]));
172 iv
.unite (Interval (o
[a
],o
[a
]));
182 Bezier::flip (Axis a
)
184 for (int i
= CONTROL_COUNT
; i
--;)
185 control_
[i
][a
] = - control_
[i
][a
];
189 Bezier::rotate (Real phi
)
191 Offset
rot (complex_exp (Offset (0, phi
)));
192 for (int i
= 0; i
< CONTROL_COUNT
; i
++)
193 control_
[i
] = complex_multiply (rot
, control_
[i
]);
197 Bezier::translate (Offset o
)
199 for (int i
= 0; i
< CONTROL_COUNT
; i
++)
204 Bezier::assert_sanity () const
206 for (int i
=0; i
< CONTROL_COUNT
; i
++)
207 assert (!isnan (control_
[i
].length ())
208 && !isinf (control_
[i
].length ()));
215 for (int i
=0; i
< CONTROL_COUNT
; i
++)
216 b2
.control_
[CONTROL_COUNT
-i
-1] = control_
[i
];