images/math.h

   1 #ifndef _IMAGES_MATH_H
   2 #define _IMAGES_MATH_H
   3
   4 extern const u32 fast_div_tab[];
   5 extern const byte fast_sqrt_tab[];
   6
   7 static inline uns
   8 isqr(int x)
   9 {
  10   return x * x;
  11 }
  12
  13 static inline uns
  14 fast_div_u32_u8(uns x, uns y)
  15 {
  16 #ifdef CPU_I386
  17   int ret, dmy;
  18   asm volatile (
  19     "mull %3"
  20     :"=d"(ret),"=a"(dmy)
  21     :"1"(x),"g"(fast_div_tab[y])
  22   );
  23   return ret;
  24 #else
  25   return ((u64)(x) * fast_div_tab[y]) >> 32;
  26 #endif
  27 }
  28
  29 static inline uns
  30 fast_sqrt_u16(uns x)
  31 {
  32   uns y;
  33   if (x < (1 << 10) - 3)
  34     y = fast_sqrt_tab[(x + 3) >> 2] >> 3;
  35   else if (x < (1 << 14) - 28)
  36     y = fast_sqrt_tab[(x + 28) >> 6] >> 1;
  37   else
  38     y = fast_sqrt_tab[x >> 8];
  39   return (x < y * y) ? y - 1 : y;
  40 }
  41
  42 static inline uns
  43 fast_sqrt_u32(uns x)
  44 {
  45   uns y;
  46   if (x < (1 << 16))
  47     {
  48       if (x < (1 << 10) - 3)
  49         y = fast_sqrt_tab[(x + 3) >> 2] >> 3;
  50       else if (x < (1 << 14) - 28)
  51         y = fast_sqrt_tab[(x + 28) >> 6] >> 1;
  52       else
  53         y = fast_sqrt_tab[x >> 8];
  54     }
  55   else
  56     {
  57       if (x < (1 << 24))
  58         {
  59           if (x < (1 << 20))
  60             {
  61               y = fast_sqrt_tab[x >> 12];
  62               y = (fast_div_u32_u8(x, y) >> 3) + (y << 1);
  63             }
  64           else
  65             {
  66               y = fast_sqrt_tab[x >> 16];
  67               y = (fast_div_u32_u8(x, y) >> 5) + (y << 3);
  68             }
  69         }
  70       else
  71         {
  72           if (x < (1 << 28))
  73             {
  74               if (x < (1 << 26))
  75                 {
  76                   y = fast_sqrt_tab[x >> 18];
  77                   y = (fast_div_u32_u8(x, y) >> 6) + (y << 4);
  78                 }
  79               else
  80                 {
  81                   y = fast_sqrt_tab[x >> 20];
  82                   y = (fast_div_u32_u8(x, y) >> 7) + (y << 5);
  83                 }
  84             }
  85           else
  86             {
  87               if (x < (1 << 30))
  88                 {
  89                   y = fast_sqrt_tab[x >> 22];
  90                   y = (fast_div_u32_u8(x, y) >> 8) + (y << 6);
  91                 }
  92               else
  93                 {
  94                   y = fast_sqrt_tab[x >> 24];
  95                   y = (fast_div_u32_u8(x, y) >> 9) + (y << 7);
  96                 }
  97             }
  98         }
  99     }
 100   return (x < y * y) ? y - 1 : y;
 101 }
 102
 103 #endif