1/*
2 * Copyright (C) 2006-2008 The Android Open Source Project
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 *      http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16
17#include "SkPoint.h"
18
19void SkIPoint::rotateCW(SkIPoint* dst) const {
20    SkASSERT(dst);
21
22    // use a tmp in case this == dst
23    int32_t tmp = fX;
24    dst->fX = -fY;
25    dst->fY = tmp;
26}
27
28void SkIPoint::rotateCCW(SkIPoint* dst) const {
29    SkASSERT(dst);
30
31    // use a tmp in case this == dst
32    int32_t tmp = fX;
33    dst->fX = fY;
34    dst->fY = -tmp;
35}
36
37///////////////////////////////////////////////////////////////////////////////
38
39void SkPoint::setIRectFan(int l, int t, int r, int b, size_t stride) {
40    SkASSERT(stride >= sizeof(SkPoint));
41
42    ((SkPoint*)((intptr_t)this + 0 * stride))->set(SkIntToScalar(l),
43                                                   SkIntToScalar(t));
44    ((SkPoint*)((intptr_t)this + 1 * stride))->set(SkIntToScalar(l),
45                                                   SkIntToScalar(b));
46    ((SkPoint*)((intptr_t)this + 2 * stride))->set(SkIntToScalar(r),
47                                                   SkIntToScalar(b));
48    ((SkPoint*)((intptr_t)this + 3 * stride))->set(SkIntToScalar(r),
49                                                   SkIntToScalar(t));
50}
51
52void SkPoint::setRectFan(SkScalar l, SkScalar t, SkScalar r, SkScalar b,
53                         size_t stride) {
54    SkASSERT(stride >= sizeof(SkPoint));
55
56    ((SkPoint*)((intptr_t)this + 0 * stride))->set(l, t);
57    ((SkPoint*)((intptr_t)this + 1 * stride))->set(l, b);
58    ((SkPoint*)((intptr_t)this + 2 * stride))->set(r, b);
59    ((SkPoint*)((intptr_t)this + 3 * stride))->set(r, t);
60}
61
62void SkPoint::rotateCW(SkPoint* dst) const {
63    SkASSERT(dst);
64
65    // use a tmp in case this == dst
66    SkScalar tmp = fX;
67    dst->fX = -fY;
68    dst->fY = tmp;
69}
70
71void SkPoint::rotateCCW(SkPoint* dst) const {
72    SkASSERT(dst);
73
74    // use a tmp in case this == dst
75    SkScalar tmp = fX;
76    dst->fX = fY;
77    dst->fY = -tmp;
78}
79
80void SkPoint::scale(SkScalar scale, SkPoint* dst) const {
81    SkASSERT(dst);
82    dst->set(SkScalarMul(fX, scale), SkScalarMul(fY, scale));
83}
84
85#define kNearlyZero     (SK_Scalar1 / 8092)
86
87bool SkPoint::normalize() {
88    return this->setLength(fX, fY, SK_Scalar1);
89}
90
91bool SkPoint::setNormalize(SkScalar x, SkScalar y) {
92    return this->setLength(x, y, SK_Scalar1);
93}
94
95bool SkPoint::setLength(SkScalar length) {
96    return this->setLength(fX, fY, length);
97}
98
99#ifdef SK_SCALAR_IS_FLOAT
100
101SkScalar SkPoint::Length(SkScalar dx, SkScalar dy) {
102    return sk_float_sqrt(dx * dx + dy * dy);
103}
104
105SkScalar SkPoint::Normalize(SkPoint* pt) {
106    float mag = SkPoint::Length(pt->fX, pt->fY);
107    if (mag > kNearlyZero) {
108        float scale = 1 / mag;
109        pt->fX *= scale;
110        pt->fY *= scale;
111        return mag;
112    }
113    return 0;
114}
115
116bool SkPoint::setLength(float x, float y, float length) {
117    float mag = sk_float_sqrt(x * x + y * y);
118    if (mag > kNearlyZero) {
119        length /= mag;
120        fX = x * length;
121        fY = y * length;
122        return true;
123    }
124    return false;
125}
126
127#else
128
129#include "Sk64.h"
130
131SkScalar SkPoint::Length(SkScalar dx, SkScalar dy) {
132    Sk64    tmp1, tmp2;
133
134    tmp1.setMul(dx, dx);
135    tmp2.setMul(dy, dy);
136    tmp1.add(tmp2);
137
138    return tmp1.getSqrt();
139}
140
141#ifdef SK_DEBUGx
142static SkFixed fixlen(SkFixed x, SkFixed y) {
143    float fx = (float)x;
144    float fy = (float)y;
145
146    return (int)floorf(sqrtf(fx*fx + fy*fy) + 0.5f);
147}
148#endif
149
150static inline uint32_t squarefixed(unsigned x) {
151    x >>= 16;
152    return x*x;
153}
154
155#if 1   // Newton iter for setLength
156
157static inline unsigned invsqrt_iter(unsigned V, unsigned U) {
158    unsigned x = V * U >> 14;
159    x = x * U >> 14;
160    x = (3 << 14) - x;
161    x = (U >> 1) * x >> 14;
162    return x;
163}
164
165static const uint16_t gInvSqrt14GuessTable[] = {
166    0x4000, 0x3c57, 0x393e, 0x3695, 0x3441, 0x3235, 0x3061,
167    0x2ebd, 0x2d41, 0x2be7, 0x2aaa, 0x2987, 0x287a, 0x2780,
168    0x2698, 0x25be, 0x24f3, 0x2434, 0x2380, 0x22d6, 0x2235,
169    0x219d, 0x210c, 0x2083, 0x2000, 0x1f82, 0x1f0b, 0x1e99,
170    0x1e2b, 0x1dc2, 0x1d5d, 0x1cfc, 0x1c9f, 0x1c45, 0x1bee,
171    0x1b9b, 0x1b4a, 0x1afc, 0x1ab0, 0x1a67, 0x1a20, 0x19dc,
172    0x1999, 0x1959, 0x191a, 0x18dd, 0x18a2, 0x1868, 0x1830,
173    0x17fa, 0x17c4, 0x1791, 0x175e, 0x172d, 0x16fd, 0x16ce
174};
175
176#define BUILD_INVSQRT_TABLEx
177#ifdef BUILD_INVSQRT_TABLE
178static void build_invsqrt14_guess_table() {
179    for (int i = 8; i <= 63; i++) {
180        unsigned x = SkToU16((1 << 28) / SkSqrt32(i << 25));
181        printf("0x%x, ", x);
182    }
183    printf("\n");
184}
185#endif
186
187static unsigned fast_invsqrt(uint32_t x) {
188#ifdef BUILD_INVSQRT_TABLE
189    unsigned top2 = x >> 25;
190    SkASSERT(top2 >= 8 && top2 <= 63);
191
192    static bool gOnce;
193    if (!gOnce) {
194        build_invsqrt14_guess_table();
195        gOnce = true;
196    }
197#endif
198
199    unsigned V = x >> 14;   // make V .14
200
201    unsigned top = x >> 25;
202    SkASSERT(top >= 8 && top <= 63);
203    SkASSERT(top - 8 < SK_ARRAY_COUNT(gInvSqrt14GuessTable));
204    unsigned U = gInvSqrt14GuessTable[top - 8];
205
206    U = invsqrt_iter(V, U);
207    return invsqrt_iter(V, U);
208}
209
210/*  We "normalize" x,y to be .14 values (so we can square them and stay 32bits.
211    Then we Newton-iterate this in .14 space to compute the invser-sqrt, and
212    scale by it at the end. The .14 space means we can execute our iterations
213    and stay in 32bits as well, making the multiplies much cheaper than calling
214    SkFixedMul.
215*/
216bool SkPoint::setLength(SkFixed ox, SkFixed oy, SkFixed length) {
217    if (ox == 0) {
218        if (oy == 0) {
219            return false;
220        }
221        this->set(0, SkApplySign(length, SkExtractSign(oy)));
222        return true;
223    }
224    if (oy == 0) {
225        this->set(SkApplySign(length, SkExtractSign(ox)), 0);
226        return true;
227    }
228
229    unsigned x = SkAbs32(ox);
230    unsigned y = SkAbs32(oy);
231    int zeros = SkCLZ(x | y);
232
233    // make x,y 1.14 values so our fast sqr won't overflow
234    if (zeros > 17) {
235        x <<= zeros - 17;
236        y <<= zeros - 17;
237    } else {
238        x >>= 17 - zeros;
239        y >>= 17 - zeros;
240    }
241    SkASSERT((x | y) <= 0x7FFF);
242
243    unsigned invrt = fast_invsqrt(x*x + y*y);
244
245    x = x * invrt >> 12;
246    y = y * invrt >> 12;
247
248    if (length != SK_Fixed1) {
249        x = SkFixedMul(x, length);
250        y = SkFixedMul(y, length);
251    }
252    this->set(SkApplySign(x, SkExtractSign(ox)),
253              SkApplySign(y, SkExtractSign(oy)));
254    return true;
255}
256#else
257/*
258    Normalize x,y, and then scale them by length.
259
260    The obvious way to do this would be the following:
261        S64 tmp1, tmp2;
262        tmp1.setMul(x,x);
263        tmp2.setMul(y,y);
264        tmp1.add(tmp2);
265        len = tmp1.getSqrt();
266        x' = SkFixedDiv(x, len);
267        y' = SkFixedDiv(y, len);
268    This is fine, but slower than what we do below.
269
270    The present technique does not compute the starting length, but
271    rather fiddles with x,y iteratively, all the while checking its
272    magnitude^2 (avoiding a sqrt).
273
274    We normalize by first shifting x,y so that at least one of them
275    has bit 31 set (after taking the abs of them).
276    Then we loop, refining x,y by squaring them and comparing
277    against a very large 1.0 (1 << 28), and then adding or subtracting
278    a delta (which itself is reduced by half each time through the loop).
279    For speed we want the squaring to be with a simple integer mul. To keep
280    that from overflowing we shift our coordinates down until we are dealing
281    with at most 15 bits (2^15-1)^2 * 2 says withing 32 bits)
282    When our square is close to 1.0, we shift x,y down into fixed range.
283*/
284bool SkPoint::setLength(SkFixed ox, SkFixed oy, SkFixed length) {
285    if (ox == 0) {
286        if (oy == 0)
287            return false;
288        this->set(0, SkApplySign(length, SkExtractSign(oy)));
289        return true;
290    }
291    if (oy == 0) {
292        this->set(SkApplySign(length, SkExtractSign(ox)), 0);
293        return true;
294    }
295
296    SkFixed x = SkAbs32(ox);
297    SkFixed y = SkAbs32(oy);
298
299    // shift x,y so that the greater of them is 15bits (1.14 fixed point)
300    {
301        int shift = SkCLZ(x | y);
302        // make them .30
303        x <<= shift - 1;
304        y <<= shift - 1;
305    }
306
307    SkFixed dx = x;
308    SkFixed dy = y;
309
310    for (int i = 0; i < 17; i++) {
311        dx >>= 1;
312        dy >>= 1;
313
314        U32 len2 = squarefixed(x) + squarefixed(y);
315        if (len2 >> 28) {
316            x -= dx;
317            y -= dy;
318        } else {
319            x += dx;
320            y += dy;
321        }
322    }
323    x >>= 14;
324    y >>= 14;
325
326#ifdef SK_DEBUGx    // measure how far we are from unit-length
327    {
328        static int gMaxError;
329        static int gMaxDiff;
330
331        SkFixed len = fixlen(x, y);
332        int err = len - SK_Fixed1;
333        err = SkAbs32(err);
334
335        if (err > gMaxError) {
336            gMaxError = err;
337            SkDebugf("gMaxError %d\n", err);
338        }
339
340        float fx = SkAbs32(ox)/65536.0f;
341        float fy = SkAbs32(oy)/65536.0f;
342        float mag = sqrtf(fx*fx + fy*fy);
343        fx /= mag;
344        fy /= mag;
345        SkFixed xx = (int)floorf(fx * 65536 + 0.5f);
346        SkFixed yy = (int)floorf(fy * 65536 + 0.5f);
347        err = SkMax32(SkAbs32(xx-x), SkAbs32(yy-y));
348        if (err > gMaxDiff) {
349            gMaxDiff = err;
350            SkDebugf("gMaxDiff %d\n", err);
351        }
352    }
353#endif
354
355    x = SkApplySign(x, SkExtractSign(ox));
356    y = SkApplySign(y, SkExtractSign(oy));
357    if (length != SK_Fixed1) {
358        x = SkFixedMul(x, length);
359        y = SkFixedMul(y, length);
360    }
361
362    this->set(x, y);
363    return true;
364}
365#endif
366
367#endif
368
369///////////////////////////////////////////////////////////////////////////////
370
371SkScalar SkPoint::distanceToLineSegmentBetweenSqd(const SkPoint& a,
372                                                  const SkPoint& b) const {
373    // See comments to distanceToLineBetweenSqd. If the projection of c onto
374    // u is between a and b then this returns the same result as that
375    // function. Otherwise, it returns the distance to the closer of a and
376    // b. Let the projection of v onto u be v'.  There are three cases:
377    //    1. v' points opposite to u. c is not between a and b and is closer
378    //       to a than b.
379    //    2. v' points along u and has magnitude less than y. c is between
380    //       a and b and the distance to the segment is the same as distance
381    //       to the line ab.
382    //    3. v' points along u and has greater magnitude than u. c is not
383    //       not between a and b and is closer to b than a.
384    // v' = (u dot v) * u / |u|. So if (u dot v)/|u| is less than zero we're
385    // in case 1. If (u dot v)/|u| is > |u| we are in case 3. Otherwise
386    // we're in case 2. We actually compare (u dot v) to 0 and |u|^2 to
387    // avoid a sqrt to compute |u|.
388
389    SkVector u = b - a;
390    SkVector v = *this - a;
391
392    SkScalar uLengthSqd = u.lengthSqd();
393    SkScalar uDotV = SkPoint::DotProduct(u, v);
394
395    if (uDotV <= 0) {
396        return v.lengthSqd();
397    } else if (uDotV > uLengthSqd) {
398        return b.distanceToSqd(*this);
399    } else {
400        SkScalar det = u.cross(v);
401        return SkScalarMulDiv(det, det, uLengthSqd);
402    }
403}
404
405