xref: /aosp_15_r20/external/freetype/src/sdf/ftsdf.c (revision 63949dbd25bcc50c4e1178497ff9e9574d44fc5a)
1 /****************************************************************************
2  *
3  * ftsdf.c
4  *
5  *   Signed Distance Field support for outline fonts (body).
6  *
7  * Copyright (C) 2020-2023 by
8  * David Turner, Robert Wilhelm, and Werner Lemberg.
9  *
10  * Written by Anuj Verma.
11  *
12  * This file is part of the FreeType project, and may only be used,
13  * modified, and distributed under the terms of the FreeType project
14  * license, LICENSE.TXT.  By continuing to use, modify, or distribute
15  * this file you indicate that you have read the license and
16  * understand and accept it fully.
17  *
18  */
19 
20 
21 #include <freetype/internal/ftobjs.h>
22 #include <freetype/internal/ftdebug.h>
23 #include <freetype/ftoutln.h>
24 #include <freetype/fttrigon.h>
25 #include <freetype/ftbitmap.h>
26 #include "ftsdf.h"
27 
28 #include "ftsdferrs.h"
29 
30 
31   /**************************************************************************
32    *
33    * A brief technical overview of how the SDF rasterizer works
34    * ----------------------------------------------------------
35    *
36    * [Notes]:
37    *   * SDF stands for Signed Distance Field everywhere.
38    *
39    *   * This renderer generates SDF directly from outlines.  There is
40    *     another renderer called 'bsdf', which converts bitmaps to SDF; see
41    *     file `ftbsdf.c` for more.
42    *
43    *   * The basic idea of generating the SDF is taken from Viktor Chlumsky's
44    *     research paper.  The paper explains both single and multi-channel
45    *     SDF, however, this implementation only generates single-channel SDF.
46    *
47    *       Chlumsky, Viktor: Shape Decomposition for Multi-channel Distance
48    *       Fields.  Master's thesis.  Czech Technical University in Prague,
49    *       Faculty of InformationTechnology, 2015.
50    *
51    *     For more information: https://github.com/Chlumsky/msdfgen
52    *
53    * ========================================================================
54    *
55    * Generating SDF from outlines is pretty straightforward.
56    *
57    * (1) We have a set of contours that make the outline of a shape/glyph.
58    *     Each contour comprises of several edges, with three types of edges.
59    *
60    *     * line segments
61    *     * conic Bezier curves
62    *     * cubic Bezier curves
63    *
64    * (2) Apart from the outlines we also have a two-dimensional grid, namely
65    *     the bitmap that is used to represent the final SDF data.
66    *
67    * (3) In order to generate SDF, our task is to find shortest signed
68    *     distance from each grid point to the outline.  The 'signed
69    *     distance' means that if the grid point is filled by any contour
70    *     then its sign is positive, otherwise it is negative.  The pseudo
71    *     code is as follows.
72    *
73    *     ```
74    *     foreach grid_point (x, y):
75    *     {
76    *       int min_dist = INT_MAX;
77    *
78    *       foreach contour in outline:
79    *       {
80    *         foreach edge in contour:
81    *         {
82    *           // get shortest distance from point (x, y) to the edge
83    *           d = get_min_dist(x, y, edge);
84    *
85    *           if (d < min_dist)
86    *             min_dist = d;
87    *         }
88    *
89    *         bitmap[x, y] = min_dist;
90    *       }
91    *     }
92    *     ```
93    *
94    * (4) After running this algorithm the bitmap contains information about
95    *     the shortest distance from each point to the outline of the shape.
96    *     Of course, while this is the most straightforward way of generating
97    *     SDF, we use various optimizations in our implementation.  See the
98    *     `sdf_generate_*' functions in this file for all details.
99    *
100    *     The optimization currently used by default is subdivision; see
101    *     function `sdf_generate_subdivision` for more.
102    *
103    *     Also, to see how we compute the shortest distance from a point to
104    *     each type of edge, check out the `get_min_distance_*' functions.
105    *
106    */
107 
108 
109   /**************************************************************************
110    *
111    * The macro FT_COMPONENT is used in trace mode.  It is an implicit
112    * parameter of the FT_TRACE() and FT_ERROR() macros, used to print/log
113    * messages during execution.
114    */
115 #undef  FT_COMPONENT
116 #define FT_COMPONENT  sdf
117 
118 
119   /**************************************************************************
120    *
121    * definitions
122    *
123    */
124 
125   /*
126    * If set to 1, the rasterizer uses Newton-Raphson's method for finding
127    * the shortest distance from a point to a conic curve.
128    *
129    * If set to 0, an analytical method gets used instead, which computes the
130    * roots of a cubic polynomial to find the shortest distance.  However,
131    * the analytical method can currently underflow; we thus use Newton's
132    * method by default.
133    */
134 #ifndef USE_NEWTON_FOR_CONIC
135 #define USE_NEWTON_FOR_CONIC  1
136 #endif
137 
138   /*
139    * The number of intervals a Bezier curve gets sampled and checked to find
140    * the shortest distance.
141    */
142 #define MAX_NEWTON_DIVISIONS  4
143 
144   /*
145    * The number of steps of Newton's iterations in each interval of the
146    * Bezier curve.  Basically, we run Newton's approximation
147    *
148    *   x -= Q(t) / Q'(t)
149    *
150    * for each division to get the shortest distance.
151    */
152 #define MAX_NEWTON_STEPS  4
153 
154   /*
155    * The epsilon distance (in 16.16 fractional units) used for corner
156    * resolving.  If the difference of two distances is less than this value
157    * they will be checked for a corner if they are ambiguous.
158    */
159 #define CORNER_CHECK_EPSILON  32
160 
161 #if 0
162   /*
163    * Coarse grid dimension.  Will probably be removed in the future because
164    * coarse grid optimization is the slowest algorithm.
165    */
166 #define CG_DIMEN  8
167 #endif
168 
169 
170   /**************************************************************************
171    *
172    * macros
173    *
174    */
175 
176 #define MUL_26D6( a, b )  ( ( ( a ) * ( b ) ) / 64 )
177 #define VEC_26D6_DOT( p, q )  ( MUL_26D6( p.x, q.x ) + \
178                                 MUL_26D6( p.y, q.y ) )
179 
180 
181   /**************************************************************************
182    *
183    * structures and enums
184    *
185    */
186 
187   /**************************************************************************
188    *
189    * @Struct:
190    *   SDF_TRaster
191    *
192    * @Description:
193    *   This struct is used in place of @FT_Raster and is stored within the
194    *   internal FreeType renderer struct.  While rasterizing it is passed to
195    *   the @FT_Raster_RenderFunc function, which then can be used however we
196    *   want.
197    *
198    * @Fields:
199    *   memory ::
200    *     Used internally to allocate intermediate memory while raterizing.
201    *
202    */
203   typedef struct  SDF_TRaster_
204   {
205     FT_Memory  memory;
206 
207   } SDF_TRaster, *SDF_PRaster;
208 
209 
210   /**************************************************************************
211    *
212    * @Enum:
213    *   SDF_Edge_Type
214    *
215    * @Description:
216    *   Enumeration of all curve types present in fonts.
217    *
218    * @Fields:
219    *   SDF_EDGE_UNDEFINED ::
220    *     Undefined edge, simply used to initialize and detect errors.
221    *
222    *   SDF_EDGE_LINE ::
223    *     Line segment with start and end point.
224    *
225    *   SDF_EDGE_CONIC ::
226    *     A conic/quadratic Bezier curve with start, end, and one control
227    *     point.
228    *
229    *   SDF_EDGE_CUBIC ::
230    *     A cubic Bezier curve with start, end, and two control points.
231    *
232    */
233   typedef enum  SDF_Edge_Type_
234   {
235     SDF_EDGE_UNDEFINED = 0,
236     SDF_EDGE_LINE      = 1,
237     SDF_EDGE_CONIC     = 2,
238     SDF_EDGE_CUBIC     = 3
239 
240   } SDF_Edge_Type;
241 
242 
243   /**************************************************************************
244    *
245    * @Enum:
246    *   SDF_Contour_Orientation
247    *
248    * @Description:
249    *   Enumeration of all orientation values of a contour.  We determine the
250    *   orientation by calculating the area covered by a contour.  Contrary
251    *   to values returned by @FT_Outline_Get_Orientation,
252    *   `SDF_Contour_Orientation` is independent of the fill rule, which can
253    *   be different for different font formats.
254    *
255    * @Fields:
256    *   SDF_ORIENTATION_NONE ::
257    *     Undefined orientation, used for initialization and error detection.
258    *
259    *   SDF_ORIENTATION_CW ::
260    *     Clockwise orientation (positive area covered).
261    *
262    *   SDF_ORIENTATION_CCW ::
263    *     Counter-clockwise orientation (negative area covered).
264    *
265    * @Note:
266    *   See @FT_Outline_Get_Orientation for more details.
267    *
268    */
269   typedef enum  SDF_Contour_Orientation_
270   {
271     SDF_ORIENTATION_NONE = 0,
272     SDF_ORIENTATION_CW   = 1,
273     SDF_ORIENTATION_CCW  = 2
274 
275   } SDF_Contour_Orientation;
276 
277 
278   /**************************************************************************
279    *
280    * @Struct:
281    *   SDF_Edge
282    *
283    * @Description:
284    *   Represent an edge of a contour.
285    *
286    * @Fields:
287    *   start_pos ::
288    *     Start position of an edge.  Valid for all types of edges.
289    *
290    *   end_pos ::
291    *     Etart position of an edge.  Valid for all types of edges.
292    *
293    *   control_a ::
294    *     A control point of the edge.  Valid only for `SDF_EDGE_CONIC`
295    *     and `SDF_EDGE_CUBIC`.
296    *
297    *   control_b ::
298    *     Another control point of the edge.  Valid only for
299    *     `SDF_EDGE_CONIC`.
300    *
301    *   edge_type ::
302    *     Type of the edge, see @SDF_Edge_Type for all possible edge types.
303    *
304    *   next ::
305    *     Used to create a singly linked list, which can be interpreted
306    *     as a contour.
307    *
308    */
309   typedef struct  SDF_Edge_
310   {
311     FT_26D6_Vec  start_pos;
312     FT_26D6_Vec  end_pos;
313     FT_26D6_Vec  control_a;
314     FT_26D6_Vec  control_b;
315 
316     SDF_Edge_Type  edge_type;
317 
318     struct SDF_Edge_*  next;
319 
320   } SDF_Edge;
321 
322 
323   /**************************************************************************
324    *
325    * @Struct:
326    *   SDF_Contour
327    *
328    * @Description:
329    *   Represent a complete contour, which contains a list of edges.
330    *
331    * @Fields:
332    *   last_pos ::
333    *     Contains the value of `end_pos' of the last edge in the list of
334    *     edges.  Useful while decomposing the outline with
335    *     @FT_Outline_Decompose.
336    *
337    *   edges ::
338    *     Linked list of all the edges that make the contour.
339    *
340    *   next ::
341    *     Used to create a singly linked list, which can be interpreted as a
342    *     complete shape or @FT_Outline.
343    *
344    */
345   typedef struct  SDF_Contour_
346   {
347     FT_26D6_Vec  last_pos;
348     SDF_Edge*    edges;
349 
350     struct SDF_Contour_*  next;
351 
352   } SDF_Contour;
353 
354 
355   /**************************************************************************
356    *
357    * @Struct:
358    *   SDF_Shape
359    *
360    * @Description:
361    *   Represent a complete shape, which is the decomposition of
362    *   @FT_Outline.
363    *
364    * @Fields:
365    *   memory ::
366    *     Used internally to allocate memory.
367    *
368    *   contours ::
369    *     Linked list of all the contours that make the shape.
370    *
371    */
372   typedef struct  SDF_Shape_
373   {
374     FT_Memory     memory;
375     SDF_Contour*  contours;
376 
377   } SDF_Shape;
378 
379 
380   /**************************************************************************
381    *
382    * @Struct:
383    *   SDF_Signed_Distance
384    *
385    * @Description:
386    *   Represent signed distance of a point, i.e., the distance of the edge
387    *   nearest to the point.
388    *
389    * @Fields:
390    *   distance ::
391    *     Distance of the point from the nearest edge.  Can be squared or
392    *     absolute depending on the `USE_SQUARED_DISTANCES` macro defined in
393    *     file `ftsdfcommon.h`.
394    *
395    *   cross ::
396    *     Cross product of the shortest distance vector (i.e., the vector
397    *     from the point to the nearest edge) and the direction of the edge
398    *     at the nearest point.  This is used to resolve ambiguities of
399    *     `sign`.
400    *
401    *   sign ::
402    *     A value used to indicate whether the distance vector is outside or
403    *     inside the contour corresponding to the edge.
404    *
405    * @Note:
406    *   `sign` may or may not be correct, therefore it must be checked
407    *   properly in case there is an ambiguity.
408    *
409    */
410   typedef struct SDF_Signed_Distance_
411   {
412     FT_16D16  distance;
413     FT_16D16  cross;
414     FT_Char   sign;
415 
416   } SDF_Signed_Distance;
417 
418 
419   /**************************************************************************
420    *
421    * @Struct:
422    *   SDF_Params
423    *
424    * @Description:
425    *   Yet another internal parameters required by the rasterizer.
426    *
427    * @Fields:
428    *   orientation ::
429    *     This is not the @SDF_Contour_Orientation value but @FT_Orientation,
430    *     which determines whether clockwise-oriented outlines are to be
431    *     filled or counter-clockwise-oriented ones.
432    *
433    *   flip_sign ::
434    *     If set to true, flip the sign.  By default the points filled by the
435    *     outline are positive.
436    *
437    *   flip_y ::
438    *     If set to true the output bitmap is upside-down.  Can be useful
439    *     because OpenGL and DirectX use different coordinate systems for
440    *     textures.
441    *
442    *   overload_sign ::
443    *     In the subdivision and bounding box optimization, the default
444    *     outside sign is taken as -1.  This parameter can be used to modify
445    *     that behaviour.  For example, while generating SDF for a single
446    *     counter-clockwise contour, the outside sign should be 1.
447    *
448    */
449   typedef struct SDF_Params_
450   {
451     FT_Orientation  orientation;
452     FT_Bool         flip_sign;
453     FT_Bool         flip_y;
454 
455     FT_Int  overload_sign;
456 
457   } SDF_Params;
458 
459 
460   /**************************************************************************
461    *
462    * constants, initializer, and destructor
463    *
464    */
465 
466   static
467   const FT_Vector  zero_vector = { 0, 0 };
468 
469   static
470   const SDF_Edge  null_edge = { { 0, 0 }, { 0, 0 },
471                                 { 0, 0 }, { 0, 0 },
472                                 SDF_EDGE_UNDEFINED, NULL };
473 
474   static
475   const SDF_Contour  null_contour = { { 0, 0 }, NULL, NULL };
476 
477   static
478   const SDF_Shape  null_shape = { NULL, NULL };
479 
480   static
481   const SDF_Signed_Distance  max_sdf = { INT_MAX, 0, 0 };
482 
483 
484   /* Create a new @SDF_Edge on the heap and assigns the `edge` */
485   /* pointer to the newly allocated memory.                    */
486   static FT_Error
sdf_edge_new(FT_Memory memory,SDF_Edge ** edge)487   sdf_edge_new( FT_Memory   memory,
488                 SDF_Edge**  edge )
489   {
490     FT_Error   error = FT_Err_Ok;
491     SDF_Edge*  ptr   = NULL;
492 
493 
494     if ( !memory || !edge )
495     {
496       error = FT_THROW( Invalid_Argument );
497       goto Exit;
498     }
499 
500     if ( !FT_QNEW( ptr ) )
501     {
502       *ptr = null_edge;
503       *edge = ptr;
504     }
505 
506   Exit:
507     return error;
508   }
509 
510 
511   /* Free the allocated `edge` variable. */
512   static void
sdf_edge_done(FT_Memory memory,SDF_Edge ** edge)513   sdf_edge_done( FT_Memory   memory,
514                  SDF_Edge**  edge )
515   {
516     if ( !memory || !edge || !*edge )
517       return;
518 
519     FT_FREE( *edge );
520   }
521 
522 
523   /* Create a new @SDF_Contour on the heap and assign     */
524   /* the `contour` pointer to the newly allocated memory. */
525   static FT_Error
sdf_contour_new(FT_Memory memory,SDF_Contour ** contour)526   sdf_contour_new( FT_Memory      memory,
527                    SDF_Contour**  contour )
528   {
529     FT_Error      error = FT_Err_Ok;
530     SDF_Contour*  ptr   = NULL;
531 
532 
533     if ( !memory || !contour )
534     {
535       error = FT_THROW( Invalid_Argument );
536       goto Exit;
537     }
538 
539     if ( !FT_QNEW( ptr ) )
540     {
541       *ptr     = null_contour;
542       *contour = ptr;
543     }
544 
545   Exit:
546     return error;
547   }
548 
549 
550   /* Free the allocated `contour` variable. */
551   /* Also free the list of edges.           */
552   static void
sdf_contour_done(FT_Memory memory,SDF_Contour ** contour)553   sdf_contour_done( FT_Memory      memory,
554                     SDF_Contour**  contour )
555   {
556     SDF_Edge*  edges;
557     SDF_Edge*  temp;
558 
559 
560     if ( !memory || !contour || !*contour )
561       return;
562 
563     edges = (*contour)->edges;
564 
565     /* release all edges */
566     while ( edges )
567     {
568       temp  = edges;
569       edges = edges->next;
570 
571       sdf_edge_done( memory, &temp );
572     }
573 
574     FT_FREE( *contour );
575   }
576 
577 
578   /* Create a new @SDF_Shape on the heap and assign     */
579   /* the `shape` pointer to the newly allocated memory. */
580   static FT_Error
sdf_shape_new(FT_Memory memory,SDF_Shape ** shape)581   sdf_shape_new( FT_Memory    memory,
582                  SDF_Shape**  shape )
583   {
584     FT_Error    error = FT_Err_Ok;
585     SDF_Shape*  ptr   = NULL;
586 
587 
588     if ( !memory || !shape )
589     {
590       error = FT_THROW( Invalid_Argument );
591       goto Exit;
592     }
593 
594     if ( !FT_QNEW( ptr ) )
595     {
596       *ptr        = null_shape;
597       ptr->memory = memory;
598       *shape      = ptr;
599     }
600 
601   Exit:
602     return error;
603   }
604 
605 
606   /* Free the allocated `shape` variable. */
607   /* Also free the list of contours.      */
608   static void
sdf_shape_done(SDF_Shape ** shape)609   sdf_shape_done( SDF_Shape**  shape )
610   {
611     FT_Memory     memory;
612     SDF_Contour*  contours;
613     SDF_Contour*  temp;
614 
615 
616     if ( !shape || !*shape )
617       return;
618 
619     memory   = (*shape)->memory;
620     contours = (*shape)->contours;
621 
622     if ( !memory )
623       return;
624 
625     /* release all contours */
626     while ( contours )
627     {
628       temp     = contours;
629       contours = contours->next;
630 
631       sdf_contour_done( memory, &temp );
632     }
633 
634     /* release the allocated shape struct  */
635     FT_FREE( *shape );
636   }
637 
638 
639   /**************************************************************************
640    *
641    * shape decomposition functions
642    *
643    */
644 
645   /* This function is called when starting a new contour at `to`, */
646   /* which gets added to the shape's list.                        */
647   static FT_Error
sdf_move_to(const FT_26D6_Vec * to,void * user)648   sdf_move_to( const FT_26D6_Vec* to,
649                void*              user )
650   {
651     SDF_Shape*    shape   = ( SDF_Shape* )user;
652     SDF_Contour*  contour = NULL;
653 
654     FT_Error   error  = FT_Err_Ok;
655     FT_Memory  memory = shape->memory;
656 
657 
658     if ( !to || !user )
659     {
660       error = FT_THROW( Invalid_Argument );
661       goto Exit;
662     }
663 
664     FT_CALL( sdf_contour_new( memory, &contour ) );
665 
666     contour->last_pos = *to;
667     contour->next     = shape->contours;
668     shape->contours   = contour;
669 
670   Exit:
671     return error;
672   }
673 
674 
675   /* This function is called when there is a line in the      */
676   /* contour.  The line starts at the previous edge point and */
677   /* stops at `to`.                                           */
678   static FT_Error
sdf_line_to(const FT_26D6_Vec * to,void * user)679   sdf_line_to( const FT_26D6_Vec*  to,
680                void*               user )
681   {
682     SDF_Shape*    shape    = ( SDF_Shape* )user;
683     SDF_Edge*     edge     = NULL;
684     SDF_Contour*  contour  = NULL;
685 
686     FT_Error      error    = FT_Err_Ok;
687     FT_Memory     memory   = shape->memory;
688 
689 
690     if ( !to || !user )
691     {
692       error = FT_THROW( Invalid_Argument );
693       goto Exit;
694     }
695 
696     contour = shape->contours;
697 
698     if ( contour->last_pos.x == to->x &&
699          contour->last_pos.y == to->y )
700       goto Exit;
701 
702     FT_CALL( sdf_edge_new( memory, &edge ) );
703 
704     edge->edge_type = SDF_EDGE_LINE;
705     edge->start_pos = contour->last_pos;
706     edge->end_pos   = *to;
707 
708     edge->next        = contour->edges;
709     contour->edges    = edge;
710     contour->last_pos = *to;
711 
712   Exit:
713     return error;
714   }
715 
716 
717   /* This function is called when there is a conic Bezier curve   */
718   /* in the contour.  The curve starts at the previous edge point */
719   /* and stops at `to`, with control point `control_1`.           */
720   static FT_Error
sdf_conic_to(const FT_26D6_Vec * control_1,const FT_26D6_Vec * to,void * user)721   sdf_conic_to( const FT_26D6_Vec*  control_1,
722                 const FT_26D6_Vec*  to,
723                 void*               user )
724   {
725     SDF_Shape*    shape    = ( SDF_Shape* )user;
726     SDF_Edge*     edge     = NULL;
727     SDF_Contour*  contour  = NULL;
728 
729     FT_Error   error  = FT_Err_Ok;
730     FT_Memory  memory = shape->memory;
731 
732 
733     if ( !control_1 || !to || !user )
734     {
735       error = FT_THROW( Invalid_Argument );
736       goto Exit;
737     }
738 
739     contour = shape->contours;
740 
741     /* If the control point coincides with any of the end points */
742     /* then it is a line and should be treated as one to avoid   */
743     /* unnecessary complexity later in the algorithm.            */
744     if ( ( contour->last_pos.x == control_1->x &&
745            contour->last_pos.y == control_1->y ) ||
746          ( control_1->x == to->x &&
747            control_1->y == to->y )               )
748     {
749       sdf_line_to( to, user );
750       goto Exit;
751     }
752 
753     FT_CALL( sdf_edge_new( memory, &edge ) );
754 
755     edge->edge_type = SDF_EDGE_CONIC;
756     edge->start_pos = contour->last_pos;
757     edge->control_a = *control_1;
758     edge->end_pos   = *to;
759 
760     edge->next        = contour->edges;
761     contour->edges    = edge;
762     contour->last_pos = *to;
763 
764   Exit:
765     return error;
766   }
767 
768 
769   /* This function is called when there is a cubic Bezier curve   */
770   /* in the contour.  The curve starts at the previous edge point */
771   /* and stops at `to`, with two control points `control_1` and   */
772   /* `control_2`.                                                 */
773   static FT_Error
sdf_cubic_to(const FT_26D6_Vec * control_1,const FT_26D6_Vec * control_2,const FT_26D6_Vec * to,void * user)774   sdf_cubic_to( const FT_26D6_Vec*  control_1,
775                 const FT_26D6_Vec*  control_2,
776                 const FT_26D6_Vec*  to,
777                 void*               user )
778   {
779     SDF_Shape*    shape   = ( SDF_Shape* )user;
780     SDF_Edge*     edge    = NULL;
781     SDF_Contour*  contour = NULL;
782 
783     FT_Error   error  = FT_Err_Ok;
784     FT_Memory  memory = shape->memory;
785 
786 
787     if ( !control_2 || !control_1 || !to || !user )
788     {
789       error = FT_THROW( Invalid_Argument );
790       goto Exit;
791     }
792 
793     contour = shape->contours;
794 
795     FT_CALL( sdf_edge_new( memory, &edge ) );
796 
797     edge->edge_type = SDF_EDGE_CUBIC;
798     edge->start_pos = contour->last_pos;
799     edge->control_a = *control_1;
800     edge->control_b = *control_2;
801     edge->end_pos   = *to;
802 
803     edge->next        = contour->edges;
804     contour->edges    = edge;
805     contour->last_pos = *to;
806 
807   Exit:
808     return error;
809   }
810 
811 
812   /* Construct the structure to hold all four outline */
813   /* decomposition functions.                         */
814   FT_DEFINE_OUTLINE_FUNCS(
815     sdf_decompose_funcs,
816 
817     (FT_Outline_MoveTo_Func) sdf_move_to,   /* move_to  */
818     (FT_Outline_LineTo_Func) sdf_line_to,   /* line_to  */
819     (FT_Outline_ConicTo_Func)sdf_conic_to,  /* conic_to */
820     (FT_Outline_CubicTo_Func)sdf_cubic_to,  /* cubic_to */
821 
822     0,                                      /* shift    */
823     0                                       /* delta    */
824   )
825 
826 
827   /* Decompose `outline` and put it into the `shape` structure.  */
828   static FT_Error
sdf_outline_decompose(FT_Outline * outline,SDF_Shape * shape)829   sdf_outline_decompose( FT_Outline*  outline,
830                          SDF_Shape*   shape )
831   {
832     FT_Error  error = FT_Err_Ok;
833 
834 
835     if ( !outline || !shape )
836     {
837       error = FT_THROW( Invalid_Argument );
838       goto Exit;
839     }
840 
841     error = FT_Outline_Decompose( outline,
842                                   &sdf_decompose_funcs,
843                                   (void*)shape );
844 
845   Exit:
846     return error;
847   }
848 
849 
850   /**************************************************************************
851    *
852    * utility functions
853    *
854    */
855 
856   /* Return the control box of an edge.  The control box is a rectangle */
857   /* in which all the control points can fit tightly.                   */
858   static FT_CBox
get_control_box(SDF_Edge edge)859   get_control_box( SDF_Edge  edge )
860   {
861     FT_CBox  cbox   = { 0, 0, 0, 0 };
862     FT_Bool  is_set = 0;
863 
864 
865     switch ( edge.edge_type )
866     {
867     case SDF_EDGE_CUBIC:
868       cbox.xMin = edge.control_b.x;
869       cbox.xMax = edge.control_b.x;
870       cbox.yMin = edge.control_b.y;
871       cbox.yMax = edge.control_b.y;
872 
873       is_set = 1;
874       FALL_THROUGH;
875 
876     case SDF_EDGE_CONIC:
877       if ( is_set )
878       {
879         cbox.xMin = edge.control_a.x < cbox.xMin
880                     ? edge.control_a.x
881                     : cbox.xMin;
882         cbox.xMax = edge.control_a.x > cbox.xMax
883                     ? edge.control_a.x
884                     : cbox.xMax;
885 
886         cbox.yMin = edge.control_a.y < cbox.yMin
887                     ? edge.control_a.y
888                     : cbox.yMin;
889         cbox.yMax = edge.control_a.y > cbox.yMax
890                     ? edge.control_a.y
891                     : cbox.yMax;
892       }
893       else
894       {
895         cbox.xMin = edge.control_a.x;
896         cbox.xMax = edge.control_a.x;
897         cbox.yMin = edge.control_a.y;
898         cbox.yMax = edge.control_a.y;
899 
900         is_set = 1;
901       }
902       FALL_THROUGH;
903 
904     case SDF_EDGE_LINE:
905       if ( is_set )
906       {
907         cbox.xMin = edge.start_pos.x < cbox.xMin
908                     ? edge.start_pos.x
909                     : cbox.xMin;
910         cbox.xMax = edge.start_pos.x > cbox.xMax
911                     ? edge.start_pos.x
912                     : cbox.xMax;
913 
914         cbox.yMin = edge.start_pos.y < cbox.yMin
915                     ? edge.start_pos.y
916                     : cbox.yMin;
917         cbox.yMax = edge.start_pos.y > cbox.yMax
918                     ? edge.start_pos.y
919                     : cbox.yMax;
920       }
921       else
922       {
923         cbox.xMin = edge.start_pos.x;
924         cbox.xMax = edge.start_pos.x;
925         cbox.yMin = edge.start_pos.y;
926         cbox.yMax = edge.start_pos.y;
927       }
928 
929       cbox.xMin = edge.end_pos.x < cbox.xMin
930                   ? edge.end_pos.x
931                   : cbox.xMin;
932       cbox.xMax = edge.end_pos.x > cbox.xMax
933                   ? edge.end_pos.x
934                   : cbox.xMax;
935 
936       cbox.yMin = edge.end_pos.y < cbox.yMin
937                   ? edge.end_pos.y
938                   : cbox.yMin;
939       cbox.yMax = edge.end_pos.y > cbox.yMax
940                   ? edge.end_pos.y
941                   : cbox.yMax;
942 
943       break;
944 
945     default:
946       break;
947     }
948 
949     return cbox;
950   }
951 
952 
953   /* Return orientation of a single contour.                    */
954   /* Note that the orientation is independent of the fill rule! */
955   /* So, for TTF a clockwise-oriented contour has to be filled  */
956   /* and the opposite for OTF fonts.                            */
957   static SDF_Contour_Orientation
get_contour_orientation(SDF_Contour * contour)958   get_contour_orientation ( SDF_Contour*  contour )
959   {
960     SDF_Edge*  head = NULL;
961     FT_26D6    area = 0;
962 
963 
964     /* return none if invalid parameters */
965     if ( !contour || !contour->edges )
966       return SDF_ORIENTATION_NONE;
967 
968     head = contour->edges;
969 
970     /* Calculate the area of the control box for all edges. */
971     while ( head )
972     {
973       switch ( head->edge_type )
974       {
975       case SDF_EDGE_LINE:
976         area += MUL_26D6( ( head->end_pos.x - head->start_pos.x ),
977                           ( head->end_pos.y + head->start_pos.y ) );
978         break;
979 
980       case SDF_EDGE_CONIC:
981         area += MUL_26D6( head->control_a.x - head->start_pos.x,
982                           head->control_a.y + head->start_pos.y );
983         area += MUL_26D6( head->end_pos.x - head->control_a.x,
984                           head->end_pos.y + head->control_a.y );
985         break;
986 
987       case SDF_EDGE_CUBIC:
988         area += MUL_26D6( head->control_a.x - head->start_pos.x,
989                           head->control_a.y + head->start_pos.y );
990         area += MUL_26D6( head->control_b.x - head->control_a.x,
991                           head->control_b.y + head->control_a.y );
992         area += MUL_26D6( head->end_pos.x - head->control_b.x,
993                           head->end_pos.y + head->control_b.y );
994         break;
995 
996       default:
997         return SDF_ORIENTATION_NONE;
998       }
999 
1000       head = head->next;
1001     }
1002 
1003     /* Clockwise contours cover a positive area, and counter-clockwise */
1004     /* contours cover a negative area.                                 */
1005     if ( area > 0 )
1006       return SDF_ORIENTATION_CW;
1007     else
1008       return SDF_ORIENTATION_CCW;
1009   }
1010 
1011 
1012   /* This function is exactly the same as the one */
1013   /* in the smooth renderer.  It splits a conic   */
1014   /* into two conics exactly half way at t = 0.5. */
1015   static void
split_conic(FT_26D6_Vec * base)1016   split_conic( FT_26D6_Vec*  base )
1017   {
1018     FT_26D6  a, b;
1019 
1020 
1021     base[4].x = base[2].x;
1022     a         = base[0].x + base[1].x;
1023     b         = base[1].x + base[2].x;
1024     base[3].x = b / 2;
1025     base[2].x = ( a + b ) / 4;
1026     base[1].x = a / 2;
1027 
1028     base[4].y = base[2].y;
1029     a         = base[0].y + base[1].y;
1030     b         = base[1].y + base[2].y;
1031     base[3].y = b / 2;
1032     base[2].y = ( a + b ) / 4;
1033     base[1].y = a / 2;
1034   }
1035 
1036 
1037   /* This function is exactly the same as the one */
1038   /* in the smooth renderer.  It splits a cubic   */
1039   /* into two cubics exactly half way at t = 0.5. */
1040   static void
split_cubic(FT_26D6_Vec * base)1041   split_cubic( FT_26D6_Vec*  base )
1042   {
1043     FT_26D6  a, b, c;
1044 
1045 
1046     base[6].x = base[3].x;
1047     a         = base[0].x + base[1].x;
1048     b         = base[1].x + base[2].x;
1049     c         = base[2].x + base[3].x;
1050     base[5].x = c / 2;
1051     c        += b;
1052     base[4].x = c / 4;
1053     base[1].x = a / 2;
1054     a        += b;
1055     base[2].x = a / 4;
1056     base[3].x = ( a + c ) / 8;
1057 
1058     base[6].y = base[3].y;
1059     a         = base[0].y + base[1].y;
1060     b         = base[1].y + base[2].y;
1061     c         = base[2].y + base[3].y;
1062     base[5].y = c / 2;
1063     c        += b;
1064     base[4].y = c / 4;
1065     base[1].y = a / 2;
1066     a        += b;
1067     base[2].y = a / 4;
1068     base[3].y = ( a + c ) / 8;
1069   }
1070 
1071 
1072   /* Split a conic Bezier curve into a number of lines */
1073   /* and add them to `out'.                            */
1074   /*                                                   */
1075   /* This function uses recursion; we thus need        */
1076   /* parameter `max_splits' for stopping.              */
1077   static FT_Error
split_sdf_conic(FT_Memory memory,FT_26D6_Vec * control_points,FT_UInt max_splits,SDF_Edge ** out)1078   split_sdf_conic( FT_Memory     memory,
1079                    FT_26D6_Vec*  control_points,
1080                    FT_UInt       max_splits,
1081                    SDF_Edge**    out )
1082   {
1083     FT_Error     error = FT_Err_Ok;
1084     FT_26D6_Vec  cpos[5];
1085     SDF_Edge*    left,*  right;
1086 
1087 
1088     if ( !memory || !out  )
1089     {
1090       error = FT_THROW( Invalid_Argument );
1091       goto Exit;
1092     }
1093 
1094     /* split conic outline */
1095     cpos[0] = control_points[0];
1096     cpos[1] = control_points[1];
1097     cpos[2] = control_points[2];
1098 
1099     split_conic( cpos );
1100 
1101     /* If max number of splits is done */
1102     /* then stop and add the lines to  */
1103     /* the list.                       */
1104     if ( max_splits <= 2 )
1105       goto Append;
1106 
1107     /* Otherwise keep splitting. */
1108     FT_CALL( split_sdf_conic( memory, &cpos[0], max_splits / 2, out ) );
1109     FT_CALL( split_sdf_conic( memory, &cpos[2], max_splits / 2, out ) );
1110 
1111     /* [NOTE]: This is not an efficient way of   */
1112     /* splitting the curve.  Check the deviation */
1113     /* instead and stop if the deviation is less */
1114     /* than a pixel.                             */
1115 
1116     goto Exit;
1117 
1118   Append:
1119     /* Do allocation and add the lines to the list. */
1120 
1121     FT_CALL( sdf_edge_new( memory, &left ) );
1122     FT_CALL( sdf_edge_new( memory, &right ) );
1123 
1124     left->start_pos  = cpos[0];
1125     left->end_pos    = cpos[2];
1126     left->edge_type  = SDF_EDGE_LINE;
1127 
1128     right->start_pos = cpos[2];
1129     right->end_pos   = cpos[4];
1130     right->edge_type = SDF_EDGE_LINE;
1131 
1132     left->next  = right;
1133     right->next = (*out);
1134     *out        = left;
1135 
1136   Exit:
1137     return error;
1138   }
1139 
1140 
1141   /* Split a cubic Bezier curve into a number of lines */
1142   /* and add them to `out`.                            */
1143   /*                                                   */
1144   /* This function uses recursion; we thus need        */
1145   /* parameter `max_splits' for stopping.              */
1146   static FT_Error
split_sdf_cubic(FT_Memory memory,FT_26D6_Vec * control_points,FT_UInt max_splits,SDF_Edge ** out)1147   split_sdf_cubic( FT_Memory     memory,
1148                    FT_26D6_Vec*  control_points,
1149                    FT_UInt       max_splits,
1150                    SDF_Edge**    out )
1151   {
1152     FT_Error       error = FT_Err_Ok;
1153     FT_26D6_Vec    cpos[7];
1154     SDF_Edge*      left, *right;
1155     const FT_26D6  threshold = ONE_PIXEL / 4;
1156 
1157 
1158     if ( !memory || !out )
1159     {
1160       error = FT_THROW( Invalid_Argument );
1161       goto Exit;
1162     }
1163 
1164     /* split the cubic */
1165     cpos[0] = control_points[0];
1166     cpos[1] = control_points[1];
1167     cpos[2] = control_points[2];
1168     cpos[3] = control_points[3];
1169 
1170     /* If the segment is flat enough we won't get any benefit by */
1171     /* splitting it further, so we can just stop splitting.      */
1172     /*                                                           */
1173     /* Check the deviation of the Bezier curve and stop if it is */
1174     /* smaller than the pre-defined `threshold` value.           */
1175     if ( FT_ABS( 2 * cpos[0].x - 3 * cpos[1].x + cpos[3].x ) < threshold &&
1176          FT_ABS( 2 * cpos[0].y - 3 * cpos[1].y + cpos[3].y ) < threshold &&
1177          FT_ABS( cpos[0].x - 3 * cpos[2].x + 2 * cpos[3].x ) < threshold &&
1178          FT_ABS( cpos[0].y - 3 * cpos[2].y + 2 * cpos[3].y ) < threshold )
1179     {
1180       split_cubic( cpos );
1181       goto Append;
1182     }
1183 
1184     split_cubic( cpos );
1185 
1186     /* If max number of splits is done */
1187     /* then stop and add the lines to  */
1188     /* the list.                       */
1189     if ( max_splits <= 2 )
1190       goto Append;
1191 
1192     /* Otherwise keep splitting. */
1193     FT_CALL( split_sdf_cubic( memory, &cpos[0], max_splits / 2, out ) );
1194     FT_CALL( split_sdf_cubic( memory, &cpos[3], max_splits / 2, out ) );
1195 
1196     /* [NOTE]: This is not an efficient way of   */
1197     /* splitting the curve.  Check the deviation */
1198     /* instead and stop if the deviation is less */
1199     /* than a pixel.                             */
1200 
1201     goto Exit;
1202 
1203   Append:
1204     /* Do allocation and add the lines to the list. */
1205 
1206     FT_CALL( sdf_edge_new( memory, &left) );
1207     FT_CALL( sdf_edge_new( memory, &right) );
1208 
1209     left->start_pos  = cpos[0];
1210     left->end_pos    = cpos[3];
1211     left->edge_type  = SDF_EDGE_LINE;
1212 
1213     right->start_pos = cpos[3];
1214     right->end_pos   = cpos[6];
1215     right->edge_type = SDF_EDGE_LINE;
1216 
1217     left->next  = right;
1218     right->next = (*out);
1219     *out        = left;
1220 
1221   Exit:
1222     return error;
1223   }
1224 
1225 
1226   /* Subdivide an entire shape into line segments */
1227   /* such that it doesn't look visually different */
1228   /* from the original curve.                     */
1229   static FT_Error
split_sdf_shape(SDF_Shape * shape)1230   split_sdf_shape( SDF_Shape*  shape )
1231   {
1232     FT_Error   error = FT_Err_Ok;
1233     FT_Memory  memory;
1234 
1235     SDF_Contour*  contours;
1236     SDF_Contour*  new_contours = NULL;
1237 
1238 
1239     if ( !shape || !shape->memory )
1240     {
1241       error = FT_THROW( Invalid_Argument );
1242       goto Exit;
1243     }
1244 
1245     contours = shape->contours;
1246     memory   = shape->memory;
1247 
1248     /* for each contour */
1249     while ( contours )
1250     {
1251       SDF_Edge*  edges     = contours->edges;
1252       SDF_Edge*  new_edges = NULL;
1253 
1254       SDF_Contour*  tempc;
1255 
1256 
1257       /* for each edge */
1258       while ( edges )
1259       {
1260         SDF_Edge*  edge = edges;
1261         SDF_Edge*  temp;
1262 
1263         switch ( edge->edge_type )
1264         {
1265         case SDF_EDGE_LINE:
1266           /* Just create a duplicate edge in case     */
1267           /* it is a line.  We can use the same edge. */
1268           FT_CALL( sdf_edge_new( memory, &temp ) );
1269 
1270           ft_memcpy( temp, edge, sizeof ( *edge ) );
1271 
1272           temp->next = new_edges;
1273           new_edges  = temp;
1274           break;
1275 
1276         case SDF_EDGE_CONIC:
1277           /* Subdivide the curve and add it to the list. */
1278           {
1279             FT_26D6_Vec  ctrls[3];
1280             FT_26D6      dx, dy;
1281             FT_UInt      num_splits;
1282 
1283 
1284             ctrls[0] = edge->start_pos;
1285             ctrls[1] = edge->control_a;
1286             ctrls[2] = edge->end_pos;
1287 
1288             dx = FT_ABS( ctrls[2].x + ctrls[0].x - 2 * ctrls[1].x );
1289             dy = FT_ABS( ctrls[2].y + ctrls[0].y - 2 * ctrls[1].y );
1290             if ( dx < dy )
1291               dx = dy;
1292 
1293             /* Calculate the number of necessary bisections.  Each      */
1294             /* bisection causes a four-fold reduction of the deviation, */
1295             /* hence we bisect the Bezier curve until the deviation     */
1296             /* becomes less than 1/8 of a pixel.  For more details      */
1297             /* check file `ftgrays.c`.                                  */
1298             num_splits = 1;
1299             while ( dx > ONE_PIXEL / 8 )
1300             {
1301               dx         >>= 2;
1302               num_splits <<= 1;
1303             }
1304 
1305             error = split_sdf_conic( memory, ctrls, num_splits, &new_edges );
1306           }
1307           break;
1308 
1309         case SDF_EDGE_CUBIC:
1310           /* Subdivide the curve and add it to the list. */
1311           {
1312             FT_26D6_Vec  ctrls[4];
1313 
1314 
1315             ctrls[0] = edge->start_pos;
1316             ctrls[1] = edge->control_a;
1317             ctrls[2] = edge->control_b;
1318             ctrls[3] = edge->end_pos;
1319 
1320             error = split_sdf_cubic( memory, ctrls, 32, &new_edges );
1321           }
1322           break;
1323 
1324         default:
1325           error = FT_THROW( Invalid_Argument );
1326         }
1327 
1328         if ( error != FT_Err_Ok )
1329           goto Exit;
1330 
1331         edges = edges->next;
1332       }
1333 
1334       /* add to the contours list */
1335       FT_CALL( sdf_contour_new( memory, &tempc ) );
1336 
1337       tempc->next  = new_contours;
1338       tempc->edges = new_edges;
1339       new_contours = tempc;
1340       new_edges    = NULL;
1341 
1342       /* deallocate the contour */
1343       tempc    = contours;
1344       contours = contours->next;
1345 
1346       sdf_contour_done( memory, &tempc );
1347     }
1348 
1349     shape->contours = new_contours;
1350 
1351   Exit:
1352     return error;
1353   }
1354 
1355 
1356   /**************************************************************************
1357    *
1358    * for debugging
1359    *
1360    */
1361 
1362 #ifdef FT_DEBUG_LEVEL_TRACE
1363 
1364   static void
sdf_shape_dump(SDF_Shape * shape)1365   sdf_shape_dump( SDF_Shape*  shape )
1366   {
1367     FT_UInt  num_contours = 0;
1368 
1369     FT_UInt  total_edges = 0;
1370     FT_UInt  total_lines = 0;
1371     FT_UInt  total_conic = 0;
1372     FT_UInt  total_cubic = 0;
1373 
1374     SDF_Contour*  contour_list;
1375 
1376 
1377     if ( !shape )
1378     {
1379       FT_TRACE5(( "sdf_shape_dump: null shape\n" ));
1380       return;
1381     }
1382 
1383     contour_list = shape->contours;
1384 
1385     FT_TRACE5(( "sdf_shape_dump (values are in 26.6 format):\n" ));
1386 
1387     while ( contour_list )
1388     {
1389       FT_UInt       num_edges = 0;
1390       SDF_Edge*     edge_list;
1391       SDF_Contour*  contour = contour_list;
1392 
1393 
1394       FT_TRACE5(( "  Contour %d\n", num_contours ));
1395 
1396       edge_list = contour->edges;
1397 
1398       while ( edge_list )
1399       {
1400         SDF_Edge*  edge = edge_list;
1401 
1402 
1403         FT_TRACE5(( "  %3d: ", num_edges ));
1404 
1405         switch ( edge->edge_type )
1406         {
1407         case SDF_EDGE_LINE:
1408           FT_TRACE5(( "Line:  (%ld, %ld) -- (%ld, %ld)\n",
1409                       edge->start_pos.x, edge->start_pos.y,
1410                       edge->end_pos.x, edge->end_pos.y ));
1411           total_lines++;
1412           break;
1413 
1414         case SDF_EDGE_CONIC:
1415           FT_TRACE5(( "Conic: (%ld, %ld) .. (%ld, %ld) .. (%ld, %ld)\n",
1416                       edge->start_pos.x, edge->start_pos.y,
1417                       edge->control_a.x, edge->control_a.y,
1418                       edge->end_pos.x, edge->end_pos.y ));
1419           total_conic++;
1420           break;
1421 
1422         case SDF_EDGE_CUBIC:
1423           FT_TRACE5(( "Cubic: (%ld, %ld) .. (%ld, %ld)"
1424                       " .. (%ld, %ld) .. (%ld %ld)\n",
1425                       edge->start_pos.x, edge->start_pos.y,
1426                       edge->control_a.x, edge->control_a.y,
1427                       edge->control_b.x, edge->control_b.y,
1428                       edge->end_pos.x, edge->end_pos.y ));
1429           total_cubic++;
1430           break;
1431 
1432         default:
1433           break;
1434         }
1435 
1436         num_edges++;
1437         total_edges++;
1438         edge_list = edge_list->next;
1439       }
1440 
1441       num_contours++;
1442       contour_list = contour_list->next;
1443     }
1444 
1445     FT_TRACE5(( "\n" ));
1446     FT_TRACE5(( "  total number of contours = %d\n", num_contours ));
1447     FT_TRACE5(( "  total number of edges    = %d\n", total_edges ));
1448     FT_TRACE5(( "    |__lines = %d\n", total_lines ));
1449     FT_TRACE5(( "    |__conic = %d\n", total_conic ));
1450     FT_TRACE5(( "    |__cubic = %d\n", total_cubic ));
1451   }
1452 
1453 #endif /* FT_DEBUG_LEVEL_TRACE */
1454 
1455 
1456   /**************************************************************************
1457    *
1458    * math functions
1459    *
1460    */
1461 
1462 #if !USE_NEWTON_FOR_CONIC
1463 
1464   /* [NOTE]: All the functions below down until rasterizer */
1465   /*         can be avoided if we decide to subdivide the  */
1466   /*         curve into lines.                             */
1467 
1468   /* This function uses Newton's iteration to find */
1469   /* the cube root of a fixed-point integer.       */
1470   static FT_16D16
cube_root(FT_16D16 val)1471   cube_root( FT_16D16  val )
1472   {
1473     /* [IMPORTANT]: This function is not good as it may */
1474     /* not break, so use a lookup table instead.  Or we */
1475     /* can use an algorithm similar to `square_root`.   */
1476 
1477     FT_Int  v, g, c;
1478 
1479 
1480     if ( val == 0                  ||
1481          val == -FT_INT_16D16( 1 ) ||
1482          val ==  FT_INT_16D16( 1 ) )
1483       return val;
1484 
1485     v = val < 0 ? -val : val;
1486     g = square_root( v );
1487     c = 0;
1488 
1489     while ( 1 )
1490     {
1491       c = FT_MulFix( FT_MulFix( g, g ), g ) - v;
1492       c = FT_DivFix( c, 3 * FT_MulFix( g, g ) );
1493 
1494       g -= c;
1495 
1496       if ( ( c < 0 ? -c : c ) < 30 )
1497         break;
1498     }
1499 
1500     return val < 0 ? -g : g;
1501   }
1502 
1503 
1504   /* Calculate the perpendicular by using '1 - base^2'. */
1505   /* Then use arctan to compute the angle.              */
1506   static FT_16D16
arc_cos(FT_16D16 val)1507   arc_cos( FT_16D16  val )
1508   {
1509     FT_16D16  p;
1510     FT_16D16  b   = val;
1511     FT_16D16  one = FT_INT_16D16( 1 );
1512 
1513 
1514     if ( b > one )
1515       b = one;
1516     if ( b < -one )
1517       b = -one;
1518 
1519     p = one - FT_MulFix( b, b );
1520     p = square_root( p );
1521 
1522     return FT_Atan2( b, p );
1523   }
1524 
1525 
1526   /* Compute roots of a quadratic polynomial, assign them to `out`, */
1527   /* and return number of real roots.                               */
1528   /*                                                                */
1529   /* The procedure can be found at                                  */
1530   /*                                                                */
1531   /*   https://mathworld.wolfram.com/QuadraticFormula.html          */
1532   static FT_UShort
solve_quadratic_equation(FT_26D6 a,FT_26D6 b,FT_26D6 c,FT_16D16 out[2])1533   solve_quadratic_equation( FT_26D6   a,
1534                             FT_26D6   b,
1535                             FT_26D6   c,
1536                             FT_16D16  out[2] )
1537   {
1538     FT_16D16  discriminant = 0;
1539 
1540 
1541     a = FT_26D6_16D16( a );
1542     b = FT_26D6_16D16( b );
1543     c = FT_26D6_16D16( c );
1544 
1545     if ( a == 0 )
1546     {
1547       if ( b == 0 )
1548         return 0;
1549       else
1550       {
1551         out[0] = FT_DivFix( -c, b );
1552 
1553         return 1;
1554       }
1555     }
1556 
1557     discriminant = FT_MulFix( b, b ) - 4 * FT_MulFix( a, c );
1558 
1559     if ( discriminant < 0 )
1560       return 0;
1561     else if ( discriminant == 0 )
1562     {
1563       out[0] = FT_DivFix( -b, 2 * a );
1564 
1565       return 1;
1566     }
1567     else
1568     {
1569       discriminant = square_root( discriminant );
1570 
1571       out[0] = FT_DivFix( -b + discriminant, 2 * a );
1572       out[1] = FT_DivFix( -b - discriminant, 2 * a );
1573 
1574       return 2;
1575     }
1576   }
1577 
1578 
1579   /* Compute roots of a cubic polynomial, assign them to `out`, */
1580   /* and return number of real roots.                           */
1581   /*                                                            */
1582   /* The procedure can be found at                              */
1583   /*                                                            */
1584   /*   https://mathworld.wolfram.com/CubicFormula.html          */
1585   static FT_UShort
solve_cubic_equation(FT_26D6 a,FT_26D6 b,FT_26D6 c,FT_26D6 d,FT_16D16 out[3])1586   solve_cubic_equation( FT_26D6   a,
1587                         FT_26D6   b,
1588                         FT_26D6   c,
1589                         FT_26D6   d,
1590                         FT_16D16  out[3] )
1591   {
1592     FT_16D16  q = 0;      /* intermediate */
1593     FT_16D16  r = 0;      /* intermediate */
1594 
1595     FT_16D16  a2 = b;     /* x^2 coefficients */
1596     FT_16D16  a1 = c;     /* x coefficients   */
1597     FT_16D16  a0 = d;     /* constant         */
1598 
1599     FT_16D16  q3   = 0;
1600     FT_16D16  r2   = 0;
1601     FT_16D16  a23  = 0;
1602     FT_16D16  a22  = 0;
1603     FT_16D16  a1x2 = 0;
1604 
1605 
1606     /* cutoff value for `a` to be a cubic, otherwise solve quadratic */
1607     if ( a == 0 || FT_ABS( a ) < 16 )
1608       return solve_quadratic_equation( b, c, d, out );
1609 
1610     if ( d == 0 )
1611     {
1612       out[0] = 0;
1613 
1614       return solve_quadratic_equation( a, b, c, out + 1 ) + 1;
1615     }
1616 
1617     /* normalize the coefficients; this also makes them 16.16 */
1618     a2 = FT_DivFix( a2, a );
1619     a1 = FT_DivFix( a1, a );
1620     a0 = FT_DivFix( a0, a );
1621 
1622     /* compute intermediates */
1623     a1x2 = FT_MulFix( a1, a2 );
1624     a22  = FT_MulFix( a2, a2 );
1625     a23  = FT_MulFix( a22, a2 );
1626 
1627     q = ( 3 * a1 - a22 ) / 9;
1628     r = ( 9 * a1x2 - 27 * a0 - 2 * a23 ) / 54;
1629 
1630     /* [BUG]: `q3` and `r2` still cause underflow. */
1631 
1632     q3 = FT_MulFix( q, q );
1633     q3 = FT_MulFix( q3, q );
1634 
1635     r2 = FT_MulFix( r, r );
1636 
1637     if ( q3 < 0 && r2 < -q3 )
1638     {
1639       FT_16D16  t = 0;
1640 
1641 
1642       q3 = square_root( -q3 );
1643       t  = FT_DivFix( r, q3 );
1644 
1645       if ( t > ( 1 << 16 ) )
1646         t =  ( 1 << 16 );
1647       if ( t < -( 1 << 16 ) )
1648         t = -( 1 << 16 );
1649 
1650       t   = arc_cos( t );
1651       a2 /= 3;
1652       q   = 2 * square_root( -q );
1653 
1654       out[0] = FT_MulFix( q, FT_Cos( t / 3 ) ) - a2;
1655       out[1] = FT_MulFix( q, FT_Cos( ( t + FT_ANGLE_PI * 2 ) / 3 ) ) - a2;
1656       out[2] = FT_MulFix( q, FT_Cos( ( t + FT_ANGLE_PI * 4 ) / 3 ) ) - a2;
1657 
1658       return 3;
1659     }
1660 
1661     else if ( r2 == -q3 )
1662     {
1663       FT_16D16  s = 0;
1664 
1665 
1666       s   = cube_root( r );
1667       a2 /= -3;
1668 
1669       out[0] = a2 + ( 2 * s );
1670       out[1] = a2 - s;
1671 
1672       return 2;
1673     }
1674 
1675     else
1676     {
1677       FT_16D16  s    = 0;
1678       FT_16D16  t    = 0;
1679       FT_16D16  dis  = 0;
1680 
1681 
1682       if ( q3 == 0 )
1683         dis = FT_ABS( r );
1684       else
1685         dis = square_root( q3 + r2 );
1686 
1687       s = cube_root( r + dis );
1688       t = cube_root( r - dis );
1689       a2 /= -3;
1690       out[0] = ( a2 + ( s + t ) );
1691 
1692       return 1;
1693     }
1694   }
1695 
1696 #endif /* !USE_NEWTON_FOR_CONIC */
1697 
1698 
1699   /*************************************************************************/
1700   /*************************************************************************/
1701   /**                                                                     **/
1702   /**  RASTERIZER                                                         **/
1703   /**                                                                     **/
1704   /*************************************************************************/
1705   /*************************************************************************/
1706 
1707   /**************************************************************************
1708    *
1709    * @Function:
1710    *   resolve_corner
1711    *
1712    * @Description:
1713    *   At some places on the grid two edges can give opposite directions;
1714    *   this happens when the closest point is on one of the endpoint.  In
1715    *   that case we need to check the proper sign.
1716    *
1717    *   This can be visualized by an example:
1718    *
1719    *   ```
1720    *                x
1721    *
1722    *                   o
1723    *                  ^ \
1724    *                 /   \
1725    *                /     \
1726    *           (a) /       \  (b)
1727    *              /         \
1728    *             /           \
1729    *            /             v
1730    *   ```
1731    *
1732    *   Suppose `x` is the point whose shortest distance from an arbitrary
1733    *   contour we want to find out.  It is clear that `o` is the nearest
1734    *   point on the contour.  Now to determine the sign we do a cross
1735    *   product of the shortest distance vector and the edge direction, i.e.,
1736    *
1737    *   ```
1738    *   => sign = cross(x - o, direction(a))
1739    *   ```
1740    *
1741    *   Using the right hand thumb rule we can see that the sign will be
1742    *   positive.
1743    *
1744    *   If we use `b', however, we have
1745    *
1746    *   ```
1747    *   => sign = cross(x - o, direction(b))
1748    *   ```
1749    *
1750    *   In this case the sign will be negative.  To determine the correct
1751    *   sign we thus divide the plane in two halves and check which plane the
1752    *   point lies in.
1753    *
1754    *   ```
1755    *                   |
1756    *                x  |
1757    *                   |
1758    *                   o
1759    *                  ^|\
1760    *                 / | \
1761    *                /  |  \
1762    *           (a) /   |   \  (b)
1763    *              /    |    \
1764    *             /           \
1765    *            /             v
1766    *   ```
1767    *
1768    *   We can see that `x` lies in the plane of `a`, so we take the sign
1769    *   determined by `a`.  This test can be easily done by calculating the
1770    *   orthogonality and taking the greater one.
1771    *
1772    *   The orthogonality is simply the sinus of the two vectors (i.e.,
1773    *   x - o) and the corresponding direction.  We efficiently pre-compute
1774    *   the orthogonality with the corresponding `get_min_distance_*`
1775    *   functions.
1776    *
1777    * @Input:
1778    *   sdf1 ::
1779    *     First signed distance (can be any of `a` or `b`).
1780    *
1781    *   sdf1 ::
1782    *     Second signed distance (can be any of `a` or `b`).
1783    *
1784    * @Return:
1785    *   The correct signed distance, which is computed by using the above
1786    *   algorithm.
1787    *
1788    * @Note:
1789    *   The function does not care about the actual distance, it simply
1790    *   returns the signed distance which has a larger cross product.  As a
1791    *   consequence, this function should not be used if the two distances
1792    *   are fairly apart.  In that case simply use the signed distance with
1793    *   a shorter absolute distance.
1794    *
1795    */
1796   static SDF_Signed_Distance
resolve_corner(SDF_Signed_Distance sdf1,SDF_Signed_Distance sdf2)1797   resolve_corner( SDF_Signed_Distance  sdf1,
1798                   SDF_Signed_Distance  sdf2 )
1799   {
1800     return FT_ABS( sdf1.cross ) > FT_ABS( sdf2.cross ) ? sdf1 : sdf2;
1801   }
1802 
1803 
1804   /**************************************************************************
1805    *
1806    * @Function:
1807    *   get_min_distance_line
1808    *
1809    * @Description:
1810    *   Find the shortest distance from the `line` segment to a given `point`
1811    *   and assign it to `out`.  Use it for line segments only.
1812    *
1813    * @Input:
1814    *   line ::
1815    *     The line segment to which the shortest distance is to be computed.
1816    *
1817    *   point ::
1818    *     Point from which the shortest distance is to be computed.
1819    *
1820    * @Output:
1821    *   out ::
1822    *     Signed distance from `point` to `line`.
1823    *
1824    * @Return:
1825    *   FreeType error, 0 means success.
1826    *
1827    * @Note:
1828    *   The `line' parameter must have an edge type of `SDF_EDGE_LINE`.
1829    *
1830    */
1831   static FT_Error
get_min_distance_line(SDF_Edge * line,FT_26D6_Vec point,SDF_Signed_Distance * out)1832   get_min_distance_line( SDF_Edge*             line,
1833                          FT_26D6_Vec           point,
1834                          SDF_Signed_Distance*  out )
1835   {
1836     /*
1837      * In order to calculate the shortest distance from a point to
1838      * a line segment, we do the following.  Let's assume that
1839      *
1840      * ```
1841      * a = start point of the line segment
1842      * b = end point of the line segment
1843      * p = point from which shortest distance is to be calculated
1844      * ```
1845      *
1846      * (1) Write the parametric equation of the line.
1847      *
1848      *     ```
1849      *     point_on_line = a + (b - a) * t   (t is the factor)
1850      *     ```
1851      *
1852      * (2) Find the projection of point `p` on the line.  The projection
1853      *     will be perpendicular to the line, which allows us to get the
1854      *     solution by making the dot product zero.
1855      *
1856      *     ```
1857      *     (point_on_line - a) . (p - point_on_line) = 0
1858      *
1859      *                (point_on_line)
1860      *      (a) x-------o----------------x (b)
1861      *                |_|
1862      *                  |
1863      *                  |
1864      *                 (p)
1865      *     ```
1866      *
1867      * (3) Simplification of the above equation yields the factor of
1868      *     `point_on_line`:
1869      *
1870      *     ```
1871      *     t = ((p - a) . (b - a)) / |b - a|^2
1872      *     ```
1873      *
1874      * (4) We clamp factor `t` between [0.0f, 1.0f] because `point_on_line`
1875      *     can be outside of the line segment:
1876      *
1877      *     ```
1878      *                                          (point_on_line)
1879      *     (a) x------------------------x (b) -----o---
1880      *                                           |_|
1881      *                                             |
1882      *                                             |
1883      *                                            (p)
1884      *     ```
1885      *
1886      * (5) Finally, the distance we are interested in is
1887      *
1888      *     ```
1889      *     |point_on_line - p|
1890      *     ```
1891      */
1892 
1893     FT_Error  error = FT_Err_Ok;
1894 
1895     FT_Vector  a;                   /* start position */
1896     FT_Vector  b;                   /* end position   */
1897     FT_Vector  p;                   /* current point  */
1898 
1899     FT_26D6_Vec  line_segment;      /* `b` - `a` */
1900     FT_26D6_Vec  p_sub_a;           /* `p` - `a` */
1901 
1902     FT_26D6   sq_line_length;       /* squared length of `line_segment` */
1903     FT_16D16  factor;               /* factor of the nearest point      */
1904     FT_26D6   cross;                /* used to determine sign           */
1905 
1906     FT_16D16_Vec  nearest_point;    /* `point_on_line`       */
1907     FT_16D16_Vec  nearest_vector;   /* `p` - `nearest_point` */
1908 
1909 
1910     if ( !line || !out )
1911     {
1912       error = FT_THROW( Invalid_Argument );
1913       goto Exit;
1914     }
1915 
1916     if ( line->edge_type != SDF_EDGE_LINE )
1917     {
1918       error = FT_THROW( Invalid_Argument );
1919       goto Exit;
1920     }
1921 
1922     a = line->start_pos;
1923     b = line->end_pos;
1924     p = point;
1925 
1926     line_segment.x = b.x - a.x;
1927     line_segment.y = b.y - a.y;
1928 
1929     p_sub_a.x = p.x - a.x;
1930     p_sub_a.y = p.y - a.y;
1931 
1932     sq_line_length = ( line_segment.x * line_segment.x ) / 64 +
1933                      ( line_segment.y * line_segment.y ) / 64;
1934 
1935     /* currently factor is 26.6 */
1936     factor = ( p_sub_a.x * line_segment.x ) / 64 +
1937              ( p_sub_a.y * line_segment.y ) / 64;
1938 
1939     /* now factor is 16.16 */
1940     factor = FT_DivFix( factor, sq_line_length );
1941 
1942     /* clamp the factor between 0.0 and 1.0 in fixed-point */
1943     if ( factor > FT_INT_16D16( 1 ) )
1944       factor = FT_INT_16D16( 1 );
1945     if ( factor < 0 )
1946       factor = 0;
1947 
1948     nearest_point.x = FT_MulFix( FT_26D6_16D16( line_segment.x ),
1949                                  factor );
1950     nearest_point.y = FT_MulFix( FT_26D6_16D16( line_segment.y ),
1951                                  factor );
1952 
1953     nearest_point.x = FT_26D6_16D16( a.x ) + nearest_point.x;
1954     nearest_point.y = FT_26D6_16D16( a.y ) + nearest_point.y;
1955 
1956     nearest_vector.x = nearest_point.x - FT_26D6_16D16( p.x );
1957     nearest_vector.y = nearest_point.y - FT_26D6_16D16( p.y );
1958 
1959     cross = FT_MulFix( nearest_vector.x, line_segment.y ) -
1960             FT_MulFix( nearest_vector.y, line_segment.x );
1961 
1962     /* assign the output */
1963     out->sign     = cross < 0 ? 1 : -1;
1964     out->distance = VECTOR_LENGTH_16D16( nearest_vector );
1965 
1966     /* Instead of finding `cross` for checking corner we */
1967     /* directly set it here.  This is more efficient     */
1968     /* because if the distance is perpendicular we can   */
1969     /* directly set it to 1.                             */
1970     if ( factor != 0 && factor != FT_INT_16D16( 1 ) )
1971       out->cross = FT_INT_16D16( 1 );
1972     else
1973     {
1974       /* [OPTIMIZATION]: Pre-compute this direction. */
1975       /* If not perpendicular then compute `cross`.  */
1976       FT_Vector_NormLen( &line_segment );
1977       FT_Vector_NormLen( &nearest_vector );
1978 
1979       out->cross = FT_MulFix( line_segment.x, nearest_vector.y ) -
1980                    FT_MulFix( line_segment.y, nearest_vector.x );
1981     }
1982 
1983   Exit:
1984     return error;
1985   }
1986 
1987 
1988   /**************************************************************************
1989    *
1990    * @Function:
1991    *   get_min_distance_conic
1992    *
1993    * @Description:
1994    *   Find the shortest distance from the `conic` Bezier curve to a given
1995    *   `point` and assign it to `out`.  Use it for conic/quadratic curves
1996    *   only.
1997    *
1998    * @Input:
1999    *   conic ::
2000    *     The conic Bezier curve to which the shortest distance is to be
2001    *     computed.
2002    *
2003    *   point ::
2004    *     Point from which the shortest distance is to be computed.
2005    *
2006    * @Output:
2007    *   out ::
2008    *     Signed distance from `point` to `conic`.
2009    *
2010    * @Return:
2011    *     FreeType error, 0 means success.
2012    *
2013    * @Note:
2014    *   The `conic` parameter must have an edge type of `SDF_EDGE_CONIC`.
2015    *
2016    */
2017 
2018 #if !USE_NEWTON_FOR_CONIC
2019 
2020   /*
2021    * The function uses an analytical method to find the shortest distance
2022    * which is faster than the Newton-Raphson method, but has underflows at
2023    * the moment.  Use Newton's method if you can see artifacts in the SDF.
2024    */
2025   static FT_Error
get_min_distance_conic(SDF_Edge * conic,FT_26D6_Vec point,SDF_Signed_Distance * out)2026   get_min_distance_conic( SDF_Edge*             conic,
2027                           FT_26D6_Vec           point,
2028                           SDF_Signed_Distance*  out )
2029   {
2030     /*
2031      * The procedure to find the shortest distance from a point to a
2032      * quadratic Bezier curve is similar to the line segment algorithm.  The
2033      * shortest distance is perpendicular to the Bezier curve; the only
2034      * difference from line is that there can be more than one
2035      * perpendicular, and we also have to check the endpoints, because the
2036      * perpendicular may not be the shortest.
2037      *
2038      * Let's assume that
2039      * ```
2040      * p0 = first endpoint
2041      * p1 = control point
2042      * p2 = second endpoint
2043      * p  = point from which shortest distance is to be calculated
2044      * ```
2045      *
2046      * (1) The equation of a quadratic Bezier curve can be written as
2047      *
2048      *     ```
2049      *     B(t) = (1 - t)^2 * p0 + 2(1 - t)t * p1 + t^2 * p2
2050      *     ```
2051      *
2052      *     with `t` a factor in the range [0.0f, 1.0f].  This equation can
2053      *     be rewritten as
2054      *
2055      *     ```
2056      *     B(t) = t^2 * (p0 - 2p1 + p2) + 2t * (p1 - p0) + p0
2057      *     ```
2058      *
2059      *     With
2060      *
2061      *     ```
2062      *     A = p0 - 2p1 + p2
2063      *     B = p1 - p0
2064      *     ```
2065      *
2066      *     we have
2067      *
2068      *     ```
2069      *     B(t) = t^2 * A + 2t * B + p0
2070      *     ```
2071      *
2072      * (2) The derivative of the last equation above is
2073      *
2074      *     ```
2075      *     B'(t) = 2 *(tA + B)
2076      *     ```
2077      *
2078      * (3) To find the shortest distance from `p` to `B(t)` we find the
2079      *     point on the curve at which the shortest distance vector (i.e.,
2080      *     `B(t) - p`) and the direction (i.e., `B'(t)`) make 90 degrees.
2081      *     In other words, we make the dot product zero.
2082      *
2083      *     ```
2084      *     (B(t) - p) . (B'(t)) = 0
2085      *     (t^2 * A + 2t * B + p0 - p) . (2 * (tA + B)) = 0
2086      *     ```
2087      *
2088      *     After simplifying we get a cubic equation
2089      *
2090      *     ```
2091      *     at^3 + bt^2 + ct + d = 0
2092      *     ```
2093      *
2094      *     with
2095      *
2096      *     ```
2097      *     a = A.A
2098      *     b = 3A.B
2099      *     c = 2B.B + A.p0 - A.p
2100      *     d = p0.B - p.B
2101      *     ```
2102      *
2103      * (4) Now the roots of the equation can be computed using 'Cardano's
2104      *     Cubic formula'; we clamp the roots in the range [0.0f, 1.0f].
2105      *
2106      * [note]: `B` and `B(t)` are different in the above equations.
2107      */
2108 
2109     FT_Error  error = FT_Err_Ok;
2110 
2111     FT_26D6_Vec  aA, bB;         /* A, B in the above comment             */
2112     FT_26D6_Vec  nearest_point = { 0, 0 };
2113                                  /* point on curve nearest to `point`     */
2114     FT_26D6_Vec  direction;      /* direction of curve at `nearest_point` */
2115 
2116     FT_26D6_Vec  p0, p1, p2;     /* control points of a conic curve       */
2117     FT_26D6_Vec  p;              /* `point` to which shortest distance    */
2118 
2119     FT_26D6  a, b, c, d;         /* cubic coefficients                    */
2120 
2121     FT_16D16  roots[3] = { 0, 0, 0 };    /* real roots of the cubic eq.   */
2122     FT_16D16  min_factor;                /* factor at `nearest_point`     */
2123     FT_16D16  cross;                     /* to determine the sign         */
2124     FT_16D16  min      = FT_INT_MAX;     /* shortest squared distance     */
2125 
2126     FT_UShort  num_roots;                /* number of real roots of cubic */
2127     FT_UShort  i;
2128 
2129 
2130     if ( !conic || !out )
2131     {
2132       error = FT_THROW( Invalid_Argument );
2133       goto Exit;
2134     }
2135 
2136     if ( conic->edge_type != SDF_EDGE_CONIC )
2137     {
2138       error = FT_THROW( Invalid_Argument );
2139       goto Exit;
2140     }
2141 
2142     p0 = conic->start_pos;
2143     p1 = conic->control_a;
2144     p2 = conic->end_pos;
2145     p  = point;
2146 
2147     /* compute substitution coefficients */
2148     aA.x = p0.x - 2 * p1.x + p2.x;
2149     aA.y = p0.y - 2 * p1.y + p2.y;
2150 
2151     bB.x = p1.x - p0.x;
2152     bB.y = p1.y - p0.y;
2153 
2154     /* compute cubic coefficients */
2155     a = VEC_26D6_DOT( aA, aA );
2156 
2157     b = 3 * VEC_26D6_DOT( aA, bB );
2158 
2159     c = 2 * VEC_26D6_DOT( bB, bB ) +
2160             VEC_26D6_DOT( aA, p0 ) -
2161             VEC_26D6_DOT( aA, p );
2162 
2163     d = VEC_26D6_DOT( p0, bB ) -
2164         VEC_26D6_DOT( p, bB );
2165 
2166     /* find the roots */
2167     num_roots = solve_cubic_equation( a, b, c, d, roots );
2168 
2169     if ( num_roots == 0 )
2170     {
2171       roots[0]  = 0;
2172       roots[1]  = FT_INT_16D16( 1 );
2173       num_roots = 2;
2174     }
2175 
2176     /* [OPTIMIZATION]: Check the roots, clamp them and discard */
2177     /*                 duplicate roots.                        */
2178 
2179     /* convert these values to 16.16 for further computation */
2180     aA.x = FT_26D6_16D16( aA.x );
2181     aA.y = FT_26D6_16D16( aA.y );
2182 
2183     bB.x = FT_26D6_16D16( bB.x );
2184     bB.y = FT_26D6_16D16( bB.y );
2185 
2186     p0.x = FT_26D6_16D16( p0.x );
2187     p0.y = FT_26D6_16D16( p0.y );
2188 
2189     p.x = FT_26D6_16D16( p.x );
2190     p.y = FT_26D6_16D16( p.y );
2191 
2192     for ( i = 0; i < num_roots; i++ )
2193     {
2194       FT_16D16  t    = roots[i];
2195       FT_16D16  t2   = 0;
2196       FT_16D16  dist = 0;
2197 
2198       FT_16D16_Vec  curve_point;
2199       FT_16D16_Vec  dist_vector;
2200 
2201       /*
2202        * Ideally we should discard the roots which are outside the range
2203        * [0.0, 1.0] and check the endpoints of the Bezier curve, but Behdad
2204        * Esfahbod proved the following lemma.
2205        *
2206        * Lemma:
2207        *
2208        * (1) If the closest point on the curve [0, 1] is to the endpoint at
2209        *     `t` = 1 and the cubic has no real roots at `t` = 1 then the
2210        *     cubic must have a real root at some `t` > 1.
2211        *
2212        * (2) Similarly, if the closest point on the curve [0, 1] is to the
2213        *     endpoint at `t` = 0 and the cubic has no real roots at `t` = 0
2214        *     then the cubic must have a real root at some `t` < 0.
2215        *
2216        * Now because of this lemma we only need to clamp the roots and that
2217        * will take care of the endpoints.
2218        *
2219        * For more details see
2220        *
2221        *   https://lists.nongnu.org/archive/html/freetype-devel/2020-06/msg00147.html
2222        */
2223 
2224       if ( t < 0 )
2225         t = 0;
2226       if ( t > FT_INT_16D16( 1 ) )
2227         t = FT_INT_16D16( 1 );
2228 
2229       t2 = FT_MulFix( t, t );
2230 
2231       /* B(t) = t^2 * A + 2t * B + p0 - p */
2232       curve_point.x = FT_MulFix( aA.x, t2 ) +
2233                       2 * FT_MulFix( bB.x, t ) + p0.x;
2234       curve_point.y = FT_MulFix( aA.y, t2 ) +
2235                       2 * FT_MulFix( bB.y, t ) + p0.y;
2236 
2237       /* `curve_point` - `p` */
2238       dist_vector.x = curve_point.x - p.x;
2239       dist_vector.y = curve_point.y - p.y;
2240 
2241       dist = VECTOR_LENGTH_16D16( dist_vector );
2242 
2243       if ( dist < min )
2244       {
2245         min           = dist;
2246         nearest_point = curve_point;
2247         min_factor    = t;
2248       }
2249     }
2250 
2251     /* B'(t) = 2 * (tA + B) */
2252     direction.x = 2 * FT_MulFix( aA.x, min_factor ) + 2 * bB.x;
2253     direction.y = 2 * FT_MulFix( aA.y, min_factor ) + 2 * bB.y;
2254 
2255     /* determine the sign */
2256     cross = FT_MulFix( nearest_point.x - p.x, direction.y ) -
2257             FT_MulFix( nearest_point.y - p.y, direction.x );
2258 
2259     /* assign the values */
2260     out->distance = min;
2261     out->sign     = cross < 0 ? 1 : -1;
2262 
2263     if ( min_factor != 0 && min_factor != FT_INT_16D16( 1 ) )
2264       out->cross = FT_INT_16D16( 1 );   /* the two are perpendicular */
2265     else
2266     {
2267       /* convert to nearest vector */
2268       nearest_point.x -= FT_26D6_16D16( p.x );
2269       nearest_point.y -= FT_26D6_16D16( p.y );
2270 
2271       /* compute `cross` if not perpendicular */
2272       FT_Vector_NormLen( &direction );
2273       FT_Vector_NormLen( &nearest_point );
2274 
2275       out->cross = FT_MulFix( direction.x, nearest_point.y ) -
2276                    FT_MulFix( direction.y, nearest_point.x );
2277     }
2278 
2279   Exit:
2280     return error;
2281   }
2282 
2283 #else /* USE_NEWTON_FOR_CONIC */
2284 
2285   /*
2286    * The function uses Newton's approximation to find the shortest distance,
2287    * which is a bit slower than the analytical method but doesn't cause
2288    * underflow.
2289    */
2290   static FT_Error
get_min_distance_conic(SDF_Edge * conic,FT_26D6_Vec point,SDF_Signed_Distance * out)2291   get_min_distance_conic( SDF_Edge*             conic,
2292                           FT_26D6_Vec           point,
2293                           SDF_Signed_Distance*  out )
2294   {
2295     /*
2296      * This method uses Newton-Raphson's approximation to find the shortest
2297      * distance from a point to a conic curve.  It does not involve solving
2298      * any cubic equation, that is why there is no risk of underflow.
2299      *
2300      * Let's assume that
2301      *
2302      * ```
2303      * p0 = first endpoint
2304      * p1 = control point
2305      * p3 = second endpoint
2306      * p  = point from which shortest distance is to be calculated
2307      * ```
2308      *
2309      * (1) The equation of a quadratic Bezier curve can be written as
2310      *
2311      *     ```
2312      *     B(t) = (1 - t)^2 * p0 + 2(1 - t)t * p1 + t^2 * p2
2313      *     ```
2314      *
2315      *     with `t` the factor in the range [0.0f, 1.0f].  The above
2316      *     equation can be rewritten as
2317      *
2318      *     ```
2319      *     B(t) = t^2 * (p0 - 2p1 + p2) + 2t * (p1 - p0) + p0
2320      *     ```
2321      *
2322      *     With
2323      *
2324      *     ```
2325      *     A = p0 - 2p1 + p2
2326      *     B = 2 * (p1 - p0)
2327      *     ```
2328      *
2329      *     we have
2330      *
2331      *     ```
2332      *     B(t) = t^2 * A + t * B + p0
2333      *     ```
2334      *
2335      * (2) The derivative of the above equation is
2336      *
2337      *     ```
2338      *     B'(t) = 2t * A + B
2339      *     ```
2340      *
2341      * (3) The second derivative of the above equation is
2342      *
2343      *     ```
2344      *     B''(t) = 2A
2345      *     ```
2346      *
2347      * (4) The equation `P(t)` of the distance from point `p` to the curve
2348      *     can be written as
2349      *
2350      *     ```
2351      *     P(t) = t^2 * A + t^2 * B + p0 - p
2352      *     ```
2353      *
2354      *     With
2355      *
2356      *     ```
2357      *     C = p0 - p
2358      *     ```
2359      *
2360      *     we have
2361      *
2362      *     ```
2363      *     P(t) = t^2 * A + t * B + C
2364      *     ```
2365      *
2366      * (5) Finally, the equation of the angle between `B(t)` and `P(t)` can
2367      *     be written as
2368      *
2369      *     ```
2370      *     Q(t) = P(t) . B'(t)
2371      *     ```
2372      *
2373      * (6) Our task is to find a value of `t` such that the above equation
2374      *     `Q(t)` becomes zero, that is, the point-to-curve vector makes
2375      *     90~degrees with the curve.  We solve this with the Newton-Raphson
2376      *     method.
2377      *
2378      * (7) We first assume an arbitrary value of factor `t`, which we then
2379      *     improve.
2380      *
2381      *     ```
2382      *     t := Q(t) / Q'(t)
2383      *     ```
2384      *
2385      *     Putting the value of `Q(t)` from the above equation gives
2386      *
2387      *     ```
2388      *     t := P(t) . B'(t) / derivative(P(t) . B'(t))
2389      *     t := P(t) . B'(t) /
2390      *            (P'(t) . B'(t) + P(t) . B''(t))
2391      *     ```
2392      *
2393      *     Note that `P'(t)` is the same as `B'(t)` because the constant is
2394      *     gone due to the derivative.
2395      *
2396      * (8) Finally we get the equation to improve the factor as
2397      *
2398      *     ```
2399      *     t := P(t) . B'(t) /
2400      *            (B'(t) . B'(t) + P(t) . B''(t))
2401      *     ```
2402      *
2403      * [note]: `B` and `B(t)` are different in the above equations.
2404      */
2405 
2406     FT_Error  error = FT_Err_Ok;
2407 
2408     FT_26D6_Vec  aA, bB, cC;     /* A, B, C in the above comment          */
2409     FT_26D6_Vec  nearest_point = { 0, 0 };
2410                                  /* point on curve nearest to `point`     */
2411     FT_26D6_Vec  direction;      /* direction of curve at `nearest_point` */
2412 
2413     FT_26D6_Vec  p0, p1, p2;     /* control points of a conic curve       */
2414     FT_26D6_Vec  p;              /* `point` to which shortest distance    */
2415 
2416     FT_16D16  min_factor = 0;            /* factor at `nearest_point'     */
2417     FT_16D16  cross;                     /* to determine the sign         */
2418     FT_16D16  min        = FT_INT_MAX;   /* shortest squared distance     */
2419 
2420     FT_UShort  iterations;
2421     FT_UShort  steps;
2422 
2423 
2424     if ( !conic || !out )
2425     {
2426       error = FT_THROW( Invalid_Argument );
2427       goto Exit;
2428     }
2429 
2430     if ( conic->edge_type != SDF_EDGE_CONIC )
2431     {
2432       error = FT_THROW( Invalid_Argument );
2433       goto Exit;
2434     }
2435 
2436     p0 = conic->start_pos;
2437     p1 = conic->control_a;
2438     p2 = conic->end_pos;
2439     p  = point;
2440 
2441     /* compute substitution coefficients */
2442     aA.x = p0.x - 2 * p1.x + p2.x;
2443     aA.y = p0.y - 2 * p1.y + p2.y;
2444 
2445     bB.x = 2 * ( p1.x - p0.x );
2446     bB.y = 2 * ( p1.y - p0.y );
2447 
2448     cC.x = p0.x;
2449     cC.y = p0.y;
2450 
2451     /* do Newton's iterations */
2452     for ( iterations = 0; iterations <= MAX_NEWTON_DIVISIONS; iterations++ )
2453     {
2454       FT_16D16  factor = FT_INT_16D16( iterations ) / MAX_NEWTON_DIVISIONS;
2455       FT_16D16  factor2;
2456       FT_16D16  length;
2457 
2458       FT_16D16_Vec  curve_point; /* point on the curve  */
2459       FT_16D16_Vec  dist_vector; /* `curve_point` - `p` */
2460 
2461       FT_26D6_Vec  d1;           /* first  derivative   */
2462       FT_26D6_Vec  d2;           /* second derivative   */
2463 
2464       FT_16D16  temp1;
2465       FT_16D16  temp2;
2466 
2467 
2468       for ( steps = 0; steps < MAX_NEWTON_STEPS; steps++ )
2469       {
2470         factor2 = FT_MulFix( factor, factor );
2471 
2472         /* B(t) = t^2 * A + t * B + p0 */
2473         curve_point.x = FT_MulFix( aA.x, factor2 ) +
2474                         FT_MulFix( bB.x, factor ) + cC.x;
2475         curve_point.y = FT_MulFix( aA.y, factor2 ) +
2476                         FT_MulFix( bB.y, factor ) + cC.y;
2477 
2478         /* convert to 16.16 */
2479         curve_point.x = FT_26D6_16D16( curve_point.x );
2480         curve_point.y = FT_26D6_16D16( curve_point.y );
2481 
2482         /* P(t) in the comment */
2483         dist_vector.x = curve_point.x - FT_26D6_16D16( p.x );
2484         dist_vector.y = curve_point.y - FT_26D6_16D16( p.y );
2485 
2486         length = VECTOR_LENGTH_16D16( dist_vector );
2487 
2488         if ( length < min )
2489         {
2490           min           = length;
2491           min_factor    = factor;
2492           nearest_point = curve_point;
2493         }
2494 
2495         /* This is Newton's approximation.          */
2496         /*                                          */
2497         /*   t := P(t) . B'(t) /                    */
2498         /*          (B'(t) . B'(t) + P(t) . B''(t)) */
2499 
2500         /* B'(t) = 2tA + B */
2501         d1.x = FT_MulFix( aA.x, 2 * factor ) + bB.x;
2502         d1.y = FT_MulFix( aA.y, 2 * factor ) + bB.y;
2503 
2504         /* B''(t) = 2A */
2505         d2.x = 2 * aA.x;
2506         d2.y = 2 * aA.y;
2507 
2508         dist_vector.x /= 1024;
2509         dist_vector.y /= 1024;
2510 
2511         /* temp1 = P(t) . B'(t) */
2512         temp1 = VEC_26D6_DOT( dist_vector, d1 );
2513 
2514         /* temp2 = B'(t) . B'(t) + P(t) . B''(t) */
2515         temp2 = VEC_26D6_DOT( d1, d1 ) +
2516                 VEC_26D6_DOT( dist_vector, d2 );
2517 
2518         factor -= FT_DivFix( temp1, temp2 );
2519 
2520         if ( factor < 0 || factor > FT_INT_16D16( 1 ) )
2521           break;
2522       }
2523     }
2524 
2525     /* B'(t) = 2t * A + B */
2526     direction.x = 2 * FT_MulFix( aA.x, min_factor ) + bB.x;
2527     direction.y = 2 * FT_MulFix( aA.y, min_factor ) + bB.y;
2528 
2529     /* determine the sign */
2530     cross = FT_MulFix( nearest_point.x - FT_26D6_16D16( p.x ),
2531                        direction.y )                           -
2532             FT_MulFix( nearest_point.y - FT_26D6_16D16( p.y ),
2533                        direction.x );
2534 
2535     /* assign the values */
2536     out->distance = min;
2537     out->sign     = cross < 0 ? 1 : -1;
2538 
2539     if ( min_factor != 0 && min_factor != FT_INT_16D16( 1 ) )
2540       out->cross = FT_INT_16D16( 1 );   /* the two are perpendicular */
2541     else
2542     {
2543       /* convert to nearest vector */
2544       nearest_point.x -= FT_26D6_16D16( p.x );
2545       nearest_point.y -= FT_26D6_16D16( p.y );
2546 
2547       /* compute `cross` if not perpendicular */
2548       FT_Vector_NormLen( &direction );
2549       FT_Vector_NormLen( &nearest_point );
2550 
2551       out->cross = FT_MulFix( direction.x, nearest_point.y ) -
2552                    FT_MulFix( direction.y, nearest_point.x );
2553     }
2554 
2555   Exit:
2556     return error;
2557   }
2558 
2559 
2560 #endif /* USE_NEWTON_FOR_CONIC */
2561 
2562 
2563   /**************************************************************************
2564    *
2565    * @Function:
2566    *   get_min_distance_cubic
2567    *
2568    * @Description:
2569    *   Find the shortest distance from the `cubic` Bezier curve to a given
2570    *   `point` and assigns it to `out`.  Use it for cubic curves only.
2571    *
2572    * @Input:
2573    *   cubic ::
2574    *     The cubic Bezier curve to which the shortest distance is to be
2575    *     computed.
2576    *
2577    *   point ::
2578    *     Point from which the shortest distance is to be computed.
2579    *
2580    * @Output:
2581    *   out ::
2582    *     Signed distance from `point` to `cubic`.
2583    *
2584    * @Return:
2585    *   FreeType error, 0 means success.
2586    *
2587    * @Note:
2588    *   The function uses Newton's approximation to find the shortest
2589    *   distance.  Another way would be to divide the cubic into conic or
2590    *   subdivide the curve into lines, but that is not implemented.
2591    *
2592    *   The `cubic` parameter must have an edge type of `SDF_EDGE_CUBIC`.
2593    *
2594    */
2595   static FT_Error
get_min_distance_cubic(SDF_Edge * cubic,FT_26D6_Vec point,SDF_Signed_Distance * out)2596   get_min_distance_cubic( SDF_Edge*             cubic,
2597                           FT_26D6_Vec           point,
2598                           SDF_Signed_Distance*  out )
2599   {
2600     /*
2601      * The procedure to find the shortest distance from a point to a cubic
2602      * Bezier curve is similar to quadratic curve algorithm.  The only
2603      * difference is that while calculating factor `t`, instead of a cubic
2604      * polynomial equation we have to find the roots of a 5th degree
2605      * polynomial equation.  Solving this would require a significant amount
2606      * of time, and still the results may not be accurate.  We are thus
2607      * going to directly approximate the value of `t` using the Newton-Raphson
2608      * method.
2609      *
2610      * Let's assume that
2611      *
2612      * ```
2613      * p0 = first endpoint
2614      * p1 = first control point
2615      * p2 = second control point
2616      * p3 = second endpoint
2617      * p  = point from which shortest distance is to be calculated
2618      * ```
2619      *
2620      * (1) The equation of a cubic Bezier curve can be written as
2621      *
2622      *     ```
2623      *     B(t) = (1 - t)^3 * p0 + 3(1 - t)^2 t * p1 +
2624      *              3(1 - t)t^2 * p2 + t^3 * p3
2625      *     ```
2626      *
2627      *     The equation can be expanded and written as
2628      *
2629      *     ```
2630      *     B(t) = t^3 * (-p0 + 3p1 - 3p2 + p3) +
2631      *              3t^2 * (p0 - 2p1 + p2) + 3t * (-p0 + p1) + p0
2632      *     ```
2633      *
2634      *     With
2635      *
2636      *     ```
2637      *     A = -p0 + 3p1 - 3p2 + p3
2638      *     B = 3(p0 - 2p1 + p2)
2639      *     C = 3(-p0 + p1)
2640      *     ```
2641      *
2642      *     we have
2643      *
2644      *     ```
2645      *     B(t) = t^3 * A + t^2 * B + t * C + p0
2646      *     ```
2647      *
2648      * (2) The derivative of the above equation is
2649      *
2650      *     ```
2651      *     B'(t) = 3t^2 * A + 2t * B + C
2652      *     ```
2653      *
2654      * (3) The second derivative of the above equation is
2655      *
2656      *     ```
2657      *     B''(t) = 6t * A + 2B
2658      *     ```
2659      *
2660      * (4) The equation `P(t)` of the distance from point `p` to the curve
2661      *     can be written as
2662      *
2663      *     ```
2664      *     P(t) = t^3 * A + t^2 * B + t * C + p0 - p
2665      *     ```
2666      *
2667      *     With
2668      *
2669      *     ```
2670      *     D = p0 - p
2671      *     ```
2672      *
2673      *     we have
2674      *
2675      *     ```
2676      *     P(t) = t^3 * A + t^2 * B + t * C + D
2677      *     ```
2678      *
2679      * (5) Finally the equation of the angle between `B(t)` and `P(t)` can
2680      *     be written as
2681      *
2682      *     ```
2683      *     Q(t) = P(t) . B'(t)
2684      *     ```
2685      *
2686      * (6) Our task is to find a value of `t` such that the above equation
2687      *     `Q(t)` becomes zero, that is, the point-to-curve vector makes
2688      *     90~degree with curve.  We solve this with the Newton-Raphson
2689      *     method.
2690      *
2691      * (7) We first assume an arbitrary value of factor `t`, which we then
2692      *     improve.
2693      *
2694      *     ```
2695      *     t := Q(t) / Q'(t)
2696      *     ```
2697      *
2698      *     Putting the value of `Q(t)` from the above equation gives
2699      *
2700      *     ```
2701      *     t := P(t) . B'(t) / derivative(P(t) . B'(t))
2702      *     t := P(t) . B'(t) /
2703      *            (P'(t) . B'(t) + P(t) . B''(t))
2704      *     ```
2705      *
2706      *     Note that `P'(t)` is the same as `B'(t)` because the constant is
2707      *     gone due to the derivative.
2708      *
2709      * (8) Finally we get the equation to improve the factor as
2710      *
2711      *     ```
2712      *     t := P(t) . B'(t) /
2713      *            (B'(t) . B'( t ) + P(t) . B''(t))
2714      *     ```
2715      *
2716      * [note]: `B` and `B(t)` are different in the above equations.
2717      */
2718 
2719     FT_Error  error = FT_Err_Ok;
2720 
2721     FT_26D6_Vec   aA, bB, cC, dD; /* A, B, C, D in the above comment       */
2722     FT_16D16_Vec  nearest_point = { 0, 0 };
2723                                   /* point on curve nearest to `point`     */
2724     FT_16D16_Vec  direction;      /* direction of curve at `nearest_point` */
2725 
2726     FT_26D6_Vec  p0, p1, p2, p3;  /* control points of a cubic curve       */
2727     FT_26D6_Vec  p;               /* `point` to which shortest distance    */
2728 
2729     FT_16D16  min_factor    = 0;            /* factor at shortest distance */
2730     FT_16D16  min_factor_sq = 0;            /* factor at shortest distance */
2731     FT_16D16  cross;                        /* to determine the sign       */
2732     FT_16D16  min           = FT_INT_MAX;   /* shortest distance           */
2733 
2734     FT_UShort  iterations;
2735     FT_UShort  steps;
2736 
2737 
2738     if ( !cubic || !out )
2739     {
2740       error = FT_THROW( Invalid_Argument );
2741       goto Exit;
2742     }
2743 
2744     if ( cubic->edge_type != SDF_EDGE_CUBIC )
2745     {
2746       error = FT_THROW( Invalid_Argument );
2747       goto Exit;
2748     }
2749 
2750     p0 = cubic->start_pos;
2751     p1 = cubic->control_a;
2752     p2 = cubic->control_b;
2753     p3 = cubic->end_pos;
2754     p  = point;
2755 
2756     /* compute substitution coefficients */
2757     aA.x = -p0.x + 3 * ( p1.x - p2.x ) + p3.x;
2758     aA.y = -p0.y + 3 * ( p1.y - p2.y ) + p3.y;
2759 
2760     bB.x = 3 * ( p0.x - 2 * p1.x + p2.x );
2761     bB.y = 3 * ( p0.y - 2 * p1.y + p2.y );
2762 
2763     cC.x = 3 * ( p1.x - p0.x );
2764     cC.y = 3 * ( p1.y - p0.y );
2765 
2766     dD.x = p0.x;
2767     dD.y = p0.y;
2768 
2769     for ( iterations = 0; iterations <= MAX_NEWTON_DIVISIONS; iterations++ )
2770     {
2771       FT_16D16  factor  = FT_INT_16D16( iterations ) / MAX_NEWTON_DIVISIONS;
2772 
2773       FT_16D16  factor2;         /* factor^2            */
2774       FT_16D16  factor3;         /* factor^3            */
2775       FT_16D16  length;
2776 
2777       FT_16D16_Vec  curve_point; /* point on the curve  */
2778       FT_16D16_Vec  dist_vector; /* `curve_point' - `p' */
2779 
2780       FT_26D6_Vec  d1;           /* first  derivative   */
2781       FT_26D6_Vec  d2;           /* second derivative   */
2782 
2783       FT_16D16  temp1;
2784       FT_16D16  temp2;
2785 
2786 
2787       for ( steps = 0; steps < MAX_NEWTON_STEPS; steps++ )
2788       {
2789         factor2 = FT_MulFix( factor, factor );
2790         factor3 = FT_MulFix( factor2, factor );
2791 
2792         /* B(t) = t^3 * A + t^2 * B + t * C + D */
2793         curve_point.x = FT_MulFix( aA.x, factor3 ) +
2794                         FT_MulFix( bB.x, factor2 ) +
2795                         FT_MulFix( cC.x, factor ) + dD.x;
2796         curve_point.y = FT_MulFix( aA.y, factor3 ) +
2797                         FT_MulFix( bB.y, factor2 ) +
2798                         FT_MulFix( cC.y, factor ) + dD.y;
2799 
2800         /* convert to 16.16 */
2801         curve_point.x = FT_26D6_16D16( curve_point.x );
2802         curve_point.y = FT_26D6_16D16( curve_point.y );
2803 
2804         /* P(t) in the comment */
2805         dist_vector.x = curve_point.x - FT_26D6_16D16( p.x );
2806         dist_vector.y = curve_point.y - FT_26D6_16D16( p.y );
2807 
2808         length = VECTOR_LENGTH_16D16( dist_vector );
2809 
2810         if ( length < min )
2811         {
2812           min           = length;
2813           min_factor    = factor;
2814           min_factor_sq = factor2;
2815           nearest_point = curve_point;
2816         }
2817 
2818         /* This the Newton's approximation.         */
2819         /*                                          */
2820         /*   t := P(t) . B'(t) /                    */
2821         /*          (B'(t) . B'(t) + P(t) . B''(t)) */
2822 
2823         /* B'(t) = 3t^2 * A + 2t * B + C */
2824         d1.x = FT_MulFix( aA.x, 3 * factor2 ) +
2825                FT_MulFix( bB.x, 2 * factor ) + cC.x;
2826         d1.y = FT_MulFix( aA.y, 3 * factor2 ) +
2827                FT_MulFix( bB.y, 2 * factor ) + cC.y;
2828 
2829         /* B''(t) = 6t * A + 2B */
2830         d2.x = FT_MulFix( aA.x, 6 * factor ) + 2 * bB.x;
2831         d2.y = FT_MulFix( aA.y, 6 * factor ) + 2 * bB.y;
2832 
2833         dist_vector.x /= 1024;
2834         dist_vector.y /= 1024;
2835 
2836         /* temp1 = P(t) . B'(t) */
2837         temp1 = VEC_26D6_DOT( dist_vector, d1 );
2838 
2839         /* temp2 = B'(t) . B'(t) + P(t) . B''(t) */
2840         temp2 = VEC_26D6_DOT( d1, d1 ) +
2841                 VEC_26D6_DOT( dist_vector, d2 );
2842 
2843         factor -= FT_DivFix( temp1, temp2 );
2844 
2845         if ( factor < 0 || factor > FT_INT_16D16( 1 ) )
2846           break;
2847       }
2848     }
2849 
2850     /* B'(t) = 3t^2 * A + 2t * B + C */
2851     direction.x = FT_MulFix( aA.x, 3 * min_factor_sq ) +
2852                   FT_MulFix( bB.x, 2 * min_factor ) + cC.x;
2853     direction.y = FT_MulFix( aA.y, 3 * min_factor_sq ) +
2854                   FT_MulFix( bB.y, 2 * min_factor ) + cC.y;
2855 
2856     /* determine the sign */
2857     cross = FT_MulFix( nearest_point.x - FT_26D6_16D16( p.x ),
2858                        direction.y )                           -
2859             FT_MulFix( nearest_point.y - FT_26D6_16D16( p.y ),
2860                        direction.x );
2861 
2862     /* assign the values */
2863     out->distance = min;
2864     out->sign     = cross < 0 ? 1 : -1;
2865 
2866     if ( min_factor != 0 && min_factor != FT_INT_16D16( 1 ) )
2867       out->cross = FT_INT_16D16( 1 );   /* the two are perpendicular */
2868     else
2869     {
2870       /* convert to nearest vector */
2871       nearest_point.x -= FT_26D6_16D16( p.x );
2872       nearest_point.y -= FT_26D6_16D16( p.y );
2873 
2874       /* compute `cross` if not perpendicular */
2875       FT_Vector_NormLen( &direction );
2876       FT_Vector_NormLen( &nearest_point );
2877 
2878       out->cross = FT_MulFix( direction.x, nearest_point.y ) -
2879                    FT_MulFix( direction.y, nearest_point.x );
2880     }
2881 
2882   Exit:
2883     return error;
2884   }
2885 
2886 
2887   /**************************************************************************
2888    *
2889    * @Function:
2890    *   sdf_edge_get_min_distance
2891    *
2892    * @Description:
2893    *   Find shortest distance from `point` to any type of `edge`.  It checks
2894    *   the edge type and then calls the relevant `get_min_distance_*`
2895    *   function.
2896    *
2897    * @Input:
2898    *   edge ::
2899    *     An edge to which the shortest distance is to be computed.
2900    *
2901    *   point ::
2902    *     Point from which the shortest distance is to be computed.
2903    *
2904    * @Output:
2905    *   out ::
2906    *     Signed distance from `point` to `edge`.
2907    *
2908    * @Return:
2909    *   FreeType error, 0 means success.
2910    *
2911    */
2912   static FT_Error
sdf_edge_get_min_distance(SDF_Edge * edge,FT_26D6_Vec point,SDF_Signed_Distance * out)2913   sdf_edge_get_min_distance( SDF_Edge*             edge,
2914                              FT_26D6_Vec           point,
2915                              SDF_Signed_Distance*  out )
2916   {
2917     FT_Error  error = FT_Err_Ok;
2918 
2919 
2920     if ( !edge || !out )
2921     {
2922       error = FT_THROW( Invalid_Argument );
2923       goto Exit;
2924     }
2925 
2926     /* edge-specific distance calculation */
2927     switch ( edge->edge_type )
2928     {
2929     case SDF_EDGE_LINE:
2930       get_min_distance_line( edge, point, out );
2931       break;
2932 
2933     case SDF_EDGE_CONIC:
2934       get_min_distance_conic( edge, point, out );
2935       break;
2936 
2937     case SDF_EDGE_CUBIC:
2938       get_min_distance_cubic( edge, point, out );
2939       break;
2940 
2941     default:
2942       error = FT_THROW( Invalid_Argument );
2943     }
2944 
2945   Exit:
2946     return error;
2947   }
2948 
2949 
2950   /* `sdf_generate' is not used at the moment */
2951 #if 0
2952 
2953   #error "DO NOT USE THIS!"
2954   #error "The function still outputs 16-bit data, which might cause memory"
2955   #error "corruption.  If required I will add this later."
2956 
2957   /**************************************************************************
2958    *
2959    * @Function:
2960    *   sdf_contour_get_min_distance
2961    *
2962    * @Description:
2963    *   Iterate over all edges that make up the contour, find the shortest
2964    *   distance from a point to this contour, and assigns result to `out`.
2965    *
2966    * @Input:
2967    *   contour ::
2968    *     A contour to which the shortest distance is to be computed.
2969    *
2970    *   point ::
2971    *     Point from which the shortest distance is to be computed.
2972    *
2973    * @Output:
2974    *   out ::
2975    *     Signed distance from the `point' to the `contour'.
2976    *
2977    * @Return:
2978    *   FreeType error, 0 means success.
2979    *
2980    * @Note:
2981    *   The function does not return a signed distance for each edge which
2982    *   makes up the contour, it simply returns the shortest of all the
2983    *   edges.
2984    *
2985    */
2986   static FT_Error
sdf_contour_get_min_distance(SDF_Contour * contour,FT_26D6_Vec point,SDF_Signed_Distance * out)2987   sdf_contour_get_min_distance( SDF_Contour*          contour,
2988                                 FT_26D6_Vec           point,
2989                                 SDF_Signed_Distance*  out )
2990   {
2991     FT_Error             error    = FT_Err_Ok;
2992     SDF_Signed_Distance  min_dist = max_sdf;
2993     SDF_Edge*            edge_list;
2994 
2995 
2996     if ( !contour || !out )
2997     {
2998       error = FT_THROW( Invalid_Argument );
2999       goto Exit;
3000     }
3001 
3002     edge_list = contour->edges;
3003 
3004     /* iterate over all the edges manually */
3005     while ( edge_list )
3006     {
3007       SDF_Signed_Distance  current_dist = max_sdf;
3008       FT_16D16             diff;
3009 
3010 
3011       FT_CALL( sdf_edge_get_min_distance( edge_list,
3012                                           point,
3013                                           &current_dist ) );
3014 
3015       if ( current_dist.distance >= 0 )
3016       {
3017         diff = current_dist.distance - min_dist.distance;
3018 
3019 
3020         if ( FT_ABS( diff ) < CORNER_CHECK_EPSILON )
3021           min_dist = resolve_corner( min_dist, current_dist );
3022         else if ( diff < 0 )
3023           min_dist = current_dist;
3024       }
3025       else
3026         FT_TRACE0(( "sdf_contour_get_min_distance: Overflow.\n" ));
3027 
3028       edge_list = edge_list->next;
3029     }
3030 
3031     *out = min_dist;
3032 
3033   Exit:
3034     return error;
3035   }
3036 
3037 
3038   /**************************************************************************
3039    *
3040    * @Function:
3041    *   sdf_generate
3042    *
3043    * @Description:
3044    *   This is the main function that is responsible for generating signed
3045    *   distance fields.  The function does not align or compute the size of
3046    *   `bitmap`; therefore the calling application must set up `bitmap`
3047    *   properly and transform the `shape' appropriately in advance.
3048    *
3049    *   Currently we check all pixels against all contours and all edges.
3050    *
3051    * @Input:
3052    *   internal_params ::
3053    *     Internal parameters and properties required by the rasterizer.  See
3054    *     @SDF_Params for more.
3055    *
3056    *   shape ::
3057    *     A complete shape which is used to generate SDF.
3058    *
3059    *   spread ::
3060    *     Maximum distances to be allowed in the output bitmap.
3061    *
3062    * @Output:
3063    *   bitmap ::
3064    *     The output bitmap which will contain the SDF information.
3065    *
3066    * @Return:
3067    *   FreeType error, 0 means success.
3068    *
3069    */
3070   static FT_Error
sdf_generate(const SDF_Params internal_params,const SDF_Shape * shape,FT_UInt spread,const FT_Bitmap * bitmap)3071   sdf_generate( const SDF_Params  internal_params,
3072                 const SDF_Shape*  shape,
3073                 FT_UInt           spread,
3074                 const FT_Bitmap*  bitmap )
3075   {
3076     FT_Error  error = FT_Err_Ok;
3077 
3078     FT_UInt  width = 0;
3079     FT_UInt  rows  = 0;
3080     FT_UInt  x     = 0;   /* used to loop in x direction, i.e., width     */
3081     FT_UInt  y     = 0;   /* used to loop in y direction, i.e., rows      */
3082     FT_UInt  sp_sq = 0;   /* `spread` [* `spread`] as a 16.16 fixed value */
3083 
3084     FT_Short*  buffer;
3085 
3086 
3087     if ( !shape || !bitmap )
3088     {
3089       error = FT_THROW( Invalid_Argument );
3090       goto Exit;
3091     }
3092 
3093     if ( spread < MIN_SPREAD || spread > MAX_SPREAD )
3094     {
3095       error = FT_THROW( Invalid_Argument );
3096       goto Exit;
3097     }
3098 
3099     width  = bitmap->width;
3100     rows   = bitmap->rows;
3101     buffer = (FT_Short*)bitmap->buffer;
3102 
3103     if ( USE_SQUARED_DISTANCES )
3104       sp_sq = FT_INT_16D16( spread * spread );
3105     else
3106       sp_sq = FT_INT_16D16( spread );
3107 
3108     if ( width == 0 || rows == 0 )
3109     {
3110       FT_TRACE0(( "sdf_generate:"
3111                   " Cannot render glyph with width/height == 0\n" ));
3112       FT_TRACE0(( "             "
3113                   " (width, height provided [%d, %d])\n",
3114                   width, rows ));
3115 
3116       error = FT_THROW( Cannot_Render_Glyph );
3117       goto Exit;
3118     }
3119 
3120     /* loop over all rows */
3121     for ( y = 0; y < rows; y++ )
3122     {
3123       /* loop over all pixels of a row */
3124       for ( x = 0; x < width; x++ )
3125       {
3126         /* `grid_point` is the current pixel position; */
3127         /* our task is to find the shortest distance   */
3128         /* from this point to the entire shape.        */
3129         FT_26D6_Vec          grid_point = zero_vector;
3130         SDF_Signed_Distance  min_dist   = max_sdf;
3131         SDF_Contour*         contour_list;
3132 
3133         FT_UInt   index;
3134         FT_Short  value;
3135 
3136 
3137         grid_point.x = FT_INT_26D6( x );
3138         grid_point.y = FT_INT_26D6( y );
3139 
3140         /* This `grid_point' is at the corner, but we */
3141         /* use the center of the pixel.               */
3142         grid_point.x += FT_INT_26D6( 1 ) / 2;
3143         grid_point.y += FT_INT_26D6( 1 ) / 2;
3144 
3145         contour_list = shape->contours;
3146 
3147         /* iterate over all contours manually */
3148         while ( contour_list )
3149         {
3150           SDF_Signed_Distance  current_dist = max_sdf;
3151 
3152 
3153           FT_CALL( sdf_contour_get_min_distance( contour_list,
3154                                                  grid_point,
3155                                                  &current_dist ) );
3156 
3157           if ( current_dist.distance < min_dist.distance )
3158             min_dist = current_dist;
3159 
3160           contour_list = contour_list->next;
3161         }
3162 
3163         /* [OPTIMIZATION]: if (min_dist > sp_sq) then simply clamp  */
3164         /*                 the value to spread to avoid square_root */
3165 
3166         /* clamp the values to spread */
3167         if ( min_dist.distance > sp_sq )
3168           min_dist.distance = sp_sq;
3169 
3170         /* square_root the values and fit in a 6.10 fixed-point */
3171         if ( USE_SQUARED_DISTANCES )
3172           min_dist.distance = square_root( min_dist.distance );
3173 
3174         if ( internal_params.orientation == FT_ORIENTATION_FILL_LEFT )
3175           min_dist.sign = -min_dist.sign;
3176         if ( internal_params.flip_sign )
3177           min_dist.sign = -min_dist.sign;
3178 
3179         min_dist.distance /= 64; /* convert from 16.16 to 22.10 */
3180 
3181         value  = min_dist.distance & 0x0000FFFF; /* truncate to 6.10 */
3182         value *= min_dist.sign;
3183 
3184         if ( internal_params.flip_y )
3185           index = y * width + x;
3186         else
3187           index = ( rows - y - 1 ) * width + x;
3188 
3189         buffer[index] = value;
3190       }
3191     }
3192 
3193   Exit:
3194     return error;
3195   }
3196 
3197 #endif /* 0 */
3198 
3199 
3200   /**************************************************************************
3201    *
3202    * @Function:
3203    *   sdf_generate_bounding_box
3204    *
3205    * @Description:
3206    *   This function does basically the same thing as `sdf_generate` above
3207    *   but more efficiently.
3208    *
3209    *   Instead of checking all pixels against all edges, we loop over all
3210    *   edges and only check pixels around the control box of the edge; the
3211    *   control box is increased by the spread in all directions.  Anything
3212    *   outside of the control box that exceeds `spread` doesn't need to be
3213    *   computed.
3214    *
3215    *   Lastly, to determine the sign of unchecked pixels, we do a single
3216    *   pass of all rows starting with a '+' sign and flipping when we come
3217    *   across a '-' sign and continue.  This also eliminates the possibility
3218    *   of overflow because we only check the proximity of the curve.
3219    *   Therefore we can use squared distanced safely.
3220    *
3221    * @Input:
3222    *   internal_params ::
3223    *     Internal parameters and properties required by the rasterizer.
3224    *     See @SDF_Params for more.
3225    *
3226    *   shape ::
3227    *     A complete shape which is used to generate SDF.
3228    *
3229    *   spread ::
3230    *     Maximum distances to be allowed in the output bitmap.
3231    *
3232    * @Output:
3233    *   bitmap ::
3234    *     The output bitmap which will contain the SDF information.
3235    *
3236    * @Return:
3237    *   FreeType error, 0 means success.
3238    *
3239    */
3240   static FT_Error
sdf_generate_bounding_box(const SDF_Params internal_params,const SDF_Shape * shape,FT_UInt spread,const FT_Bitmap * bitmap)3241   sdf_generate_bounding_box( const SDF_Params  internal_params,
3242                              const SDF_Shape*  shape,
3243                              FT_UInt           spread,
3244                              const FT_Bitmap*  bitmap )
3245   {
3246     FT_Error   error  = FT_Err_Ok;
3247     FT_Memory  memory = NULL;
3248 
3249     FT_Int  width, rows, i, j;
3250     FT_Int  sp_sq;            /* max value to check   */
3251 
3252     SDF_Contour*   contours;  /* list of all contours */
3253     FT_SDFFormat*  buffer;    /* the bitmap buffer    */
3254 
3255     /* This buffer has the same size in indices as the    */
3256     /* bitmap buffer.  When we check a pixel position for */
3257     /* a shortest distance we keep it in this buffer.     */
3258     /* This way we can find out which pixel is set,       */
3259     /* and also determine the signs properly.             */
3260     SDF_Signed_Distance*  dists = NULL;
3261 
3262     const FT_16D16  fixed_spread = (FT_16D16)FT_INT_16D16( spread );
3263 
3264 
3265     if ( !shape || !bitmap )
3266     {
3267       error = FT_THROW( Invalid_Argument );
3268       goto Exit;
3269     }
3270 
3271     if ( spread < MIN_SPREAD || spread > MAX_SPREAD )
3272     {
3273       error = FT_THROW( Invalid_Argument );
3274       goto Exit;
3275     }
3276 
3277     memory = shape->memory;
3278     if ( !memory )
3279     {
3280       error = FT_THROW( Invalid_Argument );
3281       goto Exit;
3282     }
3283 
3284     if ( FT_ALLOC( dists,
3285                    bitmap->width * bitmap->rows * sizeof ( *dists ) ) )
3286       goto Exit;
3287 
3288     contours = shape->contours;
3289     width    = (FT_Int)bitmap->width;
3290     rows     = (FT_Int)bitmap->rows;
3291     buffer   = (FT_SDFFormat*)bitmap->buffer;
3292 
3293     if ( USE_SQUARED_DISTANCES )
3294       sp_sq = FT_INT_16D16( (FT_Int)( spread * spread ) );
3295     else
3296       sp_sq = fixed_spread;
3297 
3298     if ( width == 0 || rows == 0 )
3299     {
3300       FT_TRACE0(( "sdf_generate:"
3301                   " Cannot render glyph with width/height == 0\n" ));
3302       FT_TRACE0(( "             "
3303                   " (width, height provided [%d, %d])", width, rows ));
3304 
3305       error = FT_THROW( Cannot_Render_Glyph );
3306       goto Exit;
3307     }
3308 
3309     /* loop over all contours */
3310     while ( contours )
3311     {
3312       SDF_Edge*  edges = contours->edges;
3313 
3314 
3315       /* loop over all edges */
3316       while ( edges )
3317       {
3318         FT_CBox  cbox;
3319         FT_Int   x, y;
3320 
3321 
3322         /* get the control box and increase it by `spread' */
3323         cbox = get_control_box( *edges );
3324 
3325         cbox.xMin = ( cbox.xMin - 63 ) / 64 - ( FT_Pos )spread;
3326         cbox.xMax = ( cbox.xMax + 63 ) / 64 + ( FT_Pos )spread;
3327         cbox.yMin = ( cbox.yMin - 63 ) / 64 - ( FT_Pos )spread;
3328         cbox.yMax = ( cbox.yMax + 63 ) / 64 + ( FT_Pos )spread;
3329 
3330         /* now loop over the pixels in the control box. */
3331         for ( y = cbox.yMin; y < cbox.yMax; y++ )
3332         {
3333           for ( x = cbox.xMin; x < cbox.xMax; x++ )
3334           {
3335             FT_26D6_Vec          grid_point = zero_vector;
3336             SDF_Signed_Distance  dist       = max_sdf;
3337             FT_UInt              index      = 0;
3338             FT_16D16             diff       = 0;
3339 
3340 
3341             if ( x < 0 || x >= width )
3342               continue;
3343             if ( y < 0 || y >= rows )
3344               continue;
3345 
3346             grid_point.x = FT_INT_26D6( x );
3347             grid_point.y = FT_INT_26D6( y );
3348 
3349             /* This `grid_point` is at the corner, but we */
3350             /* use the center of the pixel.               */
3351             grid_point.x += FT_INT_26D6( 1 ) / 2;
3352             grid_point.y += FT_INT_26D6( 1 ) / 2;
3353 
3354             FT_CALL( sdf_edge_get_min_distance( edges,
3355                                                 grid_point,
3356                                                 &dist ) );
3357 
3358             if ( internal_params.orientation == FT_ORIENTATION_FILL_LEFT )
3359               dist.sign = -dist.sign;
3360 
3361             /* ignore if the distance is greater than spread;       */
3362             /* otherwise it creates artifacts due to the wrong sign */
3363             if ( dist.distance > sp_sq )
3364               continue;
3365 
3366             /* take the square root of the distance if required */
3367             if ( USE_SQUARED_DISTANCES )
3368               dist.distance = square_root( dist.distance );
3369 
3370             if ( internal_params.flip_y )
3371               index = (FT_UInt)( y * width + x );
3372             else
3373               index = (FT_UInt)( ( rows - y - 1 ) * width + x );
3374 
3375             /* check whether the pixel is set or not */
3376             if ( dists[index].sign == 0 )
3377               dists[index] = dist;
3378             else
3379             {
3380               diff = FT_ABS( dists[index].distance - dist.distance );
3381 
3382               if ( diff <= CORNER_CHECK_EPSILON )
3383                 dists[index] = resolve_corner( dists[index], dist );
3384               else if ( dists[index].distance > dist.distance )
3385                 dists[index] = dist;
3386             }
3387           }
3388         }
3389 
3390         edges = edges->next;
3391       }
3392 
3393       contours = contours->next;
3394     }
3395 
3396     /* final pass */
3397     for ( j = 0; j < rows; j++ )
3398     {
3399       /* We assume the starting pixel of each row is outside. */
3400       FT_Char  current_sign = -1;
3401       FT_UInt  index;
3402 
3403 
3404       if ( internal_params.overload_sign != 0 )
3405         current_sign = internal_params.overload_sign < 0 ? -1 : 1;
3406 
3407       for ( i = 0; i < width; i++ )
3408       {
3409         index = (FT_UInt)( j * width + i );
3410 
3411         /* if the pixel is not set                     */
3412         /* its shortest distance is more than `spread` */
3413         if ( dists[index].sign == 0 )
3414           dists[index].distance = fixed_spread;
3415         else
3416           current_sign = dists[index].sign;
3417 
3418         /* clamp the values */
3419         if ( dists[index].distance > fixed_spread )
3420           dists[index].distance = fixed_spread;
3421 
3422         /* flip sign if required */
3423         dists[index].distance *= internal_params.flip_sign ? -current_sign
3424                                                            :  current_sign;
3425 
3426         /* concatenate to appropriate format */
3427         buffer[index] = map_fixed_to_sdf( dists[index].distance,
3428                                           fixed_spread );
3429       }
3430     }
3431 
3432   Exit:
3433     FT_FREE( dists );
3434     return error;
3435   }
3436 
3437 
3438   /**************************************************************************
3439    *
3440    * @Function:
3441    *   sdf_generate_subdivision
3442    *
3443    * @Description:
3444    *   Subdivide the shape into a number of straight lines, then use the
3445    *   above `sdf_generate_bounding_box` function to generate the SDF.
3446    *
3447    *   Note: After calling this function `shape` no longer has the original
3448    *         edges, it only contains lines.
3449    *
3450    * @Input:
3451    *   internal_params ::
3452    *     Internal parameters and properties required by the rasterizer.
3453    *     See @SDF_Params for more.
3454    *
3455    *   shape ::
3456    *     A complete shape which is used to generate SDF.
3457    *
3458    *   spread ::
3459    *     Maximum distances to be allowed inthe output bitmap.
3460    *
3461    * @Output:
3462    *   bitmap ::
3463    *     The output bitmap which will contain the SDF information.
3464    *
3465    * @Return:
3466    *   FreeType error, 0 means success.
3467    *
3468    */
3469   static FT_Error
sdf_generate_subdivision(const SDF_Params internal_params,SDF_Shape * shape,FT_UInt spread,const FT_Bitmap * bitmap)3470   sdf_generate_subdivision( const SDF_Params  internal_params,
3471                             SDF_Shape*        shape,
3472                             FT_UInt           spread,
3473                             const FT_Bitmap*  bitmap )
3474   {
3475     /*
3476      * Thanks to Alexei for providing the idea of this optimization.
3477      *
3478      * We take advantage of two facts.
3479      *
3480      * (1) Computing the shortest distance from a point to a line segment is
3481      *     very fast.
3482      * (2) We don't have to compute the shortest distance for the entire
3483      *     two-dimensional grid.
3484      *
3485      * Both ideas lead to the following optimization.
3486      *
3487      * (1) Split the outlines into a number of line segments.
3488      *
3489      * (2) For each line segment, only process its neighborhood.
3490      *
3491      * (3) Compute the closest distance to the line only for neighborhood
3492      *     grid points.
3493      *
3494      * This greatly reduces the number of grid points to check.
3495      */
3496 
3497     FT_Error  error = FT_Err_Ok;
3498 
3499 
3500     FT_CALL( split_sdf_shape( shape ) );
3501     FT_CALL( sdf_generate_bounding_box( internal_params,
3502                                         shape, spread, bitmap ) );
3503 
3504   Exit:
3505     return error;
3506   }
3507 
3508 
3509   /**************************************************************************
3510    *
3511    * @Function:
3512    *   sdf_generate_with_overlaps
3513    *
3514    * @Description:
3515    *   This function can be used to generate SDF for glyphs with overlapping
3516    *   contours.  The function generates SDF for contours separately on
3517    *   separate bitmaps (to generate SDF it uses
3518    *   `sdf_generate_subdivision`).  At the end it simply combines all the
3519    *   SDF into the output bitmap; this fixes all the signs and removes
3520    *   overlaps.
3521    *
3522    * @Input:
3523    *   internal_params ::
3524    *     Internal parameters and properties required by the rasterizer.  See
3525    *     @SDF_Params for more.
3526    *
3527    *   shape ::
3528    *     A complete shape which is used to generate SDF.
3529    *
3530    *   spread ::
3531    *     Maximum distances to be allowed in the output bitmap.
3532    *
3533    * @Output:
3534    *   bitmap ::
3535    *     The output bitmap which will contain the SDF information.
3536    *
3537    * @Return:
3538    *   FreeType error, 0 means success.
3539    *
3540    * @Note:
3541    *   The function cannot generate a proper SDF for glyphs with
3542    *   self-intersecting contours because we cannot separate them into two
3543    *   separate bitmaps.  In case of self-intersecting contours it is
3544    *   necessary to remove the overlaps before generating the SDF.
3545    *
3546    */
3547   static FT_Error
sdf_generate_with_overlaps(SDF_Params internal_params,SDF_Shape * shape,FT_UInt spread,const FT_Bitmap * bitmap)3548   sdf_generate_with_overlaps( SDF_Params        internal_params,
3549                               SDF_Shape*        shape,
3550                               FT_UInt           spread,
3551                               const FT_Bitmap*  bitmap )
3552   {
3553     FT_Error  error = FT_Err_Ok;
3554 
3555     FT_Int      num_contours;        /* total number of contours      */
3556     FT_Int      i, j;                /* iterators                     */
3557     FT_Int      width, rows;         /* width and rows of the bitmap  */
3558     FT_Bitmap*  bitmaps;             /* separate bitmaps for contours */
3559 
3560     SDF_Contour*  contour;           /* temporary variable to iterate */
3561     SDF_Contour*  temp_contour;      /* temporary contour             */
3562     SDF_Contour*  head;              /* head of the contour list      */
3563     SDF_Shape     temp_shape;        /* temporary shape               */
3564 
3565     FT_Memory      memory;           /* to allocate memory            */
3566     FT_SDFFormat*  t;                /* target bitmap buffer          */
3567     FT_Bool        flip_sign;        /* flip sign?                    */
3568 
3569     /* orientation of all the separate contours */
3570     SDF_Contour_Orientation*  orientations;
3571 
3572 
3573     bitmaps      = NULL;
3574     orientations = NULL;
3575     head         = NULL;
3576 
3577     if ( !shape || !bitmap || !shape->memory )
3578       return FT_THROW( Invalid_Argument );
3579 
3580     /* Disable `flip_sign` to avoid extra complication */
3581     /* during the combination phase.                   */
3582     flip_sign                 = internal_params.flip_sign;
3583     internal_params.flip_sign = 0;
3584 
3585     contour           = shape->contours;
3586     memory            = shape->memory;
3587     temp_shape.memory = memory;
3588     width             = (FT_Int)bitmap->width;
3589     rows              = (FT_Int)bitmap->rows;
3590     num_contours      = 0;
3591 
3592     /* find the number of contours in the shape */
3593     while ( contour )
3594     {
3595       num_contours++;
3596       contour = contour->next;
3597     }
3598 
3599     /* allocate the bitmaps to generate SDF for separate contours */
3600     if ( FT_ALLOC( bitmaps,
3601                    (FT_UInt)num_contours * sizeof ( *bitmaps ) ) )
3602       goto Exit;
3603 
3604     /* allocate array to hold orientation for all contours */
3605     if ( FT_ALLOC( orientations,
3606                    (FT_UInt)num_contours * sizeof ( *orientations ) ) )
3607       goto Exit;
3608 
3609     contour = shape->contours;
3610 
3611     /* Iterate over all contours and generate SDF separately. */
3612     for ( i = 0; i < num_contours; i++ )
3613     {
3614       /* initialize the corresponding bitmap */
3615       FT_Bitmap_Init( &bitmaps[i] );
3616 
3617       bitmaps[i].width      = bitmap->width;
3618       bitmaps[i].rows       = bitmap->rows;
3619       bitmaps[i].pitch      = bitmap->pitch;
3620       bitmaps[i].num_grays  = bitmap->num_grays;
3621       bitmaps[i].pixel_mode = bitmap->pixel_mode;
3622 
3623       /* allocate memory for the buffer */
3624       if ( FT_ALLOC( bitmaps[i].buffer,
3625                      bitmap->rows * (FT_UInt)bitmap->pitch ) )
3626         goto Exit;
3627 
3628       /* determine the orientation */
3629       orientations[i] = get_contour_orientation( contour );
3630 
3631       /* The `overload_sign` property is specific to  */
3632       /* `sdf_generate_bounding_box`.  This basically */
3633       /* overloads the default sign of the outside    */
3634       /* pixels, which is necessary for               */
3635       /* counter-clockwise contours.                  */
3636       if ( orientations[i] == SDF_ORIENTATION_CCW                   &&
3637            internal_params.orientation == FT_ORIENTATION_FILL_RIGHT )
3638         internal_params.overload_sign = 1;
3639       else if ( orientations[i] == SDF_ORIENTATION_CW                   &&
3640                 internal_params.orientation == FT_ORIENTATION_FILL_LEFT )
3641         internal_params.overload_sign = 1;
3642       else
3643         internal_params.overload_sign = 0;
3644 
3645       /* Make `contour->next` NULL so that there is   */
3646       /* one contour in the list.  Also hold the next */
3647       /* contour in a temporary variable so as to     */
3648       /* restore the original value.                  */
3649       temp_contour  = contour->next;
3650       contour->next = NULL;
3651 
3652       /* Use `temp_shape` to hold the new contour. */
3653       /* Now, `temp_shape` has only one contour.   */
3654       temp_shape.contours = contour;
3655 
3656       /* finally generate the SDF */
3657       FT_CALL( sdf_generate_subdivision( internal_params,
3658                                          &temp_shape,
3659                                          spread,
3660                                          &bitmaps[i] ) );
3661 
3662       /* Restore the original `next` variable. */
3663       contour->next = temp_contour;
3664 
3665       /* Since `split_sdf_shape` deallocated the original */
3666       /* contours list we need to assign the new value to */
3667       /* the shape's contour.                             */
3668       temp_shape.contours->next = head;
3669       head                      = temp_shape.contours;
3670 
3671       /* Simply flip the orientation in case of post-script fonts */
3672       /* so as to avoid modificatons in the combining phase.      */
3673       if ( internal_params.orientation == FT_ORIENTATION_FILL_LEFT )
3674       {
3675         if ( orientations[i] == SDF_ORIENTATION_CW )
3676           orientations[i] = SDF_ORIENTATION_CCW;
3677         else if ( orientations[i] == SDF_ORIENTATION_CCW )
3678           orientations[i] = SDF_ORIENTATION_CW;
3679       }
3680 
3681       contour = contour->next;
3682     }
3683 
3684     /* assign the new contour list to `shape->contours` */
3685     shape->contours = head;
3686 
3687     /* cast the output bitmap buffer */
3688     t = (FT_SDFFormat*)bitmap->buffer;
3689 
3690     /* Iterate over all pixels and combine all separate    */
3691     /* contours.  These are the rules for combining:       */
3692     /*                                                     */
3693     /* (1) For all clockwise contours, compute the largest */
3694     /*     value.  Name this as `val_c`.                   */
3695     /* (2) For all counter-clockwise contours, compute the */
3696     /*     smallest value.  Name this as `val_ac`.         */
3697     /* (3) Now, finally use the smaller value of `val_c'   */
3698     /*     and `val_ac'.                                   */
3699     for ( j = 0; j < rows; j++ )
3700     {
3701       for ( i = 0; i < width; i++ )
3702       {
3703         FT_Int  id = j * width + i;       /* index of current pixel    */
3704         FT_Int  c;                        /* contour iterator          */
3705 
3706         FT_SDFFormat  val_c  = 0;         /* max clockwise value       */
3707         FT_SDFFormat  val_ac = UCHAR_MAX; /* min counter-clockwise val */
3708 
3709 
3710         /* iterate through all the contours */
3711         for ( c = 0; c < num_contours; c++ )
3712         {
3713           /* current contour value */
3714           FT_SDFFormat  temp = ( (FT_SDFFormat*)bitmaps[c].buffer )[id];
3715 
3716 
3717           if ( orientations[c] == SDF_ORIENTATION_CW )
3718             val_c = FT_MAX( val_c, temp );   /* clockwise         */
3719           else
3720             val_ac = FT_MIN( val_ac, temp ); /* counter-clockwise */
3721         }
3722 
3723         /* Finally find the smaller of the two and assign to output. */
3724         /* Also apply `flip_sign` if set.                            */
3725         t[id] = FT_MIN( val_c, val_ac );
3726 
3727         if ( flip_sign )
3728           t[id] = invert_sign( t[id] );
3729       }
3730     }
3731 
3732   Exit:
3733     /* deallocate orientations array */
3734     if ( orientations )
3735       FT_FREE( orientations );
3736 
3737     /* deallocate temporary bitmaps */
3738     if ( bitmaps )
3739     {
3740       if ( num_contours == 0 )
3741         error = FT_THROW( Raster_Corrupted );
3742       else
3743       {
3744         for ( i = 0; i < num_contours; i++ )
3745           FT_FREE( bitmaps[i].buffer );
3746 
3747         FT_FREE( bitmaps );
3748       }
3749     }
3750 
3751     /* restore the `flip_sign` property */
3752     internal_params.flip_sign = flip_sign;
3753 
3754     return error;
3755   }
3756 
3757 
3758   /**************************************************************************
3759    *
3760    * interface functions
3761    *
3762    */
3763 
3764   static FT_Error
sdf_raster_new(void * memory_,FT_Raster * araster_)3765   sdf_raster_new( void*       memory_,   /* FT_Memory    */
3766                   FT_Raster*  araster_ ) /* SDF_PRaster* */
3767   {
3768     FT_Memory     memory  = (FT_Memory)memory_;
3769     SDF_PRaster*  araster = (SDF_PRaster*)araster_;
3770 
3771 
3772     FT_Error     error;
3773     SDF_PRaster  raster = NULL;
3774 
3775 
3776     if ( !FT_NEW( raster ) )
3777       raster->memory = memory;
3778 
3779     *araster = raster;
3780 
3781    return error;
3782   }
3783 
3784 
3785   static void
sdf_raster_reset(FT_Raster raster,unsigned char * pool_base,unsigned long pool_size)3786   sdf_raster_reset( FT_Raster       raster,
3787                     unsigned char*  pool_base,
3788                     unsigned long   pool_size )
3789   {
3790     FT_UNUSED( raster );
3791     FT_UNUSED( pool_base );
3792     FT_UNUSED( pool_size );
3793   }
3794 
3795 
3796   static FT_Error
sdf_raster_set_mode(FT_Raster raster,unsigned long mode,void * args)3797   sdf_raster_set_mode( FT_Raster      raster,
3798                        unsigned long  mode,
3799                        void*          args )
3800   {
3801     FT_UNUSED( raster );
3802     FT_UNUSED( mode );
3803     FT_UNUSED( args );
3804 
3805     return FT_Err_Ok;
3806   }
3807 
3808 
3809   static FT_Error
sdf_raster_render(FT_Raster raster,const FT_Raster_Params * params)3810   sdf_raster_render( FT_Raster                raster,
3811                      const FT_Raster_Params*  params )
3812   {
3813     FT_Error                  error      = FT_Err_Ok;
3814     SDF_TRaster*              sdf_raster = (SDF_TRaster*)raster;
3815     FT_Outline*               outline    = NULL;
3816     const SDF_Raster_Params*  sdf_params = (const SDF_Raster_Params*)params;
3817 
3818     FT_Memory   memory = NULL;
3819     SDF_Shape*  shape  = NULL;
3820     SDF_Params  internal_params;
3821 
3822 
3823     /* check for valid arguments */
3824     if ( !sdf_raster || !sdf_params )
3825     {
3826       error = FT_THROW( Invalid_Argument );
3827       goto Exit;
3828     }
3829 
3830     outline = (FT_Outline*)sdf_params->root.source;
3831 
3832     /* check whether outline is valid */
3833     if ( !outline )
3834     {
3835       error = FT_THROW( Invalid_Outline );
3836       goto Exit;
3837     }
3838 
3839     /* if the outline is empty, return */
3840     if ( outline->n_points <= 0 || outline->n_contours <= 0 )
3841       goto Exit;
3842 
3843     /* check whether the outline has valid fields */
3844     if ( !outline->contours || !outline->points )
3845     {
3846       error = FT_THROW( Invalid_Outline );
3847       goto Exit;
3848     }
3849 
3850     /* check whether spread is set properly */
3851     if ( sdf_params->spread > MAX_SPREAD ||
3852          sdf_params->spread < MIN_SPREAD )
3853     {
3854       FT_TRACE0(( "sdf_raster_render:"
3855                   " The `spread' field of `SDF_Raster_Params' is invalid,\n" ));
3856       FT_TRACE0(( "                  "
3857                   " the value of this field must be within [%d, %d].\n",
3858                   MIN_SPREAD, MAX_SPREAD ));
3859       FT_TRACE0(( "                  "
3860                   " Also, you must pass `SDF_Raster_Params' instead of\n" ));
3861       FT_TRACE0(( "                  "
3862                   " the default `FT_Raster_Params' while calling\n" ));
3863       FT_TRACE0(( "                  "
3864                   " this function and set the fields properly.\n" ));
3865 
3866       error = FT_THROW( Invalid_Argument );
3867       goto Exit;
3868     }
3869 
3870     memory = sdf_raster->memory;
3871     if ( !memory )
3872     {
3873       FT_TRACE0(( "sdf_raster_render:"
3874                   " Raster not setup properly,\n" ));
3875       FT_TRACE0(( "                  "
3876                   " unable to find memory handle.\n" ));
3877 
3878       error = FT_THROW( Invalid_Handle );
3879       goto Exit;
3880     }
3881 
3882     /* set up the parameters */
3883     internal_params.orientation   = FT_Outline_Get_Orientation( outline );
3884     internal_params.flip_sign     = sdf_params->flip_sign;
3885     internal_params.flip_y        = sdf_params->flip_y;
3886     internal_params.overload_sign = 0;
3887 
3888     FT_CALL( sdf_shape_new( memory, &shape ) );
3889 
3890     FT_CALL( sdf_outline_decompose( outline, shape ) );
3891 
3892     if ( sdf_params->overlaps )
3893       FT_CALL( sdf_generate_with_overlaps( internal_params,
3894                                            shape, sdf_params->spread,
3895                                            sdf_params->root.target ) );
3896     else
3897       FT_CALL( sdf_generate_subdivision( internal_params,
3898                                          shape, sdf_params->spread,
3899                                          sdf_params->root.target ) );
3900 
3901     if ( shape )
3902       sdf_shape_done( &shape );
3903 
3904   Exit:
3905     return error;
3906   }
3907 
3908 
3909   static void
sdf_raster_done(FT_Raster raster)3910   sdf_raster_done( FT_Raster  raster )
3911   {
3912     FT_Memory  memory = (FT_Memory)((SDF_TRaster*)raster)->memory;
3913 
3914 
3915     FT_FREE( raster );
3916   }
3917 
3918 
3919   FT_DEFINE_RASTER_FUNCS(
3920     ft_sdf_raster,
3921 
3922     FT_GLYPH_FORMAT_OUTLINE,
3923 
3924     (FT_Raster_New_Func)     sdf_raster_new,       /* raster_new      */
3925     (FT_Raster_Reset_Func)   sdf_raster_reset,     /* raster_reset    */
3926     (FT_Raster_Set_Mode_Func)sdf_raster_set_mode,  /* raster_set_mode */
3927     (FT_Raster_Render_Func)  sdf_raster_render,    /* raster_render   */
3928     (FT_Raster_Done_Func)    sdf_raster_done       /* raster_done     */
3929   )
3930 
3931 
3932 /* END */
3933