1// file  is taken from aquilax/go-perlin
2
3// Package perlin provides coherent noise function over 1, 2 or 3 dimensions
4// This code is go adaptagion based on C implementation that can be found here:
5// http://git.gnome.org/browse/gegl/tree/operations/common/perlin/perlin.c
6// (original copyright Ken Perlin)
7package perlin
8
9import (
10	"math"
11	"math/rand"
12)
13
14// General constants
15const (
16	B  = 0x100
17	N  = 0x1000
18	BM = 0xff
19)
20
21// Perlin is the noise generator
22type Perlin struct {
23	alpha float64
24	beta  float64
25	n     int
26
27	p  [B + B + 2]int
28	g3 [B + B + 2][3]float64
29	g2 [B + B + 2][2]float64
30	g1 [B + B + 2]float64
31}
32
33// NewPerlin creates new Perlin noise generator
34// In what follows "alpha" is the weight when the sum is formed.
35// Typically it is 2, As this approaches 1 the function is noisier.
36// "beta" is the harmonic scaling/spacing, typically 2, n is the
37// number of iterations and seed is the math.rand seed value to use
38func NewPerlin(alpha, beta float64, n int, seed int64) *Perlin {
39	return NewPerlinRandSource(alpha, beta, n, rand.NewSource(seed))
40}
41
42// NewPerlinRandSource creates new Perlin noise generator
43// In what follows "alpha" is the weight when the sum is formed.
44// Typically it is 2, As this approaches 1 the function is noisier.
45// "beta" is the harmonic scaling/spacing, typically 2, n is the
46// number of iterations and source is source of pseudo-random int64 values
47func NewPerlinRandSource(alpha, beta float64, n int, source rand.Source) *Perlin {
48	var p Perlin
49	var i int
50
51	p.alpha = alpha
52	p.beta = beta
53	p.n = n
54
55	r := rand.New(source)
56
57	for i = 0; i < B; i++ {
58		p.p[i] = i
59		p.g1[i] = float64((r.Int()%(B+B))-B) / B
60
61		for j := 0; j < 2; j++ {
62			p.g2[i][j] = float64((r.Int()%(B+B))-B) / B
63		}
64
65		normalize2(&p.g2[i])
66
67		for j := 0; j < 3; j++ {
68			p.g3[i][j] = float64((r.Int()%(B+B))-B) / B
69		}
70		normalize3(&p.g3[i])
71	}
72
73	for ; i > 0; i-- {
74		k := p.p[i]
75		j := r.Int() % B
76		p.p[i] = p.p[j]
77		p.p[j] = k
78	}
79
80	for i := 0; i < B+2; i++ {
81		p.p[B+i] = p.p[i]
82		p.g1[B+i] = p.g1[i]
83		for j := 0; j < 2; j++ {
84			p.g2[B+i][j] = p.g2[i][j]
85		}
86		for j := 0; j < 3; j++ {
87			p.g3[B+i][j] = p.g3[i][j]
88		}
89	}
90
91	return &p
92}
93
94func normalize2(v *[2]float64) {
95	s := math.Sqrt(v[0]*v[0] + v[1]*v[1])
96	v[0] = v[0] / s
97	v[1] = v[1] / s
98}
99
100func normalize3(v *[3]float64) {
101	s := math.Sqrt(v[0]*v[0] + v[1]*v[1] + v[2]*v[2])
102	v[0] = v[0] / s
103	v[1] = v[1] / s
104	v[2] = v[2] / s
105}
106
107func at2(rx, ry float64, q [2]float64) float64 {
108	return rx*q[0] + ry*q[1]
109}
110
111func at3(rx, ry, rz float64, q [3]float64) float64 {
112	return rx*q[0] + ry*q[1] + rz*q[2]
113}
114
115func sCurve(t float64) float64 {
116	return t * t * (3. - 2.*t)
117}
118
119func lerp(t, a, b float64) float64 {
120	return a + t*(b-a)
121}
122
123func (p *Perlin) noise1(arg float64) float64 {
124	var vec [1]float64
125	vec[0] = arg
126
127	t := vec[0] + N
128	bx0 := int(t) & BM
129	bx1 := (bx0 + 1) & BM
130	rx0 := t - float64(int(t))
131	rx1 := rx0 - 1.
132
133	sx := sCurve(rx0)
134	u := rx0 * p.g1[p.p[bx0]]
135	v := rx1 * p.g1[p.p[bx1]]
136
137	return lerp(sx, u, v)
138}
139
140func (p *Perlin) noise2(vec [2]float64) float64 {
141
142	t := vec[0] + N
143	bx0 := int(t) & BM
144	bx1 := (bx0 + 1) & BM
145	rx0 := t - float64(int(t))
146	rx1 := rx0 - 1.
147
148	t = vec[1] + N
149	by0 := int(t) & BM
150	by1 := (by0 + 1) & BM
151	ry0 := t - float64(int(t))
152	ry1 := ry0 - 1.
153
154	i := p.p[bx0]
155	j := p.p[bx1]
156
157	b00 := p.p[i+by0]
158	b10 := p.p[j+by0]
159	b01 := p.p[i+by1]
160	b11 := p.p[j+by1]
161
162	sx := sCurve(rx0)
163	sy := sCurve(ry0)
164
165	q := p.g2[b00]
166	u := at2(rx0, ry0, q)
167	q = p.g2[b10]
168	v := at2(rx1, ry0, q)
169	a := lerp(sx, u, v)
170
171	q = p.g2[b01]
172	u = at2(rx0, ry1, q)
173	q = p.g2[b11]
174	v = at2(rx1, ry1, q)
175	b := lerp(sx, u, v)
176
177	return lerp(sy, a, b)
178}
179
180func (p *Perlin) noise3(vec [3]float64) float64 {
181	t := vec[0] + N
182	bx0 := int(t) & BM
183	bx1 := (bx0 + 1) & BM
184	rx0 := t - float64(int(t))
185	rx1 := rx0 - 1.
186
187	t = vec[1] + N
188	by0 := int(t) & BM
189	by1 := (by0 + 1) & BM
190	ry0 := t - float64(int(t))
191	ry1 := ry0 - 1.
192
193	t = vec[2] + N
194	bz0 := int(t) & BM
195	bz1 := (bz0 + 1) & BM
196	rz0 := t - float64(int(t))
197	rz1 := rz0 - 1.
198
199	i := p.p[bx0]
200	j := p.p[bx1]
201
202	b00 := p.p[i+by0]
203	b10 := p.p[j+by0]
204	b01 := p.p[i+by1]
205	b11 := p.p[j+by1]
206
207	t = sCurve(rx0)
208	sy := sCurve(ry0)
209	sz := sCurve(rz0)
210
211	q := p.g3[b00+bz0]
212	u := at3(rx0, ry0, rz0, q)
213	q = p.g3[b10+bz0]
214	v := at3(rx1, ry0, rz0, q)
215	a := lerp(t, u, v)
216
217	q = p.g3[b01+bz0]
218	u = at3(rx0, ry1, rz0, q)
219	q = p.g3[b11+bz0]
220	v = at3(rx1, ry1, rz0, q)
221	b := lerp(t, u, v)
222
223	c := lerp(sy, a, b)
224
225	q = p.g3[b00+bz1]
226	u = at3(rx0, ry0, rz1, q)
227	q = p.g3[b10+bz1]
228	v = at3(rx1, ry0, rz1, q)
229	a = lerp(t, u, v)
230
231	q = p.g3[b01+bz1]
232	u = at3(rx0, ry1, rz1, q)
233	q = p.g3[b11+bz1]
234	v = at3(rx1, ry1, rz1, q)
235	b = lerp(t, u, v)
236
237	d := lerp(sy, a, b)
238
239	return lerp(sz, c, d)
240}
241
242// Noise1D generates 1-dimensional Perlin Noise value
243func (p *Perlin) Noise1D(x float64) float64 {
244	var scale float64 = 1
245	var sum float64
246	px := x
247
248	for i := 0; i < p.n; i++ {
249		val := p.noise1(px)
250		sum += val / scale
251		scale *= p.alpha
252		px *= p.beta
253	}
254	return sum
255}
256
257// Noise2D Generates 2-dimensional Perlin Noise value
258func (p *Perlin) Noise2D(x, y float64) float64 {
259	var scale float64 = 1
260	var sum float64
261	var px [2]float64
262
263	px[0] = x
264	px[1] = y
265
266	for i := 0; i < p.n; i++ {
267		val := p.noise2(px)
268		sum += val / scale
269		scale *= p.alpha
270		px[0] *= p.beta
271		px[1] *= p.beta
272	}
273	return sum
274}
275
276// Noise3D Generates 3-dimensional Perlin Noise value
277func (p *Perlin) Noise3D(x, y, z float64) float64 {
278	var scale float64 = 1
279	var sum float64
280	var px [3]float64
281
282	if z < 0.0000 {
283		return p.Noise2D(x, y)
284	}
285	px[0] = x
286	px[1] = y
287	px[2] = z
288
289	for i := 0; i < p.n; i++ {
290		val := p.noise3(px)
291		sum += val / scale
292		scale *= p.alpha
293		px[0] *= p.beta
294		px[1] *= p.beta
295		px[2] *= p.beta
296	}
297	return sum
298}
299