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 "include/private/SkFloatBits.h"
8 #include "src/core/SkArenaAlloc.h"
9 #include "src/pathops/SkOpCoincidence.h"
10 #include "src/pathops/SkPathOpsTypes.h"
11 
arguments_denormalized(float a,float b,int epsilon)12 static 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
equal_ulps(float a,float b,int epsilon,int depsilon)19 static 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 
equal_ulps_no_normal_check(float a,float b,int epsilon,int depsilon)29 static 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 
equal_ulps_pin(float a,float b,int epsilon,int depsilon)36 static 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 
d_equal_ulps(float a,float b,int epsilon)49 static 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 
not_equal_ulps(float a,float b,int epsilon)56 static 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 
not_equal_ulps_pin(float a,float b,int epsilon)66 static 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 
d_not_equal_ulps(float a,float b,int epsilon)79 static 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 
less_ulps(float a,float b,int epsilon)86 static 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 
less_or_equal_ulps(float a,float b,int epsilon)96 static 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
AlmostBequalUlps(float a,float b)107 bool AlmostBequalUlps(float a, float b) {
108     const int UlpsEpsilon = 2;
109     return equal_ulps(a, b, UlpsEpsilon, UlpsEpsilon);
110 }
111 
AlmostPequalUlps(float a,float b)112 bool AlmostPequalUlps(float a, float b) {
113     const int UlpsEpsilon = 8;
114     return equal_ulps(a, b, UlpsEpsilon, UlpsEpsilon);
115 }
116 
AlmostDequalUlps(float a,float b)117 bool AlmostDequalUlps(float a, float b) {
118     const int UlpsEpsilon = 16;
119     return d_equal_ulps(a, b, UlpsEpsilon);
120 }
121 
AlmostDequalUlps(double a,double b)122 bool 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) / std::max(fabs(a), fabs(b)) < FLT_EPSILON * 16;
127 }
128 
AlmostEqualUlps(float a,float b)129 bool AlmostEqualUlps(float a, float b) {
130     const int UlpsEpsilon = 16;
131     return equal_ulps(a, b, UlpsEpsilon, UlpsEpsilon);
132 }
133 
AlmostEqualUlpsNoNormalCheck(float a,float b)134 bool AlmostEqualUlpsNoNormalCheck(float a, float b) {
135     const int UlpsEpsilon = 16;
136     return equal_ulps_no_normal_check(a, b, UlpsEpsilon, UlpsEpsilon);
137 }
138 
AlmostEqualUlps_Pin(float a,float b)139 bool AlmostEqualUlps_Pin(float a, float b) {
140     const int UlpsEpsilon = 16;
141     return equal_ulps_pin(a, b, UlpsEpsilon, UlpsEpsilon);
142 }
143 
NotAlmostEqualUlps(float a,float b)144 bool NotAlmostEqualUlps(float a, float b) {
145     const int UlpsEpsilon = 16;
146     return not_equal_ulps(a, b, UlpsEpsilon);
147 }
148 
NotAlmostEqualUlps_Pin(float a,float b)149 bool NotAlmostEqualUlps_Pin(float a, float b) {
150     const int UlpsEpsilon = 16;
151     return not_equal_ulps_pin(a, b, UlpsEpsilon);
152 }
153 
NotAlmostDequalUlps(float a,float b)154 bool NotAlmostDequalUlps(float a, float b) {
155     const int UlpsEpsilon = 16;
156     return d_not_equal_ulps(a, b, UlpsEpsilon);
157 }
158 
RoughlyEqualUlps(float a,float b)159 bool 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 
AlmostBetweenUlps(float a,float b,float c)165 bool 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 
AlmostLessUlps(float a,float b)171 bool AlmostLessUlps(float a, float b) {
172     const int UlpsEpsilon = 16;
173     return less_ulps(a, b, UlpsEpsilon);
174 }
175 
AlmostLessOrEqualUlps(float a,float b)176 bool AlmostLessOrEqualUlps(float a, float b) {
177     const int UlpsEpsilon = 16;
178     return less_or_equal_ulps(a, b, UlpsEpsilon);
179 }
180 
UlpsDistance(float a,float b)181 int 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
cbrt_5d(double d)196 static 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)
cbrta_halleyd(const double a,const double R)206 static 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)
halley_cbrt3d(double d)213 static 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 
SkDCubeRoot(double x)220 double 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 
SkOpGlobalState(SkOpContourHead * head,SkArenaAlloc * allocator SkDEBUGPARAMS (bool debugSkipAssert)SkDEBUGPARAMS (const char * testName))231 SkOpGlobalState::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