1// Implementation of 1D, 2D, and 3D Perlin noise. Based on the
2// C code by Paul Bourke:
3// http://local.wasp.uwa.edu.au/~pbourke/texture_colour/perlin/
4class Perlin {
5  int B = 0x100;
6  int BM = 0xff;
7  int N = 0x1000;
8  int NP = 12;
9  int NM = 0xfff;
10
11  int p[];
12  float g3[][];
13  float g2[][];
14  float g1[];
15
16  void normalize2(float v[]) {
17    float s = sqrt(v[0] * v[0] + v[1] * v[1]);
18    v[0] = v[0] / s;
19    v[1] = v[1] / s;
20  }
21
22  void normalize3(float v[]) {
23    float s = sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
24    v[0] = v[0] / s;
25    v[1] = v[1] / s;
26    v[2] = v[2] / s;
27  }
28
29  float sCurve(float t) {
30    return t * t * (3.0 - 2.0 * t);
31  }
32
33  float at2(float q[], float rx, float ry) {
34    return rx * q[0] + ry * q[1];
35  }
36
37  float at3(float q[], float rx, float ry, float rz) {
38    return rx * q[0] + ry * q[1] + rz * q[2];
39  }
40
41  Perlin() {
42    p = new int[B + B + 2];
43    g3 = new float[B + B + 2][3];
44    g2 = new float[B + B + 2][2];
45    g1 = new float[B + B + 2];
46
47    init();
48  }
49
50  void init() {
51    int i, j, k;
52
53    for (i = 0 ; i < B ; i++) {
54      p[i] = i;
55      g1[i] = (random(B + B) - B) / B;
56
57      for (j = 0 ; j < 2 ; j++)
58        g2[i][j] = (random(B + B) - B) / B;
59      normalize2(g2[i]);
60
61      for (j = 0 ; j < 3 ; j++)
62        g3[i][j] = (random(B + B) - B) / B;
63      normalize3(g3[i]);
64    }
65
66    while (0 < --i) {
67      k = p[i];
68      p[i] = p[j = int(random(B))];
69      p[j] = k;
70    }
71
72    for (i = 0 ; i < B + 2 ; i++) {
73      p[B + i] = p[i];
74      g1[B + i] = g1[i];
75      for (j = 0 ; j < 2 ; j++)
76        g2[B + i][j] = g2[i][j];
77      for (j = 0 ; j < 3 ; j++)
78        g3[B + i][j] = g3[i][j];
79    }
80  }
81
82  float noise1(float[] vec) {
83    int bx0, bx1;
84    float rx0, rx1, sx, t, u, v;
85
86    t = vec[0] + N;
87    bx0 = int(t) & BM;
88    bx1 = (bx0 + 1) & BM;
89    rx0 = t - int(t);
90    rx1 = rx0 - 1.0;
91
92    sx = sCurve(rx0);
93    u = rx0 * g1[p[bx0]];
94    v = rx1 * g1[p[bx1]];
95
96    return lerp(u, v, sx);
97  }
98
99  float noise2(float[] vec) {
100    int bx0, bx1, by0, by1, b00, b10, b01, b11;
101    float rx0, rx1, ry0, ry1, sx, sy, a, b, t, u, v;
102    float[] q;
103    int i, j;
104
105    t = vec[0] + N;
106    bx0 = int(t) & BM;
107    bx1 = (bx0 + 1) & BM;
108    rx0 = t - int(t);
109    rx1 = rx0 - 1.0;
110
111    t = vec[1] + N;
112    by0 = int(t) & BM;
113    by1 = (by0 + 1) & BM;
114    ry0 = t - int(t);
115    ry1 = ry0 - 1.0;
116
117    i = p[bx0];
118    j = p[bx1];
119
120    b00 = p[i + by0];
121    b10 = p[j + by0];
122    b01 = p[i + by1];
123    b11 = p[j + by1];
124
125    sx = sCurve(rx0);
126    sy = sCurve(ry0);
127
128    q = g2[b00];
129    u = at2(q, rx0, ry0);
130    q = g2[b10];
131    v = at2(q, rx1, ry0);
132    a = lerp(u, v, sx);
133
134    q = g2[b01] ;
135    u = at2(q, rx0, ry1);
136    q = g2[b11] ;
137    v = at2(q, rx1, ry1);
138    b = lerp(u, v, sx);
139
140    return lerp(a, b, sy);
141  }
142
143  float noise3(float[] vec) {
144    int bx0, bx1, by0, by1, bz0, bz1, b00, b10, b01, b11;
145    float rx0, rx1, ry0, ry1, rz0, rz1, sy, sz, a, b, c, d, t, u, v;
146    float[] q;
147    int i, j;
148
149    t = vec[0] + N;
150    bx0 = int(t) & BM;
151    bx1 = (bx0 + 1) & BM;
152    rx0 = t - int(t);
153    rx1 = rx0 - 1.0;
154
155    t = vec[1] + N;
156    by0 = int(t) & BM;
157    by1 = (by0 + 1) & BM;
158    ry0 = t - int(t);
159    ry1 = ry0 - 1.0;
160
161    t = vec[2] + N;
162    bz0 = int(t) & BM;
163    bz1 = (bz0 + 1) & BM;
164    rz0 = t - int(t);
165    rz1 = rz0 - 1.0;
166
167    i = p[bx0];
168    j = p[bx1];
169
170    b00 = p[i + by0];
171    b10 = p[j + by0];
172    b01 = p[i + by1];
173    b11 = p[j + by1];
174
175    t  = sCurve(rx0);
176    sy = sCurve(ry0);
177    sz = sCurve(rz0);
178
179    q = g3[b00 + bz0];
180    u = at3(q, rx0, ry0, rz0);
181    q = g3[b10 + bz0];
182    v = at3(q, rx1, ry0, rz0);
183    a = lerp(u, v, t);
184
185    q = g3[b01 + bz0];
186    u = at3(q, rx0, ry1, rz0);
187    q = g3[b11 + bz0];
188    v = at3(q, rx1, ry1, rz0);
189    b = lerp(u, v, t);
190
191    c = lerp(a, b, sy);
192
193    q = g3[b00 + bz1];
194    u = at3(q, rx0, ry0, rz1);
195    q = g3[b10 + bz1];
196    v = at3(q, rx1, ry0, rz1);
197    a = lerp(u, v, t);
198
199    q = g3[b01 + bz1];
200    u = at3(q, rx0, ry1, rz1);
201    q = g3[b11 + bz1];
202    v = at3(q, rx1, ry1, rz1);
203    b = lerp(u, v, t);
204
205    d = lerp(a, b, sy);
206
207    return lerp(c, d, sz);
208  }
209
210  // In what follows "nalpha" is the weight when the sum is formed.
211  // Typically it is 2, as this approaches 1 the function is noisier.
212  // "nbeta" is the harmonic scaling/spacing, typically 2. n is the
213  // number of harmonics added up in the final result. Higher number
214  // results in more detailed noise.
215
216  float noise1D(float x, float nalpha, float nbeta, int n) {
217    float val, sum = 0;
218    float v[] = {x};
219    float nscale = 1;
220
221    for (int i = 0; i < n; i++) {
222      val = noise1(v);
223      sum += val / nscale;
224      nscale *= nalpha;
225      v[0] *= nbeta;
226    }
227    return sum;
228  }
229
230  float noise2D(float x, float y, float nalpha, float nbeta, int n) {
231   float val,sum = 0;
232   float v[] = {x, y};
233   float nscale = 1;
234
235   for (int i = 0; i < n; i++) {
236      val = noise2(v);
237      sum += val / nscale;
238      nscale *= nalpha;
239      v[0] *= nbeta;
240      v[1] *= nbeta;
241   }
242   return sum;
243  }
244
245  float noise3D(float x, float y, float z, float nalpha, float nbeta, int n) {
246    float val, sum = 0;
247    float v[] = {x, y, z};
248    float nscale = 1;
249
250    for (int i = 0 ; i < n; i++) {
251      val = noise3(v);
252      sum += val / nscale;
253      nscale *= nalpha;
254      v[0] *= nbeta;
255      v[1] *= nbeta;
256      v[2] *= nbeta;
257    }
258    return sum;
259  }
260}
261
262