diff options
Diffstat (limited to 'src/core/SkCordic.cpp')
| -rw-r--r-- | src/core/SkCordic.cpp | 301 |
1 files changed, 301 insertions, 0 deletions
diff --git a/src/core/SkCordic.cpp b/src/core/SkCordic.cpp new file mode 100644 index 0000000..539bc9b --- /dev/null +++ b/src/core/SkCordic.cpp @@ -0,0 +1,301 @@ +/* libs/corecg/SkCordic.cpp +** +** Copyright 2006, The Android Open Source Project +** +** Licensed under the Apache License, Version 2.0 (the "License"); +** you may not use this file except in compliance with the License. +** You may obtain a copy of the License at +** +** http://www.apache.org/licenses/LICENSE-2.0 +** +** Unless required by applicable law or agreed to in writing, software +** distributed under the License is distributed on an "AS IS" BASIS, +** WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +** See the License for the specific language governing permissions and +** limitations under the License. +*/ + +#include "SkCordic.h" +#include "SkMath.h" +#include "Sk64.h" + +// 0x20000000 equals pi / 4 +const int32_t kATanDegrees[] = { 0x20000000, + 0x12E4051D, 0x9FB385B, 0x51111D4, 0x28B0D43, 0x145D7E1, 0xA2F61E, 0x517C55, + 0x28BE53, 0x145F2E, 0xA2F98, 0x517CC, 0x28BE6, 0x145F3, 0xA2F9, 0x517C, + 0x28BE, 0x145F, 0xA2F, 0x517, 0x28B, 0x145, 0xA2, 0x51, 0x28, 0x14, + 0xA, 0x5, 0x2, 0x1 }; + +const int32_t kFixedInvGain1 = 0x18bde0bb; // 0.607252935 + +static void SkCircularRotation(int32_t* x0, int32_t* y0, int32_t* z0) +{ + int32_t t = 0; + int32_t x = *x0; + int32_t y = *y0; + int32_t z = *z0; + const int32_t* tanPtr = kATanDegrees; + do { + int32_t x1 = y >> t; + int32_t y1 = x >> t; + int32_t tan = *tanPtr++; + if (z >= 0) { + x -= x1; + y += y1; + z -= tan; + } else { + x += x1; + y -= y1; + z += tan; + } + } while (++t < 16); // 30); + *x0 = x; + *y0 = y; + *z0 = z; +} + +SkFixed SkCordicSinCos(SkFixed radians, SkFixed* cosp) +{ + int32_t scaledRadians = radians * 0x28be; // scale radians to 65536 / PI() + int quadrant = scaledRadians >> 30; + quadrant += 1; + if (quadrant & 2) + scaledRadians = -scaledRadians + 0x80000000; + /* |a| <= 90 degrees as a 1.31 number */ + SkFixed sin = 0; + SkFixed cos = kFixedInvGain1; + SkCircularRotation(&cos, &sin, &scaledRadians); + Sk64 scaled; + scaled.setMul(sin, 0x6488d); + sin = scaled.fHi; + scaled.setMul(cos, 0x6488d); + if (quadrant & 2) + scaled.fHi = - scaled.fHi; + *cosp = scaled.fHi; + return sin; +} + +SkFixed SkCordicTan(SkFixed a) +{ + int32_t cos; + int32_t sin = SkCordicSinCos(a, &cos); + return SkFixedDiv(sin, cos); +} + +static int32_t SkCircularVector(int32_t* y0, int32_t* x0, int32_t vecMode) +{ + int32_t x = *x0; + int32_t y = *y0; + int32_t z = 0; + int32_t t = 0; + const int32_t* tanPtr = kATanDegrees; + do { + int32_t x1 = y >> t; + int32_t y1 = x >> t; + int32_t tan = *tanPtr++; + if (y < vecMode) { + x -= x1; + y += y1; + z -= tan; + } else { + x += x1; + y -= y1; + z += tan; + } + } while (++t < 16); // 30 + Sk64 scaled; + scaled.setMul(z, 0x6488d); // scale back into the SkScalar space (0x100000000/0x28be) + return scaled.fHi; +} + +SkFixed SkCordicASin(SkFixed a) { + int32_t sign = SkExtractSign(a); + int32_t z = SkFixedAbs(a); + if (z >= SK_Fixed1) + return SkApplySign(SK_FixedPI>>1, sign); + int32_t x = kFixedInvGain1; + int32_t y = 0; + z *= 0x28be; + z = SkCircularVector(&y, &x, z); + z = SkApplySign(z, ~sign); + return z; +} + +SkFixed SkCordicACos(SkFixed a) { + int32_t z = SkCordicASin(a); + z = (SK_FixedPI>>1) - z; + return z; +} + +SkFixed SkCordicATan2(SkFixed y, SkFixed x) { + if ((x | y) == 0) + return 0; + int32_t xsign = SkExtractSign(x); + x = SkFixedAbs(x); + int32_t result = SkCircularVector(&y, &x, 0); + if (xsign) { + int32_t rsign = SkExtractSign(result); + if (y == 0) + rsign = 0; + SkFixed pi = SkApplySign(SK_FixedPI, rsign); + result = pi - result; + } + return result; +} + +const int32_t kATanHDegrees[] = { + 0x1661788D, 0xA680D61, 0x51EA6FC, 0x28CBFDD, 0x1460E34, + 0xA2FCE8, 0x517D2E, 0x28BE6E, 0x145F32, + 0xA2F98, 0x517CC, 0x28BE6, 0x145F3, 0xA2F9, 0x517C, + 0x28BE, 0x145F, 0xA2F, 0x517, 0x28B, 0x145, 0xA2, 0x51, 0x28, 0x14, + 0xA, 0x5, 0x2, 0x1 }; + +const int32_t kFixedInvGain2 = 0x31330AAA; // 1.207534495 + +static void SkHyperbolic(int32_t* x0, int32_t* y0, int32_t* z0, int mode) +{ + int32_t t = 1; + int32_t x = *x0; + int32_t y = *y0; + int32_t z = *z0; + const int32_t* tanPtr = kATanHDegrees; + int k = -3; + do { + int32_t x1 = y >> t; + int32_t y1 = x >> t; + int32_t tan = *tanPtr++; + int count = 2 + (k >> 31); + if (++k == 1) + k = -2; + do { + if (((y >> 31) & mode) | ~((z >> 31) | mode)) { + x += x1; + y += y1; + z -= tan; + } else { + x -= x1; + y -= y1; + z += tan; + } + } while (--count); + } while (++t < 30); + *x0 = x; + *y0 = y; + *z0 = z; +} + +SkFixed SkCordicLog(SkFixed a) { + a *= 0x28be; + int32_t x = a + 0x28BE60DB; // 1.0 + int32_t y = a - 0x28BE60DB; + int32_t z = 0; + SkHyperbolic(&x, &y, &z, -1); + Sk64 scaled; + scaled.setMul(z, 0x6488d); + z = scaled.fHi; + return z << 1; +} + +SkFixed SkCordicExp(SkFixed a) { + int32_t cosh = kFixedInvGain2; + int32_t sinh = 0; + SkHyperbolic(&cosh, &sinh, &a, 0); + return cosh + sinh; +} + +#ifdef SK_DEBUG + +#ifdef SK_CAN_USE_FLOAT + #include "SkFloatingPoint.h" +#endif + +void SkCordic_UnitTest() +{ +#if defined(SK_SUPPORT_UNITTEST) && defined(SK_CAN_USE_FLOAT) + float val; + for (float angle = -720; angle < 720; angle += 30) { + float radian = angle * 3.1415925358f / 180.0f; + SkFixed f_angle = (int) (radian * 65536.0f); + // sincos + float sine = sinf(radian); + float cosine = cosf(radian); + SkFixed f_cosine; + SkFixed f_sine = SkCordicSinCos(f_angle, &f_cosine); + float sine2 = (float) f_sine / 65536.0f; + float cosine2 = (float) f_cosine / 65536.0f; + float error = fabsf(sine - sine2); + if (error > 0.001) + SkDebugf("sin error : angle = %g ; sin = %g ; cordic = %g\n", angle, sine, sine2); + error = fabsf(cosine - cosine2); + if (error > 0.001) + SkDebugf("cos error : angle = %g ; cos = %g ; cordic = %g\n", angle, cosine, cosine2); + // tan + float _tan = tanf(radian); + SkFixed f_tan = SkCordicTan(f_angle); + float tan2 = (float) f_tan / 65536.0f; + error = fabsf(_tan - tan2); + if (error > 0.05 && fabsf(_tan) < 1e6) + SkDebugf("tan error : angle = %g ; tan = %g ; cordic = %g\n", angle, _tan, tan2); + } + for (val = -1; val <= 1; val += .1f) { + SkFixed f_val = (int) (val * 65536.0f); + // asin + float arcsine = asinf(val); + SkFixed f_arcsine = SkCordicASin(f_val); + float arcsine2 = (float) f_arcsine / 65536.0f; + float error = fabsf(arcsine - arcsine2); + if (error > 0.001) + SkDebugf("asin error : val = %g ; asin = %g ; cordic = %g\n", val, arcsine, arcsine2); + } +#if 1 + for (val = -1; val <= 1; val += .1f) { +#else + val = .5; { +#endif + SkFixed f_val = (int) (val * 65536.0f); + // acos + float arccos = acosf(val); + SkFixed f_arccos = SkCordicACos(f_val); + float arccos2 = (float) f_arccos / 65536.0f; + float error = fabsf(arccos - arccos2); + if (error > 0.001) + SkDebugf("acos error : val = %g ; acos = %g ; cordic = %g\n", val, arccos, arccos2); + } + // atan2 +#if 1 + for (val = -1000; val <= 1000; val += 500.f) { + for (float val2 = -1000; val2 <= 1000; val2 += 500.f) { +#else + val = 0; { + float val2 = -1000; { +#endif + SkFixed f_val = (int) (val * 65536.0f); + SkFixed f_val2 = (int) (val2 * 65536.0f); + float arctan = atan2f(val, val2); + SkFixed f_arctan = SkCordicATan2(f_val, f_val2); + float arctan2 = (float) f_arctan / 65536.0f; + float error = fabsf(arctan - arctan2); + if (error > 0.001) + SkDebugf("atan2 error : val = %g ; val2 = %g ; atan2 = %g ; cordic = %g\n", val, val2, arctan, arctan2); + } + } + // log +#if 1 + for (val = 0.125f; val <= 8.f; val *= 2.0f) { +#else + val = .5; { +#endif + SkFixed f_val = (int) (val * 65536.0f); + // acos + float log = logf(val); + SkFixed f_log = SkCordicLog(f_val); + float log2 = (float) f_log / 65536.0f; + float error = fabsf(log - log2); + if (error > 0.001) + SkDebugf("log error : val = %g ; log = %g ; cordic = %g\n", val, log, log2); + } + // exp +#endif +} + +#endif |