KiCad PCB EDA Suite
bezier_curves.cpp
Go to the documentation of this file.
1 /*
2  * This program source code file is part of KiCad, a free EDA CAD application.
3  *
4  * Copyright (C) 2014 Jean-Pierre Charras, jp.charras at wanadoo.fr
5  * Copyright (C) 2014 KiCad Developers, see CHANGELOG.TXT for contributors.
6  *
7  * This program is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU General Public License
9  * as published by the Free Software Foundation; either version 2
10  * of the License, or (at your option) any later version.
11  *
12  * This program 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 this program; if not, you may find one here:
19  * http://www.gnu.org/licenses/old-licenses/gpl-2.0.html
20  * or you may search the http://www.gnu.org website for the version 2 license,
21  * or you may write to the Free Software Foundation, Inc.,
22  * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA
23  */
24 
25 /************************************/
26 /* routines to handle bezier curves */
27 /************************************/
28 
29 #include <fctsys.h>
30 #include <bezier_curves.h>
31 
32 
33 #define add_segment(segment) if(s_bezier_Points_Buffer[s_bezier_Points_Buffer.size()-1] != segment) s_bezier_Points_Buffer.push_back(segment);
34 
35 
36 // Local variables:
37 static std::vector<wxPoint> s_bezier_Points_Buffer;
38 
39 static int bezier_recursion_limit = 12;
40 static double bezier_approximation_scale = 0.5; // 1
41 
42 static double bezier_curve_collinearity_epsilon = 1e-30;
43 static double bezier_curve_angle_tolerance_epsilon = 0.0001;
44 static double bezier_distance_tolerance_square; // derived by approximation_scale
45 static double bezier_angle_tolerance = 0.0;
46 static double bezier_cusp_limit = 0.0;
47 
48 // Local functions:
49 static void recursive_bezier( int x1, int y1, int x2, int y2, int x3, int y3, int level );
50 static void recursive_bezier( int x1,
51  int y1,
52  int x2,
53  int y2,
54  int x3,
55  int y3,
56  int x4,
57  int y4,
58  int level );
59 
60 /***********************************************************************************/
61 
62 
63 std::vector<wxPoint> Bezier2Poly( wxPoint c1, wxPoint c2, wxPoint c3, wxPoint c4 )
64 {
65  return Bezier2Poly( c1.x, c1.y, c2.x, c2.y, c3.x, c3.y, c4.x, c4.y );
66 }
67 
68 
69 std::vector<wxPoint> Bezier2Poly( wxPoint c1, wxPoint c2, wxPoint c3 )
70 {
71  return Bezier2Poly( c1.x, c1.y, c2.x, c2.y, c3.x, c3.y );
72 }
73 
74 
75 inline double calc_sq_distance( int x1, int y1, int x2, int y2 )
76 {
77  int dx = x2 - x1;
78  int dy = y2 - y1;
79 
80  return (double)dx * dx + (double)dy * dy;
81 }
82 
83 inline double sqrt_len( int dx, int dy )
84 {
85  return ((double)dx * dx) + ((double)dy * dy);
86 }
87 
88 
89 std::vector<wxPoint> Bezier2Poly( int x1, int y1, int x2, int y2, int x3, int y3 )
90 {
91  s_bezier_Points_Buffer.clear();
92 
95  s_bezier_Points_Buffer.push_back( wxPoint( x1, y1 ) );
96  recursive_bezier( x1, y1, x2, y2, x3, y3, 0 );
97  s_bezier_Points_Buffer.push_back( wxPoint( x3, y3 ) );
98 
99  wxLogDebug( wxT( "Bezier Conversion - End (%d vertex)" ), s_bezier_Points_Buffer.size() );
100  return s_bezier_Points_Buffer;
101 }
102 
103 
104 std::vector<wxPoint> Bezier2Poly( int x1, int y1, int x2, int y2, int x3, int y3, int x4, int y4 )
105 {
106  s_bezier_Points_Buffer.clear();
109 
110  s_bezier_Points_Buffer.push_back( wxPoint( x1, y1 ) );
111  recursive_bezier( x1, y1, x2, y2, x3, y3, x4, y4, 0 );
112  s_bezier_Points_Buffer.push_back( wxPoint( x4, y4 ) );
113  wxLogDebug( wxT( "Bezier Conversion - End (%d vertex)" ), s_bezier_Points_Buffer.size() );
114  return s_bezier_Points_Buffer;
115 }
116 
117 
118 void recursive_bezier( int x1, int y1, int x2, int y2, int x3, int y3, int level )
119 {
120  if( abs( level ) > bezier_recursion_limit )
121  {
122  return;
123  }
124 
125  // Calculate all the mid-points of the line segments
126  //----------------------
127  int x12 = (x1 + x2) / 2;
128  int y12 = (y1 + y2) / 2;
129  int x23 = (x2 + x3) / 2;
130  int y23 = (y2 + y3) / 2;
131  int x123 = (x12 + x23) / 2;
132  int y123 = (y12 + y23) / 2;
133 
134  int dx = x3 - x1;
135  int dy = y3 - y1;
136  double d = fabs( ((double) (x2 - x3) * dy) - ((double) (y2 - y3) * dx ) );
137  double da;
138 
140  {
141  // Regular case
142  //-----------------
143  if( d * d <= bezier_distance_tolerance_square * (dx * dx + dy * dy) )
144  {
145  // If the curvature doesn't exceed the distance_tolerance value
146  // we tend to finish subdivisions.
147  //----------------------
149  {
150  add_segment( wxPoint( x123, y123 ) );
151  return;
152  }
153 
154  // Angle & Cusp Condition
155  //----------------------
156  da = fabs( atan2( (double) ( y3 - y2 ), (double) ( x3 - x2 ) ) -
157  atan2( (double) ( y2 - y1 ), (double) ( x2 - x1 ) ) );
158  if( da >=M_PI )
159  da = 2 * M_PI - da;
160 
161  if( da < bezier_angle_tolerance )
162  {
163  // Finally we can stop the recursion
164  //----------------------
165  add_segment( wxPoint( x123, y123 ) );
166  return;
167  }
168  }
169  }
170  else
171  {
172  // Collinear case
173  //------------------
174  da = sqrt_len(dx, dy);
175  if( da == 0 )
176  {
177  d = calc_sq_distance( x1, y1, x2, y2 );
178  }
179  else
180  {
181  d = ( (double)(x2 - x1) * dx + (double)(y2 - y1) * dy ) / da;
182  if( d > 0 && d < 1 )
183  {
184  // Simple collinear case, 1---2---3
185  // We can leave just two endpoints
186  return;
187  }
188  if( d <= 0 )
189  d = calc_sq_distance( x2, y2, x1, y1 );
190  else if( d >= 1 )
191  d = calc_sq_distance( x2, y2, x3, y3 );
192  else
193  d = calc_sq_distance( x2, y2, x1 + (int) d * dx,
194  y1 + (int) d * dy );
195  }
197  {
198  add_segment( wxPoint( x2, y2 ) );
199  return;
200  }
201  }
202 
203  // Continue subdivision
204  //----------------------
205  recursive_bezier( x1, y1, x12, y12, x123, y123, level + 1 );
206  recursive_bezier( x123, y123, x23, y23, x3, y3, -(level + 1) );
207 }
208 
209 
210 void recursive_bezier( int x1, int y1, int x2, int y2, int x3, int y3, int x4, int y4, int level )
211 {
212  if( abs( level ) > bezier_recursion_limit )
213  {
214  return;
215  }
216 
217  // Calculate all the mid-points of the line segments
218  //----------------------
219  int x12 = (x1 + x2) / 2;
220  int y12 = (y1 + y2) / 2;
221  int x23 = (x2 + x3) / 2;
222  int y23 = (y2 + y3) / 2;
223  int x34 = (x3 + x4) / 2;
224  int y34 = (y3 + y4) / 2;
225  int x123 = (x12 + x23) / 2;
226  int y123 = (y12 + y23) / 2;
227  int x234 = (x23 + x34) / 2;
228  int y234 = (y23 + y34) / 2;
229  int x1234 = (x123 + x234) / 2;
230  int y1234 = (y123 + y234) / 2;
231 
232 
233  // Try to approximate the full cubic curve by a single straight line
234  //------------------
235  int dx = x4 - x1;
236  int dy = y4 - y1;
237 
238  double d2 = fabs( (double) ( (x2 - x4) * dy - (y2 - y4) * dx ) );
239  double d3 = fabs( (double) ( (x3 - x4) * dy - (y3 - y4) * dx ) );
240  double da1, da2, k;
241 
242  switch( (int(d2 > bezier_curve_collinearity_epsilon) << 1) +
244  {
245  case 0:
246 
247  // All collinear OR p1==p4
248  //----------------------
249  k = dx * dx + dy * dy;
250  if( k == 0 )
251  {
252  d2 = calc_sq_distance( x1, y1, x2, y2 );
253  d3 = calc_sq_distance( x4, y4, x3, y3 );
254  }
255  else
256  {
257  k = 1 / k;
258  da1 = x2 - x1;
259  da2 = y2 - y1;
260  d2 = k * (da1 * dx + da2 * dy);
261  da1 = x3 - x1;
262  da2 = y3 - y1;
263  d3 = k * (da1 * dx + da2 * dy);
264  if( d2 > 0 && d2 < 1 && d3 > 0 && d3 < 1 )
265  {
266  // Simple collinear case, 1---2---3---4
267  // We can leave just two endpoints
268  return;
269  }
270  if( d2 <= 0 )
271  d2 = calc_sq_distance( x2, y2, x1, y1 );
272  else if( d2 >= 1 )
273  d2 = calc_sq_distance( x2, y2, x4, y4 );
274  else
275  d2 = calc_sq_distance( x2, y2, x1 + (int) d2 * dx,
276  y1 + (int) d2 * dy );
277 
278  if( d3 <= 0 )
279  d3 = calc_sq_distance( x3, y3, x1, y1 );
280  else if( d3 >= 1 )
281  d3 = calc_sq_distance( x3, y3, x4, y4 );
282  else
283  d3 = calc_sq_distance( x3, y3, x1 + (int) d3 * dx,
284  y1 + (int) d3 * dy );
285  }
286  if( d2 > d3 )
287  {
289  {
290  add_segment( wxPoint( x2, y2 ) );
291  return;
292  }
293  }
294  else
295  {
297  {
298  add_segment( wxPoint( x3, y3 ) );
299  return;
300  }
301  }
302  break;
303 
304  case 1:
305 
306  // p1,p2,p4 are collinear, p3 is significant
307  //----------------------
308  if( d3 * d3 <= bezier_distance_tolerance_square * sqrt_len(dx, dy) )
309  {
311  {
312  add_segment( wxPoint( x23, y23 ) );
313  return;
314  }
315 
316  // Angle Condition
317  //----------------------
318  da1 = fabs( atan2( (double) ( y4 - y3 ), (double) ( x4 - x3 ) ) -
319  atan2( (double) ( y3 - y2 ), (double) ( x3 - x2 ) ) );
320  if( da1 >= M_PI )
321  da1 = 2 * M_PI - da1;
322 
323  if( da1 < bezier_angle_tolerance )
324  {
325  add_segment( wxPoint( x2, y2 ) );
326  add_segment( wxPoint( x3, y3 ) );
327  return;
328  }
329 
330  if( bezier_cusp_limit != 0.0 )
331  {
332  if( da1 > bezier_cusp_limit )
333  {
334  add_segment( wxPoint( x3, y3 ) );
335  return;
336  }
337  }
338  }
339  break;
340 
341  case 2:
342 
343  // p1,p3,p4 are collinear, p2 is significant
344  //----------------------
345  if( d2 * d2 <= bezier_distance_tolerance_square * sqrt_len(dx, dy) )
346  {
348  {
349  add_segment( wxPoint( x23, y23 ) );
350  return;
351  }
352 
353  // Angle Condition
354  //----------------------
355  da1 = fabs( atan2( (double) ( y3 - y2 ), (double) ( x3 - x2 ) ) -
356  atan2( (double) ( y2 - y1 ), (double) ( x2 - x1 ) ) );
357  if( da1 >= M_PI )
358  da1 = 2 * M_PI - da1;
359 
360  if( da1 < bezier_angle_tolerance )
361  {
362  add_segment( wxPoint( x2, y2 ) );
363  add_segment( wxPoint( x3, y3 ) );
364  return;
365  }
366 
367  if( bezier_cusp_limit != 0.0 )
368  {
369  if( da1 > bezier_cusp_limit )
370  {
371  add_segment( wxPoint( x2, y2 ) );
372  return;
373  }
374  }
375  }
376  break;
377 
378  case 3:
379 
380  // Regular case
381  //-----------------
382  if( (d2 + d3) * (d2 + d3) <= bezier_distance_tolerance_square * sqrt_len(dx, dy) )
383  {
384  // If the curvature doesn't exceed the distance_tolerance value
385  // we tend to finish subdivisions.
386  //----------------------
388  {
389  add_segment( wxPoint( x23, y23 ) );
390  return;
391  }
392 
393  // Angle & Cusp Condition
394  //----------------------
395  k = atan2( (double) ( y3 - y2 ), (double) ( x3 - x2 ) );
396  da1 = fabs( k - atan2( (double) ( y2 - y1 ),
397  (double) ( x2 - x1 ) ) );
398  da2 = fabs( atan2( (double) ( y4 - y3 ),
399  (double) ( x4 - x3 ) ) - k );
400  if( da1 >= M_PI )
401  da1 = 2 * M_PI - da1;
402  if( da2 >= M_PI )
403  da2 = 2 * M_PI - da2;
404 
405  if( da1 + da2 < bezier_angle_tolerance )
406  {
407  // Finally we can stop the recursion
408  //----------------------
409  add_segment( wxPoint( x23, y23 ) );
410  return;
411  }
412 
413  if( bezier_cusp_limit != 0.0 )
414  {
415  if( da1 > bezier_cusp_limit )
416  {
417  add_segment( wxPoint( x2, y2 ) );
418  return;
419  }
420 
421  if( da2 > bezier_cusp_limit )
422  {
423  add_segment( wxPoint( x3, y3 ) );
424  return;
425  }
426  }
427  }
428  break;
429  }
430 
431  // Continue subdivision
432  //----------------------
433  recursive_bezier( x1, y1, x12, y12, x123, y123, x1234, y1234, level + 1 );
434  recursive_bezier( x1234, y1234, x234, y234, x34, y34, x4, y4, level + 1 );
435 }
static double bezier_approximation_scale
static double bezier_cusp_limit
static double bezier_angle_tolerance
static int bezier_recursion_limit
static double bezier_curve_collinearity_epsilon
#define abs(a)
Definition: auxiliary.h:84
static double bezier_curve_angle_tolerance_epsilon
static std::vector< wxPoint > s_bezier_Points_Buffer
std::vector< wxPoint > Bezier2Poly(wxPoint c1, wxPoint c2, wxPoint c3, wxPoint c4)
Function Bezier2Poly convert a Bezier curve to a polyline.
#define add_segment(segment)
static double bezier_distance_tolerance_square
double calc_sq_distance(int x1, int y1, int x2, int y2)
double sqrt_len(int dx, int dy)
static void recursive_bezier(int x1, int y1, int x2, int y2, int x3, int y3, int level)