Home | History | Annotate | Download | only in base
      1 /***************************************************************************/
      2 /*                                                                         */
      3 /*  fttrigon.c                                                             */
      4 /*                                                                         */
      5 /*    FreeType trigonometric functions (body).                             */
      6 /*                                                                         */
      7 /*  Copyright 2001-2015 by                                                 */
      8 /*  David Turner, Robert Wilhelm, and Werner Lemberg.                      */
      9 /*                                                                         */
     10 /*  This file is part of the FreeType project, and may only be used,       */
     11 /*  modified, and distributed under the terms of the FreeType project      */
     12 /*  license, LICENSE.TXT.  By continuing to use, modify, or distribute     */
     13 /*  this file you indicate that you have read the license and              */
     14 /*  understand and accept it fully.                                        */
     15 /*                                                                         */
     16 /***************************************************************************/
     17 
     18   /*************************************************************************/
     19   /*                                                                       */
     20   /* This is a fixed-point CORDIC implementation of trigonometric          */
     21   /* functions as well as transformations between Cartesian and polar      */
     22   /* coordinates.  The angles are represented as 16.16 fixed-point values  */
     23   /* in degrees, i.e., the angular resolution is 2^-16 degrees.  Note that */
     24   /* only vectors longer than 2^16*180/pi (or at least 22 bits) on a       */
     25   /* discrete Cartesian grid can have the same or better angular           */
     26   /* resolution.  Therefore, to maintain this precision, some functions    */
     27   /* require an interim upscaling of the vectors, whereas others operate   */
     28   /* with 24-bit long vectors directly.                                    */
     29   /*                                                                       */
     30   /*************************************************************************/
     31 
     32 #include <ft2build.h>
     33 #include FT_INTERNAL_OBJECTS_H
     34 #include FT_INTERNAL_CALC_H
     35 #include FT_TRIGONOMETRY_H
     36 
     37 
     38   /* the Cordic shrink factor 0.858785336480436 * 2^32 */
     39 #define FT_TRIG_SCALE      0xDBD95B16UL
     40 
     41   /* the highest bit in overflow-safe vector components, */
     42   /* MSB of 0.858785336480436 * sqrt(0.5) * 2^30         */
     43 #define FT_TRIG_SAFE_MSB   29
     44 
     45   /* this table was generated for FT_PI = 180L << 16, i.e. degrees */
     46 #define FT_TRIG_MAX_ITERS  23
     47 
     48   static const FT_Angle
     49   ft_trig_arctan_table[] =
     50   {
     51     1740967L, 919879L, 466945L, 234379L, 117304L, 58666L, 29335L,
     52     14668L, 7334L, 3667L, 1833L, 917L, 458L, 229L, 115L,
     53     57L, 29L, 14L, 7L, 4L, 2L, 1L
     54   };
     55 
     56 
     57 #ifdef FT_LONG64
     58 
     59   /* multiply a given value by the CORDIC shrink factor */
     60   static FT_Fixed
     61   ft_trig_downscale( FT_Fixed  val )
     62   {
     63     FT_Int  s = 1;
     64 
     65 
     66     if ( val < 0 )
     67     {
     68        val = -val;
     69        s = -1;
     70     }
     71 
     72     /* 0x40000000 comes from regression analysis between true */
     73     /* and CORDIC hypotenuse, so it minimizes the error       */
     74     val = (FT_Fixed)( ( (FT_Int64)val * FT_TRIG_SCALE + 0x40000000UL ) >> 32 );
     75 
     76     return s < 0 ? -val : val;
     77   }
     78 
     79 #else /* !FT_LONG64 */
     80 
     81   /* multiply a given value by the CORDIC shrink factor */
     82   static FT_Fixed
     83   ft_trig_downscale( FT_Fixed  val )
     84   {
     85     FT_Int     s = 1;
     86     FT_UInt32  lo1, hi1, lo2, hi2, lo, hi, i1, i2;
     87 
     88 
     89     if ( val < 0 )
     90     {
     91        val = -val;
     92        s = -1;
     93     }
     94 
     95     lo1 = (FT_UInt32)val & 0x0000FFFFU;
     96     hi1 = (FT_UInt32)val >> 16;
     97     lo2 = FT_TRIG_SCALE & 0x0000FFFFU;
     98     hi2 = FT_TRIG_SCALE >> 16;
     99 
    100     lo = lo1 * lo2;
    101     i1 = lo1 * hi2;
    102     i2 = lo2 * hi1;
    103     hi = hi1 * hi2;
    104 
    105     /* Check carry overflow of i1 + i2 */
    106     i1 += i2;
    107     hi += (FT_UInt32)( i1 < i2 ) << 16;
    108 
    109     hi += i1 >> 16;
    110     i1  = i1 << 16;
    111 
    112     /* Check carry overflow of i1 + lo */
    113     lo += i1;
    114     hi += ( lo < i1 );
    115 
    116     /* 0x40000000 comes from regression analysis between true */
    117     /* and CORDIC hypotenuse, so it minimizes the error       */
    118 
    119     /* Check carry overflow of lo + 0x40000000 */
    120     lo += 0x40000000UL;
    121     hi += ( lo < 0x40000000UL );
    122 
    123     val = (FT_Fixed)hi;
    124 
    125     return s < 0 ? -val : val;
    126   }
    127 
    128 #endif /* !FT_LONG64 */
    129 
    130 
    131   /* undefined and never called for zero vector */
    132   static FT_Int
    133   ft_trig_prenorm( FT_Vector*  vec )
    134   {
    135     FT_Pos  x, y;
    136     FT_Int  shift;
    137 
    138 
    139     x = vec->x;
    140     y = vec->y;
    141 
    142     shift = FT_MSB( (FT_UInt32)( FT_ABS( x ) | FT_ABS( y ) ) );
    143 
    144     if ( shift <= FT_TRIG_SAFE_MSB )
    145     {
    146       shift  = FT_TRIG_SAFE_MSB - shift;
    147       vec->x = (FT_Pos)( (FT_ULong)x << shift );
    148       vec->y = (FT_Pos)( (FT_ULong)y << shift );
    149     }
    150     else
    151     {
    152       shift -= FT_TRIG_SAFE_MSB;
    153       vec->x = x >> shift;
    154       vec->y = y >> shift;
    155       shift  = -shift;
    156     }
    157 
    158     return shift;
    159   }
    160 
    161 
    162   static void
    163   ft_trig_pseudo_rotate( FT_Vector*  vec,
    164                          FT_Angle    theta )
    165   {
    166     FT_Int           i;
    167     FT_Fixed         x, y, xtemp, b;
    168     const FT_Angle  *arctanptr;
    169 
    170 
    171     x = vec->x;
    172     y = vec->y;
    173 
    174     /* Rotate inside [-PI/4,PI/4] sector */
    175     while ( theta < -FT_ANGLE_PI4 )
    176     {
    177       xtemp  =  y;
    178       y      = -x;
    179       x      =  xtemp;
    180       theta +=  FT_ANGLE_PI2;
    181     }
    182 
    183     while ( theta > FT_ANGLE_PI4 )
    184     {
    185       xtemp  = -y;
    186       y      =  x;
    187       x      =  xtemp;
    188       theta -=  FT_ANGLE_PI2;
    189     }
    190 
    191     arctanptr = ft_trig_arctan_table;
    192 
    193     /* Pseudorotations, with right shifts */
    194     for ( i = 1, b = 1; i < FT_TRIG_MAX_ITERS; b <<= 1, i++ )
    195     {
    196       if ( theta < 0 )
    197       {
    198         xtemp  = x + ( ( y + b ) >> i );
    199         y      = y - ( ( x + b ) >> i );
    200         x      = xtemp;
    201         theta += *arctanptr++;
    202       }
    203       else
    204       {
    205         xtemp  = x - ( ( y + b ) >> i );
    206         y      = y + ( ( x + b ) >> i );
    207         x      = xtemp;
    208         theta -= *arctanptr++;
    209       }
    210     }
    211 
    212     vec->x = x;
    213     vec->y = y;
    214   }
    215 
    216 
    217   static void
    218   ft_trig_pseudo_polarize( FT_Vector*  vec )
    219   {
    220     FT_Angle         theta;
    221     FT_Int           i;
    222     FT_Fixed         x, y, xtemp, b;
    223     const FT_Angle  *arctanptr;
    224 
    225 
    226     x = vec->x;
    227     y = vec->y;
    228 
    229     /* Get the vector into [-PI/4,PI/4] sector */
    230     if ( y > x )
    231     {
    232       if ( y > -x )
    233       {
    234         theta =  FT_ANGLE_PI2;
    235         xtemp =  y;
    236         y     = -x;
    237         x     =  xtemp;
    238       }
    239       else
    240       {
    241         theta =  y > 0 ? FT_ANGLE_PI : -FT_ANGLE_PI;
    242         x     = -x;
    243         y     = -y;
    244       }
    245     }
    246     else
    247     {
    248       if ( y < -x )
    249       {
    250         theta = -FT_ANGLE_PI2;
    251         xtemp = -y;
    252         y     =  x;
    253         x     =  xtemp;
    254       }
    255       else
    256       {
    257         theta = 0;
    258       }
    259     }
    260 
    261     arctanptr = ft_trig_arctan_table;
    262 
    263     /* Pseudorotations, with right shifts */
    264     for ( i = 1, b = 1; i < FT_TRIG_MAX_ITERS; b <<= 1, i++ )
    265     {
    266       if ( y > 0 )
    267       {
    268         xtemp  = x + ( ( y + b ) >> i );
    269         y      = y - ( ( x + b ) >> i );
    270         x      = xtemp;
    271         theta += *arctanptr++;
    272       }
    273       else
    274       {
    275         xtemp  = x - ( ( y + b ) >> i );
    276         y      = y + ( ( x + b ) >> i );
    277         x      = xtemp;
    278         theta -= *arctanptr++;
    279       }
    280     }
    281 
    282     /* round theta to acknowledge its error that mostly comes */
    283     /* from accumulated rounding errors in the arctan table   */
    284     if ( theta >= 0 )
    285       theta = FT_PAD_ROUND( theta, 16 );
    286     else
    287       theta = -FT_PAD_ROUND( -theta, 16 );
    288 
    289     vec->x = x;
    290     vec->y = theta;
    291   }
    292 
    293 
    294   /* documentation is in fttrigon.h */
    295 
    296   FT_EXPORT_DEF( FT_Fixed )
    297   FT_Cos( FT_Angle  angle )
    298   {
    299     FT_Vector  v;
    300 
    301 
    302     FT_Vector_Unit( &v, angle );
    303 
    304     return v.x;
    305   }
    306 
    307 
    308   /* documentation is in fttrigon.h */
    309 
    310   FT_EXPORT_DEF( FT_Fixed )
    311   FT_Sin( FT_Angle  angle )
    312   {
    313     FT_Vector  v;
    314 
    315 
    316     FT_Vector_Unit( &v, angle );
    317 
    318     return v.y;
    319   }
    320 
    321 
    322   /* documentation is in fttrigon.h */
    323 
    324   FT_EXPORT_DEF( FT_Fixed )
    325   FT_Tan( FT_Angle  angle )
    326   {
    327     FT_Vector  v;
    328 
    329 
    330     FT_Vector_Unit( &v, angle );
    331 
    332     return FT_DivFix( v.y, v.x );
    333   }
    334 
    335 
    336   /* documentation is in fttrigon.h */
    337 
    338   FT_EXPORT_DEF( FT_Angle )
    339   FT_Atan2( FT_Fixed  dx,
    340             FT_Fixed  dy )
    341   {
    342     FT_Vector  v;
    343 
    344 
    345     if ( dx == 0 && dy == 0 )
    346       return 0;
    347 
    348     v.x = dx;
    349     v.y = dy;
    350     ft_trig_prenorm( &v );
    351     ft_trig_pseudo_polarize( &v );
    352 
    353     return v.y;
    354   }
    355 
    356 
    357   /* documentation is in fttrigon.h */
    358 
    359   FT_EXPORT_DEF( void )
    360   FT_Vector_Unit( FT_Vector*  vec,
    361                   FT_Angle    angle )
    362   {
    363     if ( !vec )
    364       return;
    365 
    366     vec->x = FT_TRIG_SCALE >> 8;
    367     vec->y = 0;
    368     ft_trig_pseudo_rotate( vec, angle );
    369     vec->x = ( vec->x + 0x80L ) >> 8;
    370     vec->y = ( vec->y + 0x80L ) >> 8;
    371   }
    372 
    373 
    374   /* these macros return 0 for positive numbers,
    375      and -1 for negative ones */
    376 #define FT_SIGN_LONG( x )   ( (x) >> ( FT_SIZEOF_LONG * 8 - 1 ) )
    377 #define FT_SIGN_INT( x )    ( (x) >> ( FT_SIZEOF_INT * 8 - 1 ) )
    378 #define FT_SIGN_INT32( x )  ( (x) >> 31 )
    379 #define FT_SIGN_INT16( x )  ( (x) >> 15 )
    380 
    381 
    382   /* documentation is in fttrigon.h */
    383 
    384   FT_EXPORT_DEF( void )
    385   FT_Vector_Rotate( FT_Vector*  vec,
    386                     FT_Angle    angle )
    387   {
    388     FT_Int     shift;
    389     FT_Vector  v;
    390 
    391 
    392     if ( !vec || !angle )
    393       return;
    394 
    395     v = *vec;
    396 
    397     if ( v.x == 0 && v.y == 0 )
    398       return;
    399 
    400     shift = ft_trig_prenorm( &v );
    401     ft_trig_pseudo_rotate( &v, angle );
    402     v.x = ft_trig_downscale( v.x );
    403     v.y = ft_trig_downscale( v.y );
    404 
    405     if ( shift > 0 )
    406     {
    407       FT_Int32  half = (FT_Int32)1L << ( shift - 1 );
    408 
    409 
    410       vec->x = ( v.x + half + FT_SIGN_LONG( v.x ) ) >> shift;
    411       vec->y = ( v.y + half + FT_SIGN_LONG( v.y ) ) >> shift;
    412     }
    413     else
    414     {
    415       shift  = -shift;
    416       vec->x = (FT_Pos)( (FT_ULong)v.x << shift );
    417       vec->y = (FT_Pos)( (FT_ULong)v.y << shift );
    418     }
    419   }
    420 
    421 
    422   /* documentation is in fttrigon.h */
    423 
    424   FT_EXPORT_DEF( FT_Fixed )
    425   FT_Vector_Length( FT_Vector*  vec )
    426   {
    427     FT_Int     shift;
    428     FT_Vector  v;
    429 
    430 
    431     if ( !vec )
    432       return 0;
    433 
    434     v = *vec;
    435 
    436     /* handle trivial cases */
    437     if ( v.x == 0 )
    438     {
    439       return FT_ABS( v.y );
    440     }
    441     else if ( v.y == 0 )
    442     {
    443       return FT_ABS( v.x );
    444     }
    445 
    446     /* general case */
    447     shift = ft_trig_prenorm( &v );
    448     ft_trig_pseudo_polarize( &v );
    449 
    450     v.x = ft_trig_downscale( v.x );
    451 
    452     if ( shift > 0 )
    453       return ( v.x + ( 1L << ( shift - 1 ) ) ) >> shift;
    454 
    455     return (FT_Fixed)( (FT_UInt32)v.x << -shift );
    456   }
    457 
    458 
    459   /* documentation is in fttrigon.h */
    460 
    461   FT_EXPORT_DEF( void )
    462   FT_Vector_Polarize( FT_Vector*  vec,
    463                       FT_Fixed   *length,
    464                       FT_Angle   *angle )
    465   {
    466     FT_Int     shift;
    467     FT_Vector  v;
    468 
    469 
    470     if ( !vec || !length || !angle )
    471       return;
    472 
    473     v = *vec;
    474 
    475     if ( v.x == 0 && v.y == 0 )
    476       return;
    477 
    478     shift = ft_trig_prenorm( &v );
    479     ft_trig_pseudo_polarize( &v );
    480 
    481     v.x = ft_trig_downscale( v.x );
    482 
    483     *length = shift >= 0 ?                      ( v.x >>  shift )
    484                          : (FT_Fixed)( (FT_UInt32)v.x << -shift );
    485     *angle  = v.y;
    486   }
    487 
    488 
    489   /* documentation is in fttrigon.h */
    490 
    491   FT_EXPORT_DEF( void )
    492   FT_Vector_From_Polar( FT_Vector*  vec,
    493                         FT_Fixed    length,
    494                         FT_Angle    angle )
    495   {
    496     if ( !vec )
    497       return;
    498 
    499     vec->x = length;
    500     vec->y = 0;
    501 
    502     FT_Vector_Rotate( vec, angle );
    503   }
    504 
    505 
    506   /* documentation is in fttrigon.h */
    507 
    508   FT_EXPORT_DEF( FT_Angle )
    509   FT_Angle_Diff( FT_Angle  angle1,
    510                  FT_Angle  angle2 )
    511   {
    512     FT_Angle  delta = angle2 - angle1;
    513 
    514 
    515     while ( delta <= -FT_ANGLE_PI )
    516       delta += FT_ANGLE_2PI;
    517 
    518     while ( delta > FT_ANGLE_PI )
    519       delta -= FT_ANGLE_2PI;
    520 
    521     return delta;
    522   }
    523 
    524 
    525 /* END */
    526