1 /* fit.c: turn a bitmap representation of a curve into a list of splines.
2 * Some of the ideas, but not the code, comes from the Phoenix thesis.
3 * See README for the reference.
4 *
5 * Copyright (C) 1992 Free Software Foundation, Inc.
6 *
7 * This program 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, or (at your option)
10 * 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, see <https://www.gnu.org/licenses/>.
19 */
20
21 #include "config.h"
22
23 #include <string.h>
24 #include <float.h>
25 #include <math.h>
26 #include <assert.h>
27
28 #include <glib.h>
29
30 #include "global.h"
31 #include "spline.h"
32 #include "vector.h"
33
34 #include "curve.h"
35 #include "fit.h"
36 #include "pxl-outline.h"
37
38 /* If two endpoints are closer than this, they are made to be equal.
39 (-align-threshold) */
40 real align_threshold = 0.5;
41
42 /* If the angle defined by a point and its predecessors and successors
43 is smaller than this, it's a corner, even if it's within
44 `corner_surround' pixels of a point with a smaller angle.
45 (-corner-always-threshold) */
46 real corner_always_threshold = 60.0;
47
48 /* Number of points to consider when determining if a point is a corner
49 or not. (-corner-surround) */
50 unsigned corner_surround = 4;
51
52 /* If a point, its predecessors, and its successors define an angle
53 smaller than this, it's a corner. Should be in range 0..180.
54 (-corner-threshold) */
55 real corner_threshold = 100.0;
56
57 /* Amount of error at which a fitted spline is unacceptable. If any
58 pixel is further away than this from the fitted curve, we try again.
59 (-error-threshold) */
60 /* real error_threshold = .8; ALT */
61 real error_threshold = .4;
62
63 /* A second number of adjacent points to consider when filtering.
64 (-filter-alternative-surround) */
65 unsigned filter_alternative_surround = 1;
66
67 /* If the angles between the vectors produced by filter_surround and
68 filter_alternative_surround points differ by more than this, use
69 the one from filter_alternative_surround. (-filter-epsilon) */
70 real filter_epsilon = 10.0;
71
72 /* Number of times to smooth original data points. Increasing this
73 number dramatically---to 50 or so---can produce vastly better
74 results. But if any points that ``should'' be corners aren't found,
75 the curve goes to hell around that point. (-filter-iterations) */
76 /* unsigned filter_iteration_count = 4; ALT */
77 unsigned filter_iteration_count = 4;
78
79 /* To produce the new point, use the old point plus this times the
80 neighbors. (-filter-percent) */
81 real filter_percent = .33;
82
83 /* Number of adjacent points to consider if `filter_surround' points
84 defines a straight line. (-filter-secondary-surround) */
85 static unsigned filter_secondary_surround = 3;
86
87 /* Number of adjacent points to consider when filtering.
88 (-filter-surround) */
89 unsigned filter_surround = 2;
90
91 /* Says whether or not to remove ``knee'' points after finding the outline.
92 (See the comments at `remove_knee_points'.) (-remove-knees). */
93 boolean keep_knees = false;
94
95 /* If a spline is closer to a straight line than this, it remains a
96 straight line, even if it would otherwise be changed back to a curve.
97 This is weighted by the square of the curve length, to make shorter
98 curves more likely to be reverted. (-line-reversion-threshold) */
99 real line_reversion_threshold = .01;
100
101 /* How many pixels (on the average) a spline can diverge from the line
102 determined by its endpoints before it is changed to a straight line.
103 (-line-threshold) */
104 /* real line_threshold = 1.0; ALT */
105 real line_threshold = 0.5;
106
107 /* If reparameterization doesn't improve the fit by this much percent,
108 stop doing it. (-reparameterize-improve) */
109 /* real reparameterize_improvement = .10; ALT */
110 real reparameterize_improvement = .01;
111
112 /* Amount of error at which it is pointless to reparameterize. This
113 happens, for example, when we are trying to fit the outline of the
114 outside of an `O' with a single spline. The initial fit is not good
115 enough for the Newton-Raphson iteration to improve it. It may be
116 that it would be better to detect the cases where we didn't find any
117 corners. (-reparameterize-threshold) */
118 /* real reparameterize_threshold = 30.0; ALT */
119 real reparameterize_threshold = 1.0;
120
121 /* Percentage of the curve away from the worst point to look for a
122 better place to subdivide. (-subdivide-search) */
123 real subdivide_search = .1;
124
125 /* Number of points to consider when deciding whether a given point is a
126 better place to subdivide. (-subdivide-surround) */
127 unsigned subdivide_surround = 4;
128
129 /* How many pixels a point can diverge from a straight line and still be
130 considered a better place to subdivide. (-subdivide-threshold) */
131 real subdivide_threshold = .03;
132
133 /* Number of points to look at on either side of a point when computing
134 the approximation to the tangent at that point. (-tangent-surround) */
135 unsigned tangent_surround = 3;
136
137
138 /* We need to manipulate lists of array indices. */
139
140 typedef struct index_list
141 {
142 unsigned *data;
143 unsigned length;
144 } index_list_type;
145
146 /* The usual accessor macros. */
147 #define GET_INDEX(i_l, n) ((i_l).data[n])
148 #define INDEX_LIST_LENGTH(i_l) ((i_l).length)
149 #define GET_LAST_INDEX(i_l) ((i_l).data[INDEX_LIST_LENGTH (i_l) - 1])
150
151 static void append_index (index_list_type *, unsigned);
152 static void free_index_list (index_list_type *);
153 static index_list_type new_index_list (void);
154 static void remove_adjacent_corners (index_list_type *, unsigned);
155
156 static void align (spline_list_type *);
157 static void change_bad_lines (spline_list_type *);
158 static void filter (curve_type);
159 static real filter_angle (vector_type, vector_type);
160 static void find_curve_vectors
161 (unsigned, curve_type, unsigned, vector_type *, vector_type *, unsigned *);
162 static unsigned find_subdivision (curve_type, unsigned initial);
163 static void find_vectors
164 (unsigned, pixel_outline_type, vector_type *, vector_type *);
165 static index_list_type find_corners (pixel_outline_type);
166 static real find_error (curve_type, spline_type, unsigned *);
167 static vector_type find_half_tangent (curve_type, boolean start, unsigned *);
168 static void find_tangent (curve_type, boolean, boolean);
169 static spline_type fit_one_spline (curve_type);
170 static spline_list_type *fit_curve (curve_type);
171 static spline_list_type fit_curve_list (curve_list_type);
172 static spline_list_type *fit_with_least_squares (curve_type);
173 static spline_list_type *fit_with_line (curve_type);
174 static void remove_knee_points (curve_type, boolean);
175 static boolean reparameterize (curve_type, spline_type);
176 static void set_initial_parameter_values (curve_type);
177 static boolean spline_linear_enough (spline_type *, curve_type);
178 static curve_list_array_type split_at_corners (pixel_outline_list_type);
179 static boolean test_subdivision_point (curve_type, unsigned, vector_type *);
180
181 /* The top-level call that transforms the list of pixels in the outlines
182 of the original character to a list of spline lists fitted to those
183 pixels. */
184
185 spline_list_array_type
fitted_splines(pixel_outline_list_type pixel_outline_list)186 fitted_splines (pixel_outline_list_type pixel_outline_list)
187 {
188 unsigned this_list;
189 unsigned total = 0;
190 spline_list_array_type char_splines = new_spline_list_array ();
191 curve_list_array_type curve_array = split_at_corners (pixel_outline_list);
192
193 for (this_list = 0; this_list < CURVE_LIST_ARRAY_LENGTH (curve_array);
194 this_list++)
195 {
196 spline_list_type curve_list_splines;
197 curve_list_type curves = CURVE_LIST_ARRAY_ELT (curve_array, this_list);
198
199 curve_list_splines = fit_curve_list (curves);
200 append_spline_list (&char_splines, curve_list_splines);
201
202 /* REPORT ("* "); */
203 }
204
205 free_curve_list_array (&curve_array);
206
207 for (this_list = 0; this_list < SPLINE_LIST_ARRAY_LENGTH (char_splines);
208 this_list++)
209 total
210 += SPLINE_LIST_LENGTH (SPLINE_LIST_ARRAY_ELT (char_splines, this_list));
211
212 /* REPORT1 ("=%u", total); */
213
214 return char_splines;
215 }
216
217 /* Set up the internal parameters from the external ones */
218
219 void
fit_set_params(SELVALS * selVals)220 fit_set_params(SELVALS *selVals)
221 {
222 align_threshold = selVals->align_threshold;
223 corner_always_threshold = selVals->corner_always_threshold;
224 corner_surround = selVals->corner_surround;
225 corner_threshold = selVals->corner_threshold;
226 error_threshold = selVals->error_threshold;
227 filter_alternative_surround = selVals->filter_alternative_surround;
228 filter_epsilon = selVals->filter_epsilon;
229 filter_iteration_count = selVals->filter_iteration_count;
230 filter_percent = selVals->filter_percent;
231 filter_secondary_surround = selVals->filter_secondary_surround;
232 filter_surround = selVals->filter_surround;
233 keep_knees = selVals->keep_knees;
234 line_reversion_threshold = selVals->line_reversion_threshold;
235 line_threshold = selVals->line_threshold;
236 reparameterize_improvement = selVals->reparameterize_improvement;
237 reparameterize_threshold = selVals->reparameterize_threshold;
238 subdivide_search = selVals->subdivide_search;
239 subdivide_surround = selVals->subdivide_surround;
240 subdivide_threshold = selVals->subdivide_threshold;
241 tangent_surround = selVals->tangent_surround;
242 }
243
244 void
fit_set_default_params(SELVALS * selVals)245 fit_set_default_params(SELVALS *selVals)
246 {
247 selVals->align_threshold = align_threshold;
248 selVals->corner_always_threshold = corner_always_threshold;
249 selVals->corner_surround = corner_surround;
250 selVals->corner_threshold = corner_threshold;
251 selVals->error_threshold = error_threshold;
252 selVals->filter_alternative_surround = filter_alternative_surround;
253 selVals->filter_epsilon = filter_epsilon;
254 selVals->filter_iteration_count = filter_iteration_count;
255 selVals->filter_percent = filter_percent;
256 selVals->filter_secondary_surround = filter_secondary_surround;
257 selVals->filter_surround = filter_surround;
258 selVals->keep_knees = keep_knees;
259 selVals->line_reversion_threshold = line_reversion_threshold;
260 selVals->line_threshold = line_threshold;
261 selVals->reparameterize_improvement = reparameterize_improvement;
262 selVals->reparameterize_threshold = reparameterize_threshold;
263 selVals->subdivide_search = subdivide_search;
264 selVals->subdivide_surround = subdivide_surround;
265 selVals->subdivide_threshold = subdivide_threshold;
266 selVals->tangent_surround = tangent_surround;
267 }
268
269 /* Fit the list of curves CURVE_LIST to a list of splines, and return
270 it. CURVE_LIST represents a single closed paths, e.g., either the
271 inside or outside outline of an `o'. */
272
273 static spline_list_type
fit_curve_list(curve_list_type curve_list)274 fit_curve_list (curve_list_type curve_list)
275 {
276 curve_type curve;
277 unsigned this_curve, this_spline;
278 unsigned curve_list_length = CURVE_LIST_LENGTH (curve_list);
279 spline_list_type curve_list_splines = *new_spline_list ();
280
281 /* Remove the extraneous ``knee'' points before filtering. Since the
282 corners have already been found, we don't need to worry about
283 removing a point that should be a corner. */
284 if (!keep_knees)
285 {
286 /* LOG ("\nRemoving knees:\n"); */
287 for (this_curve = 0; this_curve < curve_list_length; this_curve++)
288 {
289 /* LOG1 ("#%u:", this_curve); */
290 remove_knee_points (CURVE_LIST_ELT (curve_list, this_curve),
291 CURVE_LIST_CLOCKWISE (curve_list));
292 }
293 }
294
295 /* We filter all the curves in CURVE_LIST at once; otherwise, we would
296 look at an unfiltered curve when computing tangents. */
297 /* LOG ("\nFiltering curves:\n"); */
298 for (this_curve = 0; this_curve < curve_list.length; this_curve++)
299 {
300 /* LOG1 ("#%u: ", this_curve); */
301 filter (CURVE_LIST_ELT (curve_list, this_curve));
302 /* REPORT ("f"); */
303 }
304
305 /* Make the first point in the first curve also be the last point in
306 the last curve, so the fit to the whole curve list will begin and
307 end at the same point. This may cause slight errors in computing
308 the tangents and t values, but it's worth it for the continuity.
309 Of course we don't want to do this if the two points are already
310 the same, as they are if the curve is cyclic. (We don't append it
311 earlier, in `split_at_corners', because that confuses the
312 filtering.) Finally, we can't append the point if the curve is
313 exactly three points long, because we aren't adding any more data,
314 and three points isn't enough to determine a spline. Therefore,
315 the fitting will fail. */
316 curve = CURVE_LIST_ELT (curve_list, 0);
317 if (CURVE_CYCLIC (curve) && CURVE_LENGTH (curve) != 3)
318 append_point (curve, CURVE_POINT (curve, 0));
319
320 /* Finally, fit each curve in the list to a list of splines. */
321 for (this_curve = 0; this_curve < curve_list_length; this_curve++)
322 {
323 spline_list_type *curve_splines;
324 curve_type current_curve = CURVE_LIST_ELT (curve_list, this_curve);
325
326 /* REPORT1 (" %u", this_curve); */
327 /* LOG1 ("\nFitting curve #%u:\n", this_curve); */
328
329 curve_splines = fit_curve (current_curve);
330 if (curve_splines == NULL)
331 printf("Could not fit curve #%u", this_curve);
332 else
333 {
334 /* LOG1 ("Fitted splines for curve #%u:\n", this_curve); */
335 for (this_spline = 0;
336 this_spline < SPLINE_LIST_LENGTH (*curve_splines);
337 this_spline++)
338 {
339 /* LOG1 (" %u: ", this_spline); */
340 /* if (logging) */
341 /* print_spline (log_
342 file, */
343 /* SPLINE_LIST_ELT (*curve_splines, this_spline)); */
344 }
345
346 /* After fitting, we may need to change some would-be lines
347 back to curves, because they are in a list with other
348 curves. */
349 change_bad_lines (curve_splines);
350
351 concat_spline_lists (&curve_list_splines, *curve_splines);
352 /* REPORT1 ("(%u)", SPLINE_LIST_LENGTH (*curve_splines)); */
353 }
354 }
355
356
357 /* We do this for each outline's spline list because when a point
358 is changed, it needs to be changed in both segments in which it
359 appears---and the segments might be in different curves. */
360 align (&curve_list_splines);
361
362 return curve_list_splines;
363 }
364
365
366 /* Transform a set of locations to a list of splines (the fewer the
367 better). We are guaranteed that CURVE does not contain any corners.
368 We return NULL if we cannot fit the points at all. */
369
370 static spline_list_type *
fit_curve(curve_type curve)371 fit_curve (curve_type curve)
372 {
373 spline_list_type *fitted_splines;
374
375 if (CURVE_LENGTH (curve) < 2)
376 {
377 printf ("Tried to fit curve with less than two points");
378 return NULL;
379 }
380
381 /* Do we have enough points to fit with a spline? */
382 fitted_splines = CURVE_LENGTH (curve) < 4
383 ? fit_with_line (curve)
384 : fit_with_least_squares (curve);
385
386 return fitted_splines;
387 }
388
389 /* As mentioned above, the first step is to find the corners in
390 PIXEL_LIST, the list of points. (Presumably we can't fit a single
391 spline around a corner.) The general strategy is to look through all
392 the points, remembering which we want to consider corners. Then go
393 through that list, producing the curve_list. This is dictated by the
394 fact that PIXEL_LIST does not necessarily start on a corner---it just
395 starts at the character's first outline pixel, going left-to-right,
396 top-to-bottom. But we want all our splines to start and end on real
397 corners.
398
399 For example, consider the top of a capital `C' (this is in cmss20):
400 x
401 ***********
402 ******************
403
404 PIXEL_LIST will start at the pixel below the `x'. If we considered
405 this pixel a corner, we would wind up matching a very small segment
406 from there to the end of the line, probably as a straight line, which
407 is certainly not what we want.
408
409 PIXEL_LIST has one element for each closed outline on the character.
410 To preserve this information, we return an array of curve_lists, one
411 element (which in turn consists of several curves, one between each
412 pair of corners) for each element in PIXEL_LIST. */
413
414 static curve_list_array_type
split_at_corners(pixel_outline_list_type pixel_list)415 split_at_corners (pixel_outline_list_type pixel_list)
416 {
417 unsigned this_pixel_o;
418 curve_list_array_type curve_array = new_curve_list_array ();
419
420 /* LOG ("\nFinding corners:\n"); */
421
422 for (this_pixel_o = 0; this_pixel_o < O_LIST_LENGTH (pixel_list);
423 this_pixel_o++)
424 {
425 curve_type curve, first_curve;
426 index_list_type corner_list;
427 unsigned p, this_corner;
428 curve_list_type curve_list = new_curve_list ();
429 pixel_outline_type pixel_o = O_LIST_OUTLINE (pixel_list, this_pixel_o);
430
431 CURVE_LIST_CLOCKWISE (curve_list) = O_CLOCKWISE (pixel_o);
432
433 /* LOG1 ("#%u:", this_pixel_o); */
434
435 /* If the outline does not have enough points, we can't do
436 anything. The endpoints of the outlines are automatically
437 corners. We need at least `corner_surround' more pixels on
438 either side of a point before it is conceivable that we might
439 want another corner. */
440 if (O_LENGTH (pixel_o) > corner_surround * 2 + 2)
441 {
442 corner_list = find_corners (pixel_o);
443 }
444 else
445 {
446 corner_list.data = NULL;
447 corner_list.length = 0;
448 }
449
450 /* Remember the first curve so we can make it be the `next' of the
451 last one. (And vice versa.) */
452 first_curve = new_curve ();
453
454 curve = first_curve;
455
456 if (corner_list.length == 0)
457 { /* No corners. Use all of the pixel outline as the curve. */
458 for (p = 0; p < O_LENGTH (pixel_o); p++)
459 append_pixel (curve, O_COORDINATE (pixel_o, p));
460
461 /* This curve is cyclic. */
462 CURVE_CYCLIC (curve) = true;
463 }
464 else
465 { /* Each curve consists of the points between (inclusive) each pair
466 of corners. */
467 for (this_corner = 0; this_corner < corner_list.length - 1;
468 this_corner++)
469 {
470 curve_type previous_curve = curve;
471 unsigned corner = GET_INDEX (corner_list, this_corner);
472 unsigned next_corner = GET_INDEX (corner_list, this_corner + 1);
473
474 for (p = corner; p <= next_corner; p++)
475 append_pixel (curve, O_COORDINATE (pixel_o, p));
476
477 append_curve (&curve_list, curve);
478 curve = new_curve ();
479 NEXT_CURVE (previous_curve) = curve;
480 PREVIOUS_CURVE (curve) = previous_curve;
481 }
482
483 /* The last curve is different. It consists of the points
484 (inclusive) between the last corner and the end of the list,
485 and the beginning of the list and the first corner. */
486 for (p = GET_LAST_INDEX (corner_list); p < O_LENGTH (pixel_o);
487 p++)
488 append_pixel (curve, O_COORDINATE (pixel_o, p));
489
490 for (p = 0; p <= GET_INDEX (corner_list, 0); p++)
491 append_pixel (curve, O_COORDINATE (pixel_o, p));
492 }
493
494 /* LOG1 (" [%u].\n", corner_list.length); */
495
496 /* Add `curve' to the end of the list, updating the pointers in
497 the chain. */
498 append_curve (&curve_list, curve);
499 NEXT_CURVE (curve) = first_curve;
500 PREVIOUS_CURVE (first_curve) = curve;
501
502 /* And now add the just-completed curve list to the array. */
503 append_curve_list (&curve_array, curve_list);
504 } /* End of considering each pixel outline. */
505
506 return curve_array;
507 }
508
509
510 /* We consider a point to be a corner if (1) the angle defined by the
511 `corner_surround' points coming into it and going out from it is less
512 than `corner_threshold' degrees, and no point within
513 `corner_surround' points has a smaller angle; or (2) the angle is less
514 than `corner_always_threshold' degrees.
515
516 Because of the different cases, it is convenient to have the
517 following macro to append a corner on to the list we return. The
518 character argument C is simply so that the different cases can be
519 distinguished in the log file. */
520
521 #define APPEND_CORNER(index, angle, c) \
522 do \
523 { \
524 append_index (&corner_list, index); \
525 /*LOG4 (" (%d,%d)%c%.3f", */ \
526 /* O_COORDINATE (pixel_outline, index).x,*/ \
527 /* O_COORDINATE (pixel_outline, index).y,*/ \
528 /* c, angle);*/ \
529 } \
530 while (0)
531
532 static index_list_type
find_corners(pixel_outline_type pixel_outline)533 find_corners (pixel_outline_type pixel_outline)
534 {
535 unsigned p;
536 index_list_type corner_list = new_index_list ();
537
538 /* Consider each pixel on the outline in turn. */
539 for (p = 0; p < O_LENGTH (pixel_outline); p++)
540 {
541 real corner_angle;
542 vector_type in_vector, out_vector;
543
544 /* Check if the angle is small enough. */
545 find_vectors (p, pixel_outline, &in_vector, &out_vector);
546 corner_angle = Vangle (in_vector, out_vector);
547
548 if (fabs (corner_angle) <= corner_threshold)
549 {
550 /* We want to keep looking, instead of just appending the
551 first pixel we find with a small enough angle, since there
552 might be another corner within `corner_surround' pixels, with
553 a smaller angle. If that is the case, we want that one. */
554 real best_corner_angle = corner_angle;
555 unsigned best_corner_index = p;
556 index_list_type equally_good_list = new_index_list ();
557 /* As we come into the loop, `p' is the index of the point
558 that has an angle less than `corner_angle'. We use `i' to
559 move through the pixels next to that, and `q' for moving
560 through the adjacent pixels to each `p'. */
561 unsigned q = p;
562 unsigned i = p + 1;
563
564 while (true)
565 {
566 /* Perhaps the angle is sufficiently small that we want to
567 consider this a corner, even if it's not the best
568 (unless we've already wrapped around in the search,
569 i.e., `q<i', in which case we have already added the
570 corner, and we don't want to add it again). We want to
571 do this check on the first candidate we find, as well
572 as the others in the loop, hence this comes before the
573 stopping condition. */
574 if (corner_angle <= corner_always_threshold && q >= p)
575 APPEND_CORNER (q, corner_angle, '\\');
576
577 /* Exit the loop if we've looked at `corner_surround'
578 pixels past the best one we found, or if we've looked
579 at all the pixels. */
580 if (i >= best_corner_index + corner_surround
581 || i >= O_LENGTH (pixel_outline))
582 break;
583
584 /* Check the angle. */
585 q = i % O_LENGTH (pixel_outline);
586 find_vectors (q, pixel_outline, &in_vector, &out_vector);
587 corner_angle = Vangle (in_vector, out_vector);
588
589 /* If we come across a corner that is just as good as the
590 best one, we should make it a corner, too. This
591 happens, for example, at the points on the `W' in some
592 typefaces, where the ``points'' are flat. */
593 if (epsilon_equal (corner_angle, best_corner_angle))
594 append_index (&equally_good_list, q);
595
596 else if (corner_angle < best_corner_angle)
597 {
598 best_corner_angle = corner_angle;
599 /* We want to check `corner_surround' pixels beyond the
600 new best corner. */
601 i = best_corner_index = q;
602 free_index_list (&equally_good_list);
603 equally_good_list = new_index_list ();
604 }
605
606 i++;
607 }
608
609 /* After we exit the loop, `q' is the index of the last point
610 we checked. We have already added the corner if
611 `best_corner_angle' is less than `corner_always_threshold'.
612 Again, if we've already wrapped around, we don't want to
613 add the corner again. */
614 if (best_corner_angle > corner_always_threshold
615 && best_corner_index >= p)
616 {
617 unsigned i;
618
619 APPEND_CORNER (best_corner_index, best_corner_angle, '/');
620
621 for (i = 0; i < INDEX_LIST_LENGTH (equally_good_list); i++)
622 APPEND_CORNER (GET_INDEX (equally_good_list, i),
623 best_corner_angle, '@');
624 free_index_list (&equally_good_list);
625 }
626
627 /* If we wrapped around in our search, we're done; otherwise,
628 we don't want the outer loop to look at the pixels that we
629 already looked at in searching for the best corner. */
630 p = (q < p) ? O_LENGTH (pixel_outline) : q;
631 } /* End of searching for the best corner. */
632 } /* End of considering each pixel. */
633
634 if (INDEX_LIST_LENGTH (corner_list) > 0)
635 /* We never want two corners next to each other, since the
636 only way to fit such a ``curve'' would be with a straight
637 line, which usually interrupts the continuity dreadfully. */
638 remove_adjacent_corners (&corner_list, O_LENGTH (pixel_outline) - 1);
639
640 return corner_list;
641 }
642
643
644 /* Return the difference vectors coming in and going out of the outline
645 OUTLINE at the point whose index is TEST_INDEX. In Phoenix,
646 Schneider looks at a single point on either side of the point we're
647 considering. That works for him because his points are not touching.
648 But our points *are* touching, and so we have to look at
649 `corner_surround' points on either side, to get a better picture of
650 the outline's shape. */
651
652 static void
find_vectors(unsigned test_index,pixel_outline_type outline,vector_type * in,vector_type * out)653 find_vectors (unsigned test_index, pixel_outline_type outline,
654 vector_type *in, vector_type *out)
655 {
656 int i;
657 unsigned n_done;
658 coordinate_type candidate = O_COORDINATE (outline, test_index);
659
660 in->dx = 0.0;
661 in->dy = 0.0;
662 out->dx = 0.0;
663 out->dy = 0.0;
664
665 /* Add up the differences from p of the `corner_surround' points
666 before p. */
667 for (i = O_PREV (outline, test_index), n_done = 0; n_done < corner_surround;
668 i = O_PREV (outline, i), n_done++)
669 *in = Vadd (*in, IPsubtract (O_COORDINATE (outline, i), candidate));
670
671 #if 0
672 /* We don't need this code any more, because now we create the pixel
673 outline from the corners of the pixels, rather than the edges. */
674
675 /* To see why we need this test, consider the following
676 case: four pixels stacked vertically with no other adjacent pixels,
677 i.e., *
678 *x
679 *
680 *
681 *** (etc.) We are considering the pixel marked `x' for cornerhood.
682 The out vector at this point is going to be the zero vector (if
683 `corner_surround' is 3), because the first
684 pixel on the outline is the one above the x, the second pixel x
685 itself, and the third the one below x. (Remember that we go
686 around the edges of the pixels to find the outlines, not the
687 pixels themselves.) */
688 if (magnitude (*in) == 0.0)
689 {
690 WARNING ("Zero magnitude in");
691 return corner_threshold + 1.0;
692 }
693 #endif /* 0 */
694
695 /* And the points after p. */
696 for (i = O_NEXT (outline, test_index), n_done = 0; n_done < corner_surround;
697 i = O_NEXT (outline, i), n_done++)
698 *out = Vadd (*out, IPsubtract (O_COORDINATE (outline, i), candidate));
699
700 #if 0
701 /* As with the test for the in vector, we don't need this any more. */
702 if (magnitude (*out) == 0.0)
703 {
704 WARNING ("Zero magnitude out");
705 return corner_threshold + 1.0;
706 }
707 #endif /* 0 */
708 }
709
710
711 /* Remove adjacent points from the index list LIST. We do this by first
712 sorting the list and then running through it. Since these lists are
713 quite short, a straight selection sort (e.g., p.139 of the Art of
714 Computer Programming, vol.3) is good enough. LAST_INDEX is the index
715 of the last pixel on the outline, i.e., the next one is the first
716 pixel. We need this for checking the adjacency of the last corner.
717
718 We need to do this because the adjacent corners turn into
719 two-pixel-long curves, which can only be fit by straight lines. */
720
721 static void
remove_adjacent_corners(index_list_type * list,unsigned last_index)722 remove_adjacent_corners (index_list_type *list, unsigned last_index)
723 {
724 unsigned j;
725 unsigned last;
726 index_list_type new = new_index_list ();
727
728 for (j = INDEX_LIST_LENGTH (*list) - 1; j > 0; j--)
729 {
730 unsigned search;
731 unsigned temp;
732 /* Find maximal element below `j'. */
733 unsigned max_index = j;
734
735 for (search = 0; search < j; search++)
736 if (GET_INDEX (*list, search) > GET_INDEX (*list, max_index))
737 max_index = search;
738
739 if (max_index != j)
740 {
741 temp = GET_INDEX (*list, j);
742 GET_INDEX (*list, j) = GET_INDEX (*list, max_index);
743 GET_INDEX (*list, max_index) = temp;
744 printf ("needed exchange"); /* xx -- really have to sort? */
745 }
746 }
747
748 /* The list is sorted. Now look for adjacent entries. Each time
749 through the loop we insert the current entry and, if appropriate,
750 the next entry. */
751 for (j = 0; j < INDEX_LIST_LENGTH (*list) - 1; j++)
752 {
753 unsigned current = GET_INDEX (*list, j);
754 unsigned next = GET_INDEX (*list, j + 1);
755
756 /* We should never have inserted the same element twice. */
757 assert (current != next);
758
759 append_index (&new, current);
760 if (next == current + 1)
761 j++;
762 }
763
764 /* Don't append the last element if it is 1) adjacent to the previous
765 one; or 2) adjacent to the very first one. */
766 last = GET_LAST_INDEX (*list);
767 if (INDEX_LIST_LENGTH (new) == 0
768 || !(last == GET_LAST_INDEX (new) + 1
769 || (last == last_index && GET_INDEX (*list, 0) == 0)))
770 append_index (&new, last);
771
772 free_index_list (list);
773 *list = new;
774 }
775
776 /* A ``knee'' is a point which forms a ``right angle'' with its
777 predecessor and successor. See the documentation (the `Removing
778 knees' section) for an example and more details.
779
780 The argument CLOCKWISE tells us which direction we're moving. (We
781 can't figure that information out from just the single segment with
782 which we are given to work.)
783
784 We should never find two consecutive knees.
785
786 Since the first and last points are corners (unless the curve is
787 cyclic), it doesn't make sense to remove those. */
788
789 /* This evaluates to true if the vector V is zero in one direction and
790 nonzero in the other. */
791 #define ONLY_ONE_ZERO(v) \
792 (((v).dx == 0.0 && (v).dy != 0.0) || ((v).dy == 0.0 && (v).dx != 0.0))
793
794
795 /* There are four possible cases for knees, one for each of the four
796 corners of a rectangle; and then the cases differ depending on which
797 direction we are going around the curve. The tests are listed here
798 in the order of upper left, upper right, lower right, lower left.
799 Perhaps there is some simple pattern to the
800 clockwise/counterclockwise differences, but I don't see one. */
801 #define CLOCKWISE_KNEE(prev_delta, next_delta) \
802 ((prev_delta.dx == -1.0 && next_delta.dy == 1.0) \
803 || (prev_delta.dy == 1.0 && next_delta.dx == 1.0) \
804 || (prev_delta.dx == 1.0 && next_delta.dy == -1.0) \
805 || (prev_delta.dy == -1.0 && next_delta.dx == -1.0))
806
807 #define COUNTERCLOCKWISE_KNEE(prev_delta, next_delta) \
808 ((prev_delta.dy == 1.0 && next_delta.dx == -1.0) \
809 || (prev_delta.dx == 1.0 && next_delta.dy == 1.0) \
810 || (prev_delta.dy == -1.0 && next_delta.dx == 1.0) \
811 || (prev_delta.dx == -1.0 && next_delta.dy == -1.0))
812
813 static void
remove_knee_points(curve_type curve,boolean clockwise)814 remove_knee_points (curve_type curve, boolean clockwise)
815 {
816 int i;
817 unsigned offset = CURVE_CYCLIC (curve) ? 0 : 1;
818 coordinate_type previous
819 = real_to_int_coord (CURVE_POINT (curve, CURVE_PREV (curve, offset)));
820 curve_type trimmed_curve = copy_most_of_curve (curve);
821
822 if (!CURVE_CYCLIC (curve))
823 append_pixel (trimmed_curve, real_to_int_coord (CURVE_POINT (curve, 0)));
824
825 for (i = offset; i < CURVE_LENGTH (curve) - offset; i++)
826 {
827 coordinate_type current
828 = real_to_int_coord (CURVE_POINT (curve, i));
829 coordinate_type next
830 = real_to_int_coord (CURVE_POINT (curve, CURVE_NEXT (curve, i)));
831 vector_type prev_delta = IPsubtract (previous, current);
832 vector_type next_delta = IPsubtract (next, current);
833
834 if (ONLY_ONE_ZERO (prev_delta) && ONLY_ONE_ZERO (next_delta)
835 && ((clockwise && CLOCKWISE_KNEE (prev_delta, next_delta))
836 || (!clockwise
837 && COUNTERCLOCKWISE_KNEE (prev_delta, next_delta))))
838 {
839 /* LOG2 (" (%d,%d)", current.x, current.y); */
840 }
841 else
842 {
843 previous = current;
844 append_pixel (trimmed_curve, current);
845 }
846 }
847
848 if (!CURVE_CYCLIC (curve))
849 append_pixel (trimmed_curve, real_to_int_coord (LAST_CURVE_POINT (curve)));
850
851 /* if (CURVE_LENGTH (trimmed_curve) == CURVE_LENGTH (curve)) */
852 /* LOG (" (none)"); */
853
854 /* LOG (".\n"); */
855
856 free_curve (curve);
857 *curve = *trimmed_curve;
858 }
859
860 /* Smooth the curve by adding in neighboring points. Do this
861 `filter_iteration_count' times. But don't change the corners. */
862
863 #if 0
864 /* Computing the new point based on a single neighboring point and with
865 constant percentages, as the `SMOOTH' macro did, isn't quite good
866 enough. For example, at the top of an `o' the curve might well have
867 three consecutive horizontal pixels, even though there isn't really a
868 straight there. With this code, the middle point would remain
869 unfiltered. */
870
871 #define SMOOTH(axis) \
872 CURVE_POINT (curve, this_point).axis \
873 = ((1.0 - center_percent) / 2) \
874 * CURVE_POINT (curve, CURVE_PREV (curve, this_point)).axis \
875 + center_percent * CURVE_POINT (curve, this_point).axis \
876 + ((1.0 - center_percent) / 2) \
877 * CURVE_POINT (curve, CURVE_NEXT (curve, this_point)).axis
878 #endif /* 0 */
879
880 /* We sometimes need to change the information about the filtered point.
881 This macro assigns to the relevant variables. */
882 #define FILTER_ASSIGN(new) \
883 do \
884 { \
885 in = in_##new; \
886 out = out_##new; \
887 count = new##_count; \
888 angle = angle_##new; \
889 } \
890 while (0)
891
892 static void
filter(curve_type curve)893 filter (curve_type curve)
894 {
895 unsigned iteration, this_point;
896 unsigned offset = CURVE_CYCLIC (curve) ? 0 : 1;
897
898 /* We must have at least three points---the previous one, the current
899 one, and the next one. But if we don't have at least five, we will
900 probably collapse the curve down onto a single point, which means
901 we won't be able to fit it with a spline. */
902 if (CURVE_LENGTH (curve) < 5)
903 {
904 /* LOG1 ("Length is %u, not enough to filter.\n", CURVE_LENGTH (curve)); */
905 return;
906 }
907
908 for (iteration = 0; iteration < filter_iteration_count; iteration++)
909 {
910 curve_type new_curve = copy_most_of_curve (curve);
911
912 /* Keep the first point on the curve. */
913 if (offset)
914 append_point (new_curve, CURVE_POINT (curve, 0));
915
916 for (this_point = offset; this_point < CURVE_LENGTH (curve) - offset;
917 this_point++)
918 {
919 real angle, angle_alt;
920 vector_type in, in_alt, out, out_alt, sum;
921 unsigned count, alt_count;
922 real_coordinate_type new_point;
923
924 /* Find the angle using the usual number of surrounding points
925 on the curve. */
926 find_curve_vectors (this_point, curve, filter_surround,
927 &in, &out, &count);
928 angle = filter_angle (in, out);
929
930 /* Find the angle using the alternative (presumably smaller)
931 number. */
932 find_curve_vectors (this_point, curve, filter_alternative_surround,
933 &in_alt, &out_alt, &alt_count);
934 angle_alt = filter_angle (in_alt, out_alt);
935
936 /* If the alternate angle is enough larger than the usual one
937 and neither of the components of the sum are zero, use it.
938 (We don't use absolute value here, since if the alternate
939 angle is smaller, we don't particularly care, since that
940 means the curve is pretty flat around the current point,
941 anyway.) This helps keep small features from disappearing
942 into the surrounding curve. */
943 sum = Vadd (in_alt, out_alt);
944 if (angle_alt - angle >= filter_epsilon
945 && sum.dx != 0 && sum.dy != 0)
946 FILTER_ASSIGN (alt);
947
948 #if 0
949 /* This code isn't needed anymore, since we do the filtering in a
950 somewhat more general way. */
951 /* If we're left with an angle of zero, don't stop yet; we
952 might be at a straight which really isn't one (as in the `o'
953 discussed above). */
954 if (epsilon_equal (angle, 0.0))
955 {
956 real angle_secondary;
957 vector_type in_secondary, out_secondary;
958 unsigned in_secondary_count, out_secondary_count;
959
960 find_curve_vectors (this_point, curve, filter_secondary_surround,
961 &in_secondary, &out_secondary,
962 &in_secondary_count, &out_secondary_count);
963 angle_secondary = filter_angle (in_secondary, out_secondary);
964 if (!epsilon_equal (angle_secondary, 0.0))
965 FILTER_ASSIGN (secondary);
966 }
967 #endif /* 0 */
968
969 /* Start with the old point. */
970 new_point = CURVE_POINT (curve, this_point);
971 sum = Vadd (in, out);
972 new_point.x += sum.dx * filter_percent / count;
973 new_point.y += sum.dy * filter_percent / count;
974
975 /* Put the newly computed point into a separate curve, so it
976 doesn't affect future computation (on this iteration). */
977 append_point (new_curve, new_point);
978 }
979
980 /* Just as with the first point, we have to keep the last point. */
981 if (offset)
982 append_point (new_curve, LAST_CURVE_POINT (curve));
983
984 /* Set the original curve to the newly filtered one, and go again. */
985 free_curve (curve);
986 *curve = *new_curve;
987 }
988
989 /* log_curve (curve, false); */
990 /* display_curve (curve); */
991 }
992
993
994 /* Return the vectors IN and OUT, computed by looking at SURROUND points
995 on either side of TEST_INDEX. Also return the number of points in
996 the vectors in COUNT (we make sure they are the same). */
997
998 static void
find_curve_vectors(unsigned test_index,curve_type curve,unsigned surround,vector_type * in,vector_type * out,unsigned * count)999 find_curve_vectors (unsigned test_index, curve_type curve,
1000 unsigned surround,
1001 vector_type *in, vector_type *out, unsigned *count)
1002 {
1003 int i;
1004 unsigned in_count, out_count;
1005 unsigned n_done;
1006 real_coordinate_type candidate = CURVE_POINT (curve, test_index);
1007
1008 /* Add up the differences from p of the `surround' points
1009 before p. */
1010
1011 in->dx = 0.0;
1012 in->dy = 0.0;
1013
1014 for (i = CURVE_PREV (curve, test_index), n_done = 0;
1015 i >= 0 && n_done < surround; /* Do not wrap around. */
1016 i = CURVE_PREV (curve, i), n_done++)
1017 *in = Vadd (*in, Psubtract (CURVE_POINT (curve, i), candidate));
1018 in_count = n_done;
1019
1020 /* And the points after p. Don't use more points after p than we
1021 ended up with before it. */
1022 out->dx = 0.0;
1023 out->dy = 0.0;
1024
1025 for (i = CURVE_NEXT (curve, test_index), n_done = 0;
1026 i < CURVE_LENGTH (curve) && n_done < surround && n_done < in_count;
1027 i = CURVE_NEXT (curve, i), n_done++)
1028 *out = Vadd (*out, Psubtract (CURVE_POINT (curve, i), candidate));
1029 out_count = n_done;
1030
1031 /* If we used more points before p than after p, we have to go back
1032 and redo it. (We could just subtract the ones that were missing,
1033 but for this few of points, efficiency doesn't matter.) */
1034 if (out_count < in_count)
1035 {
1036 in->dx = 0.0;
1037 in->dy = 0.0;
1038
1039 for (i = CURVE_PREV (curve, test_index), n_done = 0;
1040 i >= 0 && n_done < out_count;
1041 i = CURVE_PREV (curve, i), n_done++)
1042 *in = Vadd (*in, Psubtract (CURVE_POINT (curve, i), candidate));
1043 in_count = n_done;
1044 }
1045
1046 assert (in_count == out_count);
1047 *count = in_count;
1048 }
1049
1050
1051 /* Find the angle between the vectors IN and OUT, and bring it into the
1052 range [0,45]. */
1053
1054 static real
filter_angle(vector_type in,vector_type out)1055 filter_angle (vector_type in, vector_type out)
1056 {
1057 real angle = Vangle (in, out);
1058
1059 /* What we want to do between 90 and 180 is the same as what we
1060 want to do between 0 and 90. */
1061 angle = fmod (angle, 1990.0);
1062
1063 /* And what we want to do between 45 and 90 is the same as
1064 between 0 and 45, only reversed. */
1065 if (angle > 45.0) angle = 90.0 - angle;
1066
1067 return angle;
1068 }
1069
1070 /* This routine returns the curve fitted to a straight line in a very
1071 simple way: make the first and last points on the curve be the
1072 endpoints of the line. This simplicity is justified because we are
1073 called only on very short curves. */
1074
1075 static spline_list_type *
fit_with_line(curve_type curve)1076 fit_with_line (curve_type curve)
1077 {
1078 spline_type line = new_spline ();
1079
1080 /* LOG ("Fitting with straight line:\n"); */
1081 /* REPORT ("l"); */
1082
1083 SPLINE_DEGREE (line) = LINEAR;
1084 START_POINT (line) = CONTROL1 (line) = CURVE_POINT (curve, 0);
1085 END_POINT (line) = CONTROL2 (line) = LAST_CURVE_POINT (curve);
1086
1087 /* Make sure that this line is never changed to a cubic. */
1088 SPLINE_LINEARITY (line) = 0;
1089
1090 /* if (logging) */
1091 /* { */
1092 /* LOG (" "); */
1093 /* print_spline (log_file, line); */
1094 /* } */
1095
1096 return init_spline_list (line);
1097 }
1098
1099 /* The least squares method is well described in Schneider's thesis.
1100 Briefly, we try to fit the entire curve with one spline. If that fails,
1101 we try reparameterizing, i.e., finding new, and supposedly better,
1102 t values. If that still fails, we subdivide the curve. */
1103
1104 static spline_list_type *
fit_with_least_squares(curve_type curve)1105 fit_with_least_squares (curve_type curve)
1106 {
1107 real error, best_error = FLT_MAX;
1108 spline_type spline, best_spline;
1109 spline_list_type *spline_list;
1110 unsigned worst_point;
1111 unsigned iteration = 0;
1112 real previous_error = FLT_MAX;
1113 real improvement = FLT_MAX;
1114
1115 /* FIXME: Initialize best_spline to zeroes. This is strictly not
1116 necessary as best_spline is always set in the loop below. But the
1117 compiler thinks it isn't and warns. Ideally, the code should be
1118 rewritten such that best_spline and best_error are initialized with
1119 the first values before the loop begins. */
1120 memset (&best_spline, 0, sizeof best_spline);
1121
1122 /* LOG ("\nFitting with least squares:\n"); */
1123
1124 /* Phoenix reduces the number of points with a ``linear spline
1125 technique''. But for fitting letterforms, that is
1126 inappropriate. We want all the points we can get. */
1127
1128 /* It makes no difference whether we first set the `t' values or
1129 find the tangents. This order makes the documentation a little
1130 more coherent. */
1131
1132 /* LOG ("Finding tangents:\n"); */
1133 find_tangent (curve, /* to_start */ true, /* cross_curve */ false);
1134 find_tangent (curve, /* to_start */ false, /* cross_curve */ false);
1135
1136 set_initial_parameter_values (curve);
1137
1138 /* Now we loop, reparameterizing and/or subdividing, until CURVE has
1139 been fit. */
1140 while (true)
1141 {
1142 /* LOG (" fitted to spline:\n"); */
1143
1144 spline = fit_one_spline (curve);
1145
1146 /* if (logging) */
1147 /* { */
1148 /* LOG (" "); */
1149 /* print_spline (log_file, spline); */
1150 /* } */
1151
1152 error = find_error (curve, spline, &worst_point);
1153 if (error > previous_error)
1154 {
1155 /* LOG ("Reparameterization made it worse.\n"); */
1156 /* Just fall through; we will subdivide. */
1157 }
1158 else
1159 {
1160 best_error = error;
1161 best_spline = spline;
1162 }
1163
1164 improvement = 1.0 - error / previous_error;
1165
1166 /* Don't exit, even if the error is less than `error_threshold',
1167 since we might be able to improve the fit with further
1168 reparameterization. If the reparameterization made it worse,
1169 we will exit here, since `improvement' will be negative. */
1170 if (improvement < reparameterize_improvement
1171 || error > reparameterize_threshold)
1172 break;
1173
1174 iteration++;
1175 /* LOG3 ("Reparameterization #%u (error was %.3f, a %u%% improvement):\n", */
1176 /* iteration, error, ((unsigned) (improvement * 100.0))); */
1177
1178 /* The reparameterization might fail, if the initial fit was
1179 better than `reparameterize_threshold', yet worse than the
1180 Newton-Raphson algorithm could handle. */
1181 if (!reparameterize (curve, spline))
1182 break;
1183
1184 previous_error = error;
1185 }
1186
1187 /* Go back to the best fit. */
1188 spline = best_spline;
1189 error = best_error;
1190
1191 if (error < error_threshold)
1192 {
1193 /* The points were fitted with a (possibly reparameterized)
1194 spline. We end up here whenever a fit is accepted. We have
1195 one more job: see if the ``curve'' that was fit should really
1196 be a straight line. */
1197 if (spline_linear_enough (&spline, curve))
1198 {
1199 SPLINE_DEGREE (spline) = LINEAR;
1200 /* LOG ("Changed to line.\n"); */
1201 }
1202 spline_list = init_spline_list (spline);
1203 /* LOG1 ("Accepted error of %.3f.\n", error); */
1204 }
1205 else
1206 {
1207 /* We couldn't fit the curve acceptably, so subdivide. */
1208 unsigned subdivision_index;
1209 spline_list_type *left_spline_list;
1210 spline_list_type *right_spline_list;
1211 curve_type left_curve = new_curve ();
1212 curve_type right_curve = new_curve ();
1213
1214 /* Keep the linked list of curves intact. */
1215 NEXT_CURVE (right_curve) = NEXT_CURVE (curve);
1216 PREVIOUS_CURVE (right_curve) = left_curve;
1217 NEXT_CURVE (left_curve) = right_curve;
1218 PREVIOUS_CURVE (left_curve) = curve;
1219 NEXT_CURVE (curve) = left_curve;
1220
1221 /* REPORT ("s"); */
1222 /* LOG1 ("\nSubdividing (error %.3f):\n", error); */
1223 /* LOG3 (" Original point: (%.3f,%.3f), #%u.\n", */
1224 /* CURVE_POINT (curve, worst_point).x, */
1225 /* CURVE_POINT (curve, worst_point).y, worst_point); */
1226 subdivision_index = find_subdivision (curve, worst_point);
1227 /* LOG3 (" Final point: (%.3f,%.3f), #%u.\n", */
1228 /* CURVE_POINT (curve, subdivision_index).x, */
1229 /* CURVE_POINT (curve, subdivision_index).y, subdivision_index); */
1230 /* display_subdivision (CURVE_POINT (curve, subdivision_index)); */
1231
1232 /* The last point of the left-hand curve will also be the first
1233 point of the right-hand curve. */
1234 CURVE_LENGTH (left_curve) = subdivision_index + 1;
1235 CURVE_LENGTH (right_curve) = CURVE_LENGTH (curve) - subdivision_index;
1236 left_curve->point_list = curve->point_list;
1237 right_curve->point_list = curve->point_list + subdivision_index;
1238
1239 /* We want to use the tangents of the curve which we are
1240 subdividing for the start tangent for left_curve and the
1241 end tangent for right_curve. */
1242 CURVE_START_TANGENT (left_curve) = CURVE_START_TANGENT (curve);
1243 CURVE_END_TANGENT (right_curve) = CURVE_END_TANGENT (curve);
1244
1245 /* We have to set up the two curves before finding the tangent at
1246 the subdivision point. The tangent at that point must be the
1247 same for both curves, or noticeable bumps will occur in the
1248 character. But we want to use information on both sides of the
1249 point to compute the tangent, hence cross_curve = true. */
1250 find_tangent (left_curve, /* to_start_point: */ false,
1251 /* cross_curve: */ true);
1252 CURVE_START_TANGENT (right_curve) = CURVE_END_TANGENT (left_curve);
1253
1254 /* Now that we've set up the curves, we can fit them. */
1255 left_spline_list = fit_curve (left_curve);
1256 right_spline_list = fit_curve (right_curve);
1257
1258 /* Neither of the subdivided curves could be fit, so fail. */
1259 if (left_spline_list == NULL && right_spline_list == NULL)
1260 return NULL;
1261
1262 /* Put the two together (or whichever of them exist). */
1263 spline_list = new_spline_list ();
1264
1265 if (left_spline_list == NULL)
1266 {
1267 WARNING ("could not fit left spline list");
1268 /* LOG1 ("Could not fit spline to left curve (%x).\n", */
1269 /* (unsigned) left_curve); */
1270 }
1271 else
1272 concat_spline_lists (spline_list, *left_spline_list);
1273
1274 if (right_spline_list == NULL)
1275 {
1276 WARNING ("could not fit right spline list");
1277 /* LOG1 ("Could not fit spline to right curve (%x).\n", */
1278 /* (unsigned) right_curve); */
1279 }
1280 else
1281 concat_spline_lists (spline_list, *right_spline_list);
1282 }
1283
1284 return spline_list;
1285 }
1286
1287
1288 /* Our job here is to find alpha1 (and alpha2), where t1_hat (t2_hat) is
1289 the tangent to CURVE at the starting (ending) point, such that:
1290
1291 control1 = alpha1*t1_hat + starting point
1292 control2 = alpha2*t2_hat + ending_point
1293
1294 and the resulting spline (starting_point .. control1 and control2 ..
1295 ending_point) minimizes the least-square error from CURVE.
1296
1297 See pp.57--59 of the Phoenix thesis.
1298
1299 The B?(t) here corresponds to B_i^3(U_i) there.
1300 The Bernshte\u in polynomials of degree n are defined by
1301 B_i^n(t) = { n \choose i } t^i (1-t)^{n-i}, i = 0..n */
1302
1303 #define B0(t) CUBE (1 - (t))
1304 #define B1(t) (3.0 * (t) * SQUARE (1 - (t)))
1305 #define B2(t) (3.0 * SQUARE (t) * (1 - (t)))
1306 #define B3(t) CUBE (t)
1307
1308 #define U(i) CURVE_T (curve, i)
1309
1310 static spline_type
fit_one_spline(curve_type curve)1311 fit_one_spline (curve_type curve)
1312 {
1313 /* Since our arrays are zero-based, the `C0' and `C1' here correspond
1314 to `C1' and `C2' in the paper. */
1315 real X_C1_det, C0_X_det, C0_C1_det;
1316 real alpha1, alpha2;
1317 spline_type spline;
1318 vector_type start_vector, end_vector;
1319 unsigned i;
1320 vector_type t1_hat = *CURVE_START_TANGENT (curve);
1321 vector_type t2_hat = *CURVE_END_TANGENT (curve);
1322 real C[2][2] = { { 0.0, 0.0 }, { 0.0, 0.0 } };
1323 real X[2] = { 0.0, 0.0 };
1324 int Alen = CURVE_LENGTH (curve);
1325 vector_type *A;
1326
1327 A = g_new0 (vector_type, Alen * 2);
1328
1329 START_POINT (spline) = CURVE_POINT (curve, 0);
1330 END_POINT (spline) = LAST_CURVE_POINT (curve);
1331 SPLINE_LINEARITY (spline) = 0;
1332 start_vector = make_vector (START_POINT (spline));
1333 end_vector = make_vector (END_POINT (spline));
1334
1335 for (i = 0; i < CURVE_LENGTH (curve); i++)
1336 {
1337 A[i*2+0] = Vmult_scalar (t1_hat, B1 (U (i)));
1338 A[i*2+1] = Vmult_scalar (t2_hat, B2 (U (i)));
1339 }
1340
1341 for (i = 0; i < CURVE_LENGTH (curve); i++)
1342 {
1343 vector_type temp, temp0, temp1, temp2, temp3;
1344 vector_type *Ai = &A[i*2];
1345
1346 C[0][0] += Vdot (Ai[0], Ai[0]);
1347 C[0][1] += Vdot (Ai[0], Ai[1]);
1348 /* C[1][0] = C[0][1] (this is assigned outside the loop) */
1349 C[1][1] += Vdot (Ai[1], Ai[1]);
1350
1351 /* Now the right-hand side of the equation in the paper. */
1352 temp0 = Vmult_scalar (start_vector, B0 (U (i)));
1353 temp1 = Vmult_scalar (start_vector, B1 (U (i)));
1354 temp2 = Vmult_scalar (end_vector, B2 (U (i)));
1355 temp3 = Vmult_scalar (end_vector, B3 (U (i)));
1356
1357 temp = make_vector (Vsubtract_point (CURVE_POINT (curve, i),
1358 Vadd (temp0, Vadd (temp1, Vadd (temp2, temp3)))));
1359
1360 X[0] += Vdot (temp, Ai[0]);
1361 X[1] += Vdot (temp, Ai[1]);
1362 }
1363
1364 C[1][0] = C[0][1];
1365
1366 X_C1_det = X[0] * C[1][1] - X[1] * C[0][1];
1367 C0_X_det = C[0][0] * X[1] - C[0][1] * X[0];
1368 C0_C1_det = C[0][0] * C[1][1] - C[1][0] * C[0][1];
1369 if (C0_C1_det == 0.0)
1370 FATAL ("zero determinant of C0*C1");
1371
1372 alpha1 = X_C1_det / C0_C1_det;
1373 alpha2 = C0_X_det / C0_C1_det;
1374
1375 CONTROL1 (spline) = Vadd_point (START_POINT (spline),
1376 Vmult_scalar (t1_hat, alpha1));
1377 CONTROL2 (spline) = Vadd_point (END_POINT (spline),
1378 Vmult_scalar (t2_hat, alpha2));
1379 SPLINE_DEGREE (spline) = CUBIC;
1380
1381 g_free (A);
1382
1383 return spline;
1384 }
1385
1386 /* Find closer-to-optimal t values for the given x,y's and control
1387 points, using Newton-Raphson iteration. A good description of this
1388 is in Plass & Stone. This routine performs one step in the
1389 iteration, not the whole thing. */
1390
1391 static boolean
reparameterize(curve_type curve,spline_type S)1392 reparameterize (curve_type curve, spline_type S)
1393 {
1394 unsigned p, i;
1395 spline_type S1, S2; /* S' and S''. */
1396
1397 /* REPORT ("r"); */
1398
1399 /* Find the first and second derivatives of S. To make
1400 `evaluate_spline' work, we fill the beginning points (i.e., the first
1401 two for a linear spline and the first three for a quadratic one),
1402 even though this is at odds with the rest of the program. */
1403 for (i = 0; i < 3; i++)
1404 {
1405 S1.v[i].x = 3.0 * (S.v[i + 1].x - S.v[i].x);
1406 S1.v[i].y = 3.0 * (S.v[i + 1].y - S.v[i].y);
1407 }
1408 S1.v[i].x = S1.v[i].y = -1.0; /* These will never be accessed. */
1409 SPLINE_DEGREE (S1) = QUADRATIC;
1410
1411 for (i = 0; i < 2; i++)
1412 {
1413 S2.v[i].x = 2.0 * (S1.v[i + 1].x - S1.v[i].x);
1414 S2.v[i].y = 2.0 * (S1.v[i + 1].y - S1.v[i].y);
1415 }
1416 S2.v[2].x = S2.v[2].y = S2.v[3].x = S2.v[3].y = -1.0;
1417 SPLINE_DEGREE (S2) = LINEAR;
1418
1419 for (p = 0; p < CURVE_LENGTH (curve); p++)
1420 {
1421 real new_distance, old_distance;
1422 real_coordinate_type new_point;
1423 real_coordinate_type point = CURVE_POINT (curve, p);
1424 real t = CURVE_T (curve, p);
1425
1426 /* Find the points at this t on S, S', and S''. */
1427 real_coordinate_type S_t = evaluate_spline (S, t);
1428 real_coordinate_type S1_t = evaluate_spline (S1, t);
1429 real_coordinate_type S2_t = evaluate_spline (S2, t);
1430
1431 /* The differences in x and y (Q1(t) on p.62 of Schneider's thesis). */
1432 real_coordinate_type d;
1433 real numerator;
1434 real denominator;
1435
1436 d.x = S_t.x - point.x;
1437 d.y = S_t.y - point.y;
1438
1439 /* The step size, f(t)/f'(t). */
1440 numerator = d.x * S1_t.x + d.y * S1_t.y;
1441 denominator = (SQUARE (S1_t.x) + d.x * S2_t.x
1442 + SQUARE (S1_t.y) + d.y * S2_t.y);
1443
1444 /* We compute the distances only to be able to check that we
1445 really are moving closer. I don't know how this condition can
1446 be reliably checked for in advance, but I know what it means in
1447 practice: the original fit was not good enough for the
1448 Newton-Raphson iteration to converge. Therefore, we need to
1449 abort the reparameterization, and subdivide. */
1450 old_distance = distance (S_t, point);
1451 CURVE_T (curve, p) -= numerator / denominator;
1452 new_point = evaluate_spline (S, CURVE_T (curve, p));
1453 new_distance = distance (new_point, point);
1454
1455 if (new_distance > old_distance)
1456 {
1457 /* REPORT ("!"); */
1458 /* LOG3 (" Stopped reparameterizing; %.3f > %.3f at point %u.\n", */
1459 /* new_distance, old_distance, p); */
1460 return false;
1461 }
1462
1463 /* The t value might be negative or > 1, if the choice of control
1464 points wasn't so great. We have no difficulty in evaluating
1465 a spline at a t outside the range zero to one, though, so it
1466 doesn't matter. (Although it is a little unconventional.) */
1467 }
1468 /* LOG (" reparameterized curve:\n "); */
1469 /* log_curve (curve, true); */
1470
1471 return true;
1472 }
1473
1474 /* This routine finds the best place to subdivide the curve CURVE,
1475 somewhere near to the point whose index is INITIAL. Originally,
1476 curves were always subdivided at the point of worst error, which is
1477 intuitively appealing, but doesn't always give the best results. For
1478 example, at the end of a serif that tapers into the stem, the best
1479 subdivision point is at the point where they join, even if the worst
1480 point is a little ways into the serif.
1481
1482 We return the index of the point at which to subdivide. */
1483
1484 static unsigned
find_subdivision(curve_type curve,unsigned initial)1485 find_subdivision (curve_type curve, unsigned initial)
1486 {
1487 unsigned i, n_done;
1488 int best_point = -1;
1489 vector_type best = { FLT_MAX, FLT_MAX };
1490 unsigned search = subdivide_search * CURVE_LENGTH (curve);
1491
1492 /* LOG1 (" Number of points to search: %u.\n", search); */
1493
1494 /* We don't want to find the first (or last) point in the curve. */
1495 for (i = initial, n_done = 0; i > 0 && n_done < search;
1496 i = CURVE_PREV (curve, i), n_done++)
1497 {
1498 if (test_subdivision_point (curve, i, &best))
1499 {
1500 best_point = i;
1501 /* LOG3 (" Better point: (%.3f,%.3f), #%u.\n", */
1502 /* CURVE_POINT (curve, i).x, CURVE_POINT (curve, i).y, i); */
1503 }
1504 }
1505
1506 /* If we found a good one, let's take it. */
1507 if (best_point != -1)
1508 return best_point;
1509
1510 for (i = CURVE_NEXT (curve, initial), n_done = 0;
1511 i < CURVE_LENGTH (curve) - 1 && n_done < search;
1512 i = CURVE_NEXT (curve, i), n_done++)
1513 {
1514 if (test_subdivision_point (curve, i, &best))
1515 {
1516 best_point = i;
1517 /* LOG3 (" Better point at (%.3f,%.3f), #%u.\n", */
1518 /* CURVE_POINT (curve, i).x, CURVE_POINT (curve, i).y, i); */
1519 }
1520 }
1521
1522 /* If we didn't find any better point, return the original. */
1523 return best_point == -1 ? initial : best_point;
1524 }
1525
1526
1527 /* Here are some macros that decide whether or not we're at a
1528 ``join point'', as described above. */
1529 #define ONLY_ONE_LESS(v) \
1530 (((v).dx < subdivide_threshold && (v).dy > subdivide_threshold) \
1531 || ((v).dy < subdivide_threshold && (v).dx > subdivide_threshold))
1532
1533 #define BOTH_GREATER(v) \
1534 ((v).dx > subdivide_threshold && (v).dy > subdivide_threshold)
1535
1536 /* We assume that the vectors V1 and V2 are nonnegative. */
1537 #define JOIN(v1, v2) \
1538 ((ONLY_ONE_LESS (v1) && BOTH_GREATER (v2)) \
1539 || (ONLY_ONE_LESS (v2) && BOTH_GREATER (v1)))
1540
1541 /* If the component D of the vector V is smaller than the best so far,
1542 update the best point. */
1543 #define UPDATE_BEST(v, d) \
1544 do \
1545 { \
1546 if ((v).d < subdivide_threshold && (v).d < best->d) \
1547 best->d = (v).d; \
1548 } \
1549 while (0)
1550
1551
1552 /* If the point INDEX in the curve CURVE is the best subdivision point
1553 we've found so far, update the vector BEST. */
1554
1555 static boolean
test_subdivision_point(curve_type curve,unsigned index,vector_type * best)1556 test_subdivision_point (curve_type curve, unsigned index, vector_type *best)
1557 {
1558 unsigned count;
1559 vector_type in, out;
1560 boolean join = false;
1561
1562 find_curve_vectors (index, curve, subdivide_surround, &in, &out, &count);
1563
1564 /* We don't want to subdivide at points which are very close to the
1565 endpoints, so if the vectors aren't computed from as many points as
1566 we asked for, don't bother checking this point. */
1567 if (count == subdivide_surround)
1568 {
1569 in = Vabs (in);
1570 out = Vabs (out);
1571
1572 join = JOIN (in, out);
1573 if (join)
1574 {
1575 UPDATE_BEST (in, dx);
1576 UPDATE_BEST (in, dy);
1577 UPDATE_BEST (out, dx);
1578 UPDATE_BEST (out, dy);
1579 }
1580 }
1581
1582 return join;
1583 }
1584
1585 /* Find reasonable values for t for each point on CURVE. The method is
1586 called chord-length parameterization, which is described in Plass &
1587 Stone. The basic idea is just to use the distance from one point to
1588 the next as the t value, normalized to produce values that increase
1589 from zero for the first point to one for the last point. */
1590
1591 static void
set_initial_parameter_values(curve_type curve)1592 set_initial_parameter_values (curve_type curve)
1593 {
1594 unsigned p;
1595
1596 /* LOG ("\nAssigning initial t values:\n "); */
1597
1598 CURVE_T (curve, 0) = 0.0;
1599
1600 for (p = 1; p < CURVE_LENGTH (curve); p++)
1601 {
1602 real_coordinate_type point = CURVE_POINT (curve, p),
1603 previous_p = CURVE_POINT (curve, p - 1);
1604 real d = distance (point, previous_p);
1605 CURVE_T (curve, p) = CURVE_T (curve, p - 1) + d;
1606 }
1607
1608 assert (LAST_CURVE_T (curve) != 0.0);
1609
1610 for (p = 1; p < CURVE_LENGTH (curve); p++)
1611 CURVE_T (curve, p) = CURVE_T (curve, p) / LAST_CURVE_T (curve);
1612
1613 /* log_entire_curve (curve); */
1614 }
1615
1616 /* Find an approximation to the tangent to an endpoint of CURVE (to the
1617 first point if TO_START_POINT is true, else the last). If
1618 CROSS_CURVE is true, consider points on the adjacent curve to CURVE.
1619
1620 It is important to compute an accurate approximation, because the
1621 control points that we eventually decide upon to fit the curve will
1622 be placed on the half-lines defined by the tangents and
1623 endpoints...and we never recompute the tangent after this. */
1624
1625 static void
find_tangent(curve_type curve,boolean to_start_point,boolean cross_curve)1626 find_tangent (curve_type curve, boolean to_start_point, boolean cross_curve)
1627 {
1628 vector_type tangent;
1629 vector_type **curve_tangent = to_start_point ? &(CURVE_START_TANGENT (curve))
1630 : &(CURVE_END_TANGENT (curve));
1631 unsigned n_points = 0;
1632
1633 /* LOG1 (" tangent to %s: ", to_start_point ? "start" : "end"); */
1634
1635 if (*curve_tangent == NULL)
1636 {
1637 *curve_tangent = g_new (vector_type, 1);
1638 tangent = find_half_tangent (curve, to_start_point, &n_points);
1639
1640 if (cross_curve)
1641 {
1642 curve_type adjacent_curve
1643 = to_start_point ? PREVIOUS_CURVE (curve) : NEXT_CURVE (curve);
1644 vector_type tangent2
1645 = find_half_tangent (adjacent_curve, !to_start_point, &n_points);
1646
1647 /* LOG2 ("(adjacent curve half tangent (%.3f,%.3f)) ", */
1648 /* tangent2.dx, tangent2.dy); */
1649 tangent = Vadd (tangent, tangent2);
1650 }
1651
1652 assert (n_points > 0);
1653 **curve_tangent = Vmult_scalar (tangent, 1.0 / n_points);
1654 }
1655 else
1656 {
1657 /* LOG ("(already computed) "); */
1658 }
1659
1660 /* LOG2 ("(%.3f,%.3f).\n", (*curve_tangent)->dx, (*curve_tangent)->dy); */
1661 }
1662
1663
1664 /* Find the change in y and change in x for `tangent_surround' (a global)
1665 points along CURVE. Increment N_POINTS by the number of points we
1666 actually look at. */
1667
1668 static vector_type
find_half_tangent(curve_type c,boolean to_start_point,unsigned * n_points)1669 find_half_tangent (curve_type c, boolean to_start_point, unsigned *n_points)
1670 {
1671 unsigned p;
1672 int factor = to_start_point ? 1 : -1;
1673 unsigned tangent_index = to_start_point ? 0 : c->length - 1;
1674 real_coordinate_type tangent_point = CURVE_POINT (c, tangent_index);
1675 vector_type tangent;
1676
1677 tangent.dx = 0.0;
1678 tangent.dy = 0.0;
1679
1680 for (p = 1; p <= tangent_surround; p++)
1681 {
1682 int this_index = p * factor + tangent_index;
1683 real_coordinate_type this_point;
1684
1685 if (this_index < 0 || this_index >= c->length)
1686 break;
1687
1688 this_point = CURVE_POINT (c, p * factor + tangent_index);
1689
1690 /* Perhaps we should weight the tangent from `this_point' by some
1691 factor dependent on the distance from the tangent point. */
1692 tangent = Vadd (tangent,
1693 Vmult_scalar (Psubtract (this_point, tangent_point),
1694 factor));
1695 (*n_points)++;
1696 }
1697
1698 return tangent;
1699 }
1700
1701 /* When this routine is called, we have computed a spline representation
1702 for the digitized curve. The question is, how good is it? If the
1703 fit is very good indeed, we might have an error of zero on each
1704 point, and then WORST_POINT becomes irrelevant. But normally, we
1705 return the error at the worst point, and the index of that point in
1706 WORST_POINT. The error computation itself is the Euclidean distance
1707 from the original curve CURVE to the fitted spline SPLINE. */
1708
1709 static real
find_error(curve_type curve,spline_type spline,unsigned * worst_point)1710 find_error (curve_type curve, spline_type spline, unsigned *worst_point)
1711 {
1712 unsigned this_point;
1713 real total_error = 0.0;
1714 real worst_error = FLT_MIN;
1715
1716 *worst_point = CURVE_LENGTH (curve) + 1; /* A sentinel value. */
1717
1718 for (this_point = 0; this_point < CURVE_LENGTH (curve); this_point++)
1719 {
1720 real_coordinate_type curve_point = CURVE_POINT (curve, this_point);
1721 real t = CURVE_T (curve, this_point);
1722 real_coordinate_type spline_point = evaluate_spline (spline, t);
1723 real this_error = distance (curve_point, spline_point);
1724
1725 if (this_error > worst_error)
1726 {
1727 *worst_point = this_point;
1728 worst_error = this_error;
1729 }
1730 total_error += this_error;
1731 }
1732
1733 if (*worst_point == CURVE_LENGTH (curve) + 1)
1734 { /* Didn't have any ``worst point''; the error should be zero. */
1735 if (epsilon_equal (total_error, 0.0))
1736 {
1737 /* LOG (" Every point fit perfectly.\n"); */
1738 }
1739 else
1740 printf ("No worst point found; something is wrong");
1741 }
1742 else
1743 {
1744 /* LOG4 (" Worst error (at (%.3f,%.3f), point #%u) was %.3f.\n", */
1745 /* CURVE_POINT (curve, *worst_point).x, */
1746 /* CURVE_POINT (curve, *worst_point).y, *worst_point, worst_error); */
1747 /* LOG1 (" Total error was %.3f.\n", total_error); */
1748 /* LOG2 (" Average error (over %u points) was %.3f.\n", */
1749 /* CURVE_LENGTH (curve), total_error / CURVE_LENGTH (curve)); */
1750 }
1751
1752 return worst_error;
1753 }
1754
1755 /* Supposing that we have accepted the error, another question arises:
1756 would we be better off just using a straight line? */
1757
1758 static boolean
spline_linear_enough(spline_type * spline,curve_type curve)1759 spline_linear_enough (spline_type *spline, curve_type curve)
1760 {
1761 real A, B, C, slope;
1762 unsigned this_point;
1763 real distance = 0.0;
1764
1765 /* LOG ("Checking linearity:\n"); */
1766
1767 /* For a line described by Ax + By + C = 0, the distance d from a
1768 point (x0,y0) to that line is:
1769
1770 d = | Ax0 + By0 + C | / sqrt (A^2 + B^2).
1771
1772 We can find A, B, and C from the starting and ending points,
1773 unless the line defined by those two points is vertical. */
1774
1775 if (epsilon_equal (START_POINT (*spline).x, END_POINT (*spline).x))
1776 {
1777 A = 1;
1778 B = 0;
1779 C = -START_POINT (*spline).x;
1780 }
1781 else
1782 {
1783 /* Plug the spline's starting and ending points into the two-point
1784 equation for a line, that is,
1785
1786 (y - y1) = ((y2 - y1)/(x2 - x1)) * (x - x1)
1787
1788 to get the values for A, B, and C. */
1789
1790 slope = ((END_POINT (*spline).y - START_POINT (*spline).y)
1791 / (END_POINT (*spline).x - START_POINT (*spline).x));
1792 A = -slope;
1793 B = 1;
1794 C = slope * START_POINT (*spline).x - START_POINT (*spline).y;
1795 }
1796 /* LOG3 (" Line is %.3fx + %.3fy + %.3f = 0.\n", A, B, C); */
1797
1798 for (this_point = 0; this_point < CURVE_LENGTH (curve); this_point++)
1799 {
1800 real t = CURVE_T (curve, this_point);
1801 real_coordinate_type spline_point = evaluate_spline (*spline, t);
1802
1803 distance += fabs (A * spline_point.x + B * spline_point.y + C)
1804 / sqrt (A * A + B * B);
1805 }
1806 /* LOG1 (" Total distance is %.3f, ", distance); */
1807
1808 distance /= CURVE_LENGTH (curve);
1809 /* LOG1 ("which is %.3f normalized.\n", distance); */
1810
1811 /* We want reversion of short curves to splines to be more likely than
1812 reversion of long curves, hence the second division by the curve
1813 length, for use in `change_bad_lines'. */
1814 SPLINE_LINEARITY (*spline) = distance / CURVE_LENGTH (curve);
1815 /* LOG1 (" Final linearity: %.3f.\n", SPLINE_LINEARITY (*spline)); */
1816
1817 return distance < line_threshold;
1818 }
1819
1820
1821 /* Unfortunately, we cannot tell in isolation whether a given spline
1822 should be changed to a line or not. That can only be known after the
1823 entire curve has been fit to a list of splines. (The curve is the
1824 pixel outline between two corners.) After subdividing the curve, a
1825 line may very well fit a portion of the curve just as well as the
1826 spline---but unless a spline is truly close to being a line, it
1827 should not be combined with other lines. */
1828
1829 static void
change_bad_lines(spline_list_type * spline_list)1830 change_bad_lines (spline_list_type *spline_list)
1831 {
1832 unsigned this_spline;
1833 boolean found_cubic = false;
1834 unsigned length = SPLINE_LIST_LENGTH (*spline_list);
1835
1836 /* LOG1 ("\nChecking for bad lines (length %u):\n", length); */
1837
1838 /* First see if there are any splines in the fitted shape. */
1839 for (this_spline = 0; this_spline < length; this_spline++)
1840 {
1841 if (SPLINE_DEGREE (SPLINE_LIST_ELT (*spline_list, this_spline)) == CUBIC)
1842 {
1843 found_cubic = true;
1844 break;
1845 }
1846 }
1847
1848 /* If so, change lines back to splines (we haven't done anything to
1849 their control points, so we only have to change the degree) unless
1850 the spline is close enough to being a line. */
1851 if (found_cubic)
1852 for (this_spline = 0; this_spline < length; this_spline++)
1853 {
1854 spline_type s = SPLINE_LIST_ELT (*spline_list, this_spline);
1855
1856 if (SPLINE_DEGREE (s) == LINEAR)
1857 {
1858 /* LOG1 (" #%u: ", this_spline); */
1859 if (SPLINE_LINEARITY (s) > line_reversion_threshold)
1860 {
1861 /* LOG ("reverted, "); */
1862 SPLINE_DEGREE (SPLINE_LIST_ELT (*spline_list, this_spline))
1863 = CUBIC;
1864 }
1865 /* LOG1 ("linearity %.3f.\n", SPLINE_LINEARITY (s)); */
1866 }
1867 }
1868 else
1869 {
1870 /* LOG (" No lines.\n"); */
1871 }
1872 }
1873
1874 /* When we have finished fitting an entire pixel outline to a spline
1875 list L, we go through L to ensure that any endpoints that are ``close
1876 enough'' (i.e., within `align_threshold') to being the same really
1877 are the same. */
1878
1879 /* This macro adjusts the AXIS axis on the starting and ending points on
1880 a particular spline if necessary. */
1881 #define TRY_AXIS(axis) \
1882 do \
1883 { \
1884 real delta = fabs (end.axis - start.axis); \
1885 \
1886 if (!epsilon_equal (delta, 0.0) && delta <= align_threshold) \
1887 { \
1888 spline_type *next = &NEXT_SPLINE_LIST_ELT (*l, this_spline); \
1889 spline_type *prev = &PREV_SPLINE_LIST_ELT (*l, this_spline); \
1890 \
1891 START_POINT (*s).axis = END_POINT (*s).axis \
1892 = END_POINT (*prev).axis = START_POINT (*next).axis \
1893 = (start.axis + end.axis) / 2; \
1894 spline_change = true; \
1895 } \
1896 } \
1897 while (0)
1898
1899 static void
align(spline_list_type * l)1900 align (spline_list_type *l)
1901 {
1902 boolean change;
1903 unsigned this_spline;
1904 unsigned length = SPLINE_LIST_LENGTH (*l);
1905
1906 /* LOG1 ("\nAligning spline list (length %u):\n", length); */
1907
1908 do
1909 {
1910 change = false;
1911
1912 /* LOG (" "); */
1913
1914 for (this_spline = 0; this_spline < length; this_spline++)
1915 {
1916 boolean spline_change = false;
1917 spline_type *s = &SPLINE_LIST_ELT (*l, this_spline);
1918 real_coordinate_type start = START_POINT (*s);
1919 real_coordinate_type end = END_POINT (*s);
1920
1921 TRY_AXIS (x);
1922 TRY_AXIS (y);
1923 if (spline_change)
1924 {
1925 /* LOG1 ("%u ", this_spline); */
1926 change |= spline_change;
1927 }
1928 }
1929 /* LOG ("\n"); */
1930 }
1931 while (change);
1932 }
1933
1934 /* Lists of array indices (well, that is what we use it for). */
1935
1936 static index_list_type
new_index_list(void)1937 new_index_list (void)
1938 {
1939 index_list_type index_list;
1940
1941 index_list.data = NULL;
1942 INDEX_LIST_LENGTH (index_list) = 0;
1943
1944 return index_list;
1945 }
1946
1947
1948 static void
free_index_list(index_list_type * index_list)1949 free_index_list (index_list_type *index_list)
1950 {
1951 if (INDEX_LIST_LENGTH (*index_list) > 0)
1952 {
1953 g_free (index_list->data);
1954 index_list->data = NULL;
1955 INDEX_LIST_LENGTH (*index_list) = 0;
1956 }
1957 }
1958
1959
1960 static void
append_index(index_list_type * list,unsigned new_index)1961 append_index (index_list_type *list, unsigned new_index)
1962 {
1963 INDEX_LIST_LENGTH (*list)++;
1964 list->data = (unsigned *)g_realloc(list->data,(INDEX_LIST_LENGTH (*list)) * sizeof(unsigned));
1965 /* XRETALLOC (list->data, INDEX_LIST_LENGTH (*list), unsigned); */
1966 list->data[INDEX_LIST_LENGTH (*list) - 1] = new_index;
1967 }
1968