1/*
2 * Copyright 2012 Google Inc.
3 *
4 * Use of this source code is governed by a BSD-style license that can be
5 * found in the LICENSE file.
6 */
7#include "SkArenaAlloc.h"
8#include "SkFloatBits.h"
9#include "SkOpCoincidence.h"
10#include "SkPathOpsTypes.h"
11
12static bool arguments_denormalized(float a, float b, int epsilon) {
13    float denormalizedCheck = FLT_EPSILON * epsilon / 2;
14    return fabsf(a) <= denormalizedCheck && fabsf(b) <= denormalizedCheck;
15}
16
17// from http://randomascii.wordpress.com/2012/02/25/comparing-floating-point-numbers-2012-edition/
18// FIXME: move to SkFloatBits.h
19static bool equal_ulps(float a, float b, int epsilon, int depsilon) {
20    if (arguments_denormalized(a, b, depsilon)) {
21        return true;
22    }
23    int aBits = SkFloatAs2sCompliment(a);
24    int bBits = SkFloatAs2sCompliment(b);
25    // Find the difference in ULPs.
26    return aBits < bBits + epsilon && bBits < aBits + epsilon;
27}
28
29static bool equal_ulps_no_normal_check(float a, float b, int epsilon, int depsilon) {
30    int aBits = SkFloatAs2sCompliment(a);
31    int bBits = SkFloatAs2sCompliment(b);
32    // Find the difference in ULPs.
33    return aBits < bBits + epsilon && bBits < aBits + epsilon;
34}
35
36static bool equal_ulps_pin(float a, float b, int epsilon, int depsilon) {
37    if (!SkScalarIsFinite(a) || !SkScalarIsFinite(b)) {
38        return false;
39    }
40    if (arguments_denormalized(a, b, depsilon)) {
41        return true;
42    }
43    int aBits = SkFloatAs2sCompliment(a);
44    int bBits = SkFloatAs2sCompliment(b);
45    // Find the difference in ULPs.
46    return aBits < bBits + epsilon && bBits < aBits + epsilon;
47}
48
49static bool d_equal_ulps(float a, float b, int epsilon) {
50    int aBits = SkFloatAs2sCompliment(a);
51    int bBits = SkFloatAs2sCompliment(b);
52    // Find the difference in ULPs.
53    return aBits < bBits + epsilon && bBits < aBits + epsilon;
54}
55
56static bool not_equal_ulps(float a, float b, int epsilon) {
57    if (arguments_denormalized(a, b, epsilon)) {
58        return false;
59    }
60    int aBits = SkFloatAs2sCompliment(a);
61    int bBits = SkFloatAs2sCompliment(b);
62    // Find the difference in ULPs.
63    return aBits >= bBits + epsilon || bBits >= aBits + epsilon;
64}
65
66static bool not_equal_ulps_pin(float a, float b, int epsilon) {
67    if (!SkScalarIsFinite(a) || !SkScalarIsFinite(b)) {
68        return false;
69    }
70    if (arguments_denormalized(a, b, epsilon)) {
71        return false;
72    }
73    int aBits = SkFloatAs2sCompliment(a);
74    int bBits = SkFloatAs2sCompliment(b);
75    // Find the difference in ULPs.
76    return aBits >= bBits + epsilon || bBits >= aBits + epsilon;
77}
78
79static bool d_not_equal_ulps(float a, float b, int epsilon) {
80    int aBits = SkFloatAs2sCompliment(a);
81    int bBits = SkFloatAs2sCompliment(b);
82    // Find the difference in ULPs.
83    return aBits >= bBits + epsilon || bBits >= aBits + epsilon;
84}
85
86static bool less_ulps(float a, float b, int epsilon) {
87    if (arguments_denormalized(a, b, epsilon)) {
88        return a <= b - FLT_EPSILON * epsilon;
89    }
90    int aBits = SkFloatAs2sCompliment(a);
91    int bBits = SkFloatAs2sCompliment(b);
92    // Find the difference in ULPs.
93    return aBits <= bBits - epsilon;
94}
95
96static bool less_or_equal_ulps(float a, float b, int epsilon) {
97    if (arguments_denormalized(a, b, epsilon)) {
98        return a < b + FLT_EPSILON * epsilon;
99    }
100    int aBits = SkFloatAs2sCompliment(a);
101    int bBits = SkFloatAs2sCompliment(b);
102    // Find the difference in ULPs.
103    return aBits < bBits + epsilon;
104}
105
106// equality using the same error term as between
107bool AlmostBequalUlps(float a, float b) {
108    const int UlpsEpsilon = 2;
109    return equal_ulps(a, b, UlpsEpsilon, UlpsEpsilon);
110}
111
112bool AlmostPequalUlps(float a, float b) {
113    const int UlpsEpsilon = 8;
114    return equal_ulps(a, b, UlpsEpsilon, UlpsEpsilon);
115}
116
117bool AlmostDequalUlps(float a, float b) {
118    const int UlpsEpsilon = 16;
119    return d_equal_ulps(a, b, UlpsEpsilon);
120}
121
122bool AlmostDequalUlps(double a, double b) {
123    if (fabs(a) < SK_ScalarMax && fabs(b) < SK_ScalarMax) {
124        return AlmostDequalUlps(SkDoubleToScalar(a), SkDoubleToScalar(b));
125    }
126    return fabs(a - b) / SkTMax(fabs(a), fabs(b)) < FLT_EPSILON * 16;
127}
128
129bool AlmostEqualUlps(float a, float b) {
130    const int UlpsEpsilon = 16;
131    return equal_ulps(a, b, UlpsEpsilon, UlpsEpsilon);
132}
133
134bool AlmostEqualUlpsNoNormalCheck(float a, float b) {
135    const int UlpsEpsilon = 16;
136    return equal_ulps_no_normal_check(a, b, UlpsEpsilon, UlpsEpsilon);
137}
138
139bool AlmostEqualUlps_Pin(float a, float b) {
140    const int UlpsEpsilon = 16;
141    return equal_ulps_pin(a, b, UlpsEpsilon, UlpsEpsilon);
142}
143
144bool NotAlmostEqualUlps(float a, float b) {
145    const int UlpsEpsilon = 16;
146    return not_equal_ulps(a, b, UlpsEpsilon);
147}
148
149bool NotAlmostEqualUlps_Pin(float a, float b) {
150    const int UlpsEpsilon = 16;
151    return not_equal_ulps_pin(a, b, UlpsEpsilon);
152}
153
154bool NotAlmostDequalUlps(float a, float b) {
155    const int UlpsEpsilon = 16;
156    return d_not_equal_ulps(a, b, UlpsEpsilon);
157}
158
159bool RoughlyEqualUlps(float a, float b) {
160    const int UlpsEpsilon = 256;
161    const int DUlpsEpsilon = 1024;
162    return equal_ulps(a, b, UlpsEpsilon, DUlpsEpsilon);
163}
164
165bool AlmostBetweenUlps(float a, float b, float c) {
166    const int UlpsEpsilon = 2;
167    return a <= c ? less_or_equal_ulps(a, b, UlpsEpsilon) && less_or_equal_ulps(b, c, UlpsEpsilon)
168        : less_or_equal_ulps(b, a, UlpsEpsilon) && less_or_equal_ulps(c, b, UlpsEpsilon);
169}
170
171bool AlmostLessUlps(float a, float b) {
172    const int UlpsEpsilon = 16;
173    return less_ulps(a, b, UlpsEpsilon);
174}
175
176bool AlmostLessOrEqualUlps(float a, float b) {
177    const int UlpsEpsilon = 16;
178    return less_or_equal_ulps(a, b, UlpsEpsilon);
179}
180
181int UlpsDistance(float a, float b) {
182    SkFloatIntUnion floatIntA, floatIntB;
183    floatIntA.fFloat = a;
184    floatIntB.fFloat = b;
185    // Different signs means they do not match.
186    if ((floatIntA.fSignBitInt < 0) != (floatIntB.fSignBitInt < 0)) {
187        // Check for equality to make sure +0 == -0
188        return a == b ? 0 : SK_MaxS32;
189    }
190    // Find the difference in ULPs.
191    return SkTAbs(floatIntA.fSignBitInt - floatIntB.fSignBitInt);
192}
193
194// cube root approximation using bit hack for 64-bit float
195// adapted from Kahan's cbrt
196static double cbrt_5d(double d) {
197    const unsigned int B1 = 715094163;
198    double t = 0.0;
199    unsigned int* pt = (unsigned int*) &t;
200    unsigned int* px = (unsigned int*) &d;
201    pt[1] = px[1] / 3 + B1;
202    return t;
203}
204
205// iterative cube root approximation using Halley's method (double)
206static double cbrta_halleyd(const double a, const double R) {
207    const double a3 = a * a * a;
208    const double b = a * (a3 + R + R) / (a3 + a3 + R);
209    return b;
210}
211
212// cube root approximation using 3 iterations of Halley's method (double)
213static double halley_cbrt3d(double d) {
214    double a = cbrt_5d(d);
215    a = cbrta_halleyd(a, d);
216    a = cbrta_halleyd(a, d);
217    return cbrta_halleyd(a, d);
218}
219
220double SkDCubeRoot(double x) {
221    if (approximately_zero_cubed(x)) {
222        return 0;
223    }
224    double result = halley_cbrt3d(fabs(x));
225    if (x < 0) {
226        result = -result;
227    }
228    return result;
229}
230
231SkOpGlobalState::SkOpGlobalState(SkOpContourHead* head,
232                                 SkArenaAlloc* allocator
233                                 SkDEBUGPARAMS(bool debugSkipAssert)
234                                 SkDEBUGPARAMS(const char* testName))
235    : fAllocator(allocator)
236    , fCoincidence(nullptr)
237    , fContourHead(head)
238    , fNested(0)
239    , fWindingFailed(false)
240    , fPhase(SkOpPhase::kIntersecting)
241    SkDEBUGPARAMS(fDebugTestName(testName))
242    SkDEBUGPARAMS(fAngleID(0))
243    SkDEBUGPARAMS(fCoinID(0))
244    SkDEBUGPARAMS(fContourID(0))
245    SkDEBUGPARAMS(fPtTID(0))
246    SkDEBUGPARAMS(fSegmentID(0))
247    SkDEBUGPARAMS(fSpanID(0))
248    SkDEBUGPARAMS(fDebugSkipAssert(debugSkipAssert)) {
249#if DEBUG_T_SECT_LOOP_COUNT
250    debugResetLoopCounts();
251#endif
252#if DEBUG_COIN
253    fPreviousFuncName = nullptr;
254#endif
255}
256