1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
|
(*
coranac.com's awesome atan2 implementation
Very fast and very small.
*)
unit cearn_atan;
{$mode objfpc}
{$H+}
interface
uses
ctypes, nds9;
const
__qran_seed: cint = 42;
function atan2Lookup(x, y: cint): cuint32;
function atan2Lerp(x, y: cint): cuint32;
function sqran(seed: cint): cint;
function qran(): cint; inline;
function qran_range(min, max: cint): cint; inline;
const
BRAD_PI_SHIFT = 14;
BRAD_PI = 1 shl BRAD_PI_SHIFT;
BRAD_HPI = BRAD_PI div 2;
BRAD_2PI = BRAD_PI * 2;
ATAN_ONE = $1000;
ATAN_FP = 12;
// Some constants for dealing with atanLUT.
ATANLUT_STRIDE = ATAN_ONE div $80;
ATANLUT_STRIDE_SHIFT = 5;
// Arctangents LUT. Interval: [0, 1] (one=128); PI=0x20000
atanLUT: array [0..130-1] of cushort = (
$0000, $0146, $028C, $03D2, $0517, $065D, $07A2, $08E7,
$0A2C, $0B71, $0CB5, $0DF9, $0F3C, $107F, $11C1, $1303,
$1444, $1585, $16C5, $1804, $1943, $1A80, $1BBD, $1CFA,
$1E35, $1F6F, $20A9, $21E1, $2319, $2450, $2585, $26BA,
$27ED, $291F, $2A50, $2B80, $2CAF, $2DDC, $2F08, $3033,
$315D, $3285, $33AC, $34D2, $35F6, $3719, $383A, $395A,
$3A78, $3B95, $3CB1, $3DCB, $3EE4, $3FFB, $4110, $4224,
$4336, $4447, $4556, $4664, $4770, $487A, $4983, $4A8B,
// 64
$4B90, $4C94, $4D96, $4E97, $4F96, $5093, $518F, $5289,
$5382, $5478, $556E, $5661, $5753, $5843, $5932, $5A1E,
$5B0A, $5BF3, $5CDB, $5DC1, $5EA6, $5F89, $606A, $614A,
$6228, $6305, $63E0, $64B9, $6591, $6667, $673B, $680E,
$68E0, $69B0, $6A7E, $6B4B, $6C16, $6CDF, $6DA8, $6E6E,
$6F33, $6FF7, $70B9, $717A, $7239, $72F6, $73B3, $746D,
$7527, $75DF, $7695, $774A, $77FE, $78B0, $7961, $7A10,
$7ABF, $7B6B, $7C17, $7CC1, $7D6A, $7E11, $7EB7, $7F5C,
// 128
$8000, $80A2);
implementation
// Quick (and very dirty) pseudo-random number generator
// return random in range [0,8000h>
function qran(): cint; inline;
begin
__qran_seed := 1664525 * __qran_seed + 1013904223;
result := (__qran_seed shr 16) and $7FFF;
end;
function qran_range(min, max: cint): cint; inline;
begin
result := (qran() * (max - min) shr 15) + min;
end;
// Get the octant a coordinate pair is in.
procedure OCTANTIFY(var _x, _y, _o: cint); inline;
var
_t: cint;
begin
repeat
_o := 0;
if (_y < 0) then
begin
_x := -_x;
_y := -_y;
_o := _o + 4;
end;
if (_x <= 0) then
begin
_t := _x;
_x := _y;
_y := -_t;
_o := _o + 2;
end;
if (_x <= _y) then
begin
_t := _y - _x;
_x := _x + _y;
_y := _t;
_o := _o + 1;
end;
until true;
end;
function QDIV(num, den: cint; const bits: cint): cint; inline;
begin
while (REG_DIVCNT^ and DIV_BUSY) <> 0 do;
REG_DIVCNT^ := DIV_64_32;
REG_DIV_NUMER^ := cint64(num) shl bits;
REG_DIV_DENOM_L^ := den;
while (REG_DIVCNT^ and DIV_BUSY) <> 0 do;
result := REG_DIV_RESULT_L^;
end;
function atan2Lerp(x, y: cint): cuint32;
var
phi: cint;
t, fa, fb, h: cuint32;
begin
if (y =0) then
begin
if x >= 0 then
result := 0
else
result := BRAD_PI;
exit;
end;
OCTANTIFY(x, y, phi);
phi := phi * BRAD_PI div 4;
t := QDIV(y, x, ATAN_FP);
h := t mod ATANLUT_STRIDE;
fa := atanLUT[t div ATANLUT_STRIDE];
fb := atanLUT[t div ATANLUT_STRIDE + 1];
result := phi + (fa + SarLongint((fb - fa) * h, ATANLUT_STRIDE_SHIFT)) div 8;
end;
function atan2Lookup(x, y: cint): cuint32;
var
phi: cint;
t: cuint32;
begin
if (y = 0) then
begin
if x >= 0 then
result := 0
else
result := BRAD_PI;
exit;
end;
OCTANTIFY(x, y, phi);
phi := phi * BRAD_PI div 4;
t := QDIV(y, x, ATAN_FP);
result := phi + atanLUT[t div ATANLUT_STRIDE] div 8;
end;
function sqran(seed: cint): cint;
var
old: cint;
begin
old := __qran_seed;
__qran_seed := seed;
result := old;
end;
end.
|
