-
Notifications
You must be signed in to change notification settings - Fork 1
/
ext_makegrains.go
123 lines (105 loc) · 3 KB
/
ext_makegrains.go
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
package engine
import (
"math"
"math/rand"
)
func init() {
DeclFunc("ext_makegrains", Voronoi, "Voronoi tesselation (grain size, num regions)")
}
func Voronoi(grainsize float64, numRegions, seed int) {
Refer("Lel2014")
SetBusy(true)
defer SetBusy(false)
t := newTesselation(grainsize, numRegions, int64(seed))
regions.hist = append(regions.hist, t.RegionOf)
regions.render(t.RegionOf)
}
type tesselation struct {
grainsize float64
tilesize float64
maxRegion int
cache map[int2][]center
seed int64
rnd *rand.Rand
}
// integer tile coordinate
type int2 struct{ x, y int }
// Voronoi center info
type center struct {
x, y float64 // center position (m)
region byte // region for all cells near center
}
// nRegion exclusive
func newTesselation(grainsize float64, nRegion int, seed int64) *tesselation {
return &tesselation{grainsize,
float64(float32(grainsize * TILE)), // expect 4 grains/block, 36 per 3x3 blocks = safe, relatively round number
nRegion,
make(map[int2][]center),
seed,
rand.New(rand.NewSource(0))}
}
const (
TILE = 2 // tile size in grains
LAMBDA = TILE * TILE // expected grains per tile
)
// Returns the region of the grain where cell at x,y,z belongs to
func (t *tesselation) RegionOf(x, y, z float64) int {
tile := t.tileOf(x, y) // tile containing x,y
// look for nearest center in tile + neighbors
nearest := center{x, y, 0} // dummy initial value, but safe should the infinite impossibility strike.
mindist := math.Inf(1)
for tx := tile.x - 1; tx <= tile.x+1; tx++ {
for ty := tile.y - 1; ty <= tile.y+1; ty++ {
centers := t.centersInTile(tx, ty)
for _, c := range centers {
dist := sqr(x-c.x) + sqr(y-c.y)
if dist < mindist {
nearest = c
mindist = dist
}
}
}
}
//fmt.Println("nearest", x, y, ":", nearest)
return int(nearest.region)
}
// Returns the list of Voronoi centers in tile(ix, iy), using only ix,iy to seed the random generator
func (t *tesselation) centersInTile(tx, ty int) []center {
pos := int2{tx, ty}
if c, ok := t.cache[pos]; ok {
return c
} else {
// tile-specific seed that works for positive and negative tx, ty
seed := (int64(ty)+(1<<24))*(1<<24) + (int64(tx) + (1 << 24))
t.rnd.Seed(seed ^ t.seed)
N := t.poisson(LAMBDA)
c := make([]center, N)
// absolute position of tile (m)
x0, y0 := float64(tx)*t.tilesize, float64(ty)*t.tilesize
for i := range c {
// random position inside tile
c[i].x = x0 + t.rnd.Float64()*t.tilesize
c[i].y = y0 + t.rnd.Float64()*t.tilesize
c[i].region = byte(t.rnd.Intn(t.maxRegion))
}
t.cache[pos] = c
return c
}
}
func sqr(x float64) float64 { return x * x }
func (t *tesselation) tileOf(x, y float64) int2 {
ix := int(math.Floor(x / t.tilesize))
iy := int(math.Floor(y / t.tilesize))
return int2{ix, iy}
}
// Generate poisson distributed numbers (according to Knuth)
func (t *tesselation) poisson(lambda float64) int {
L := math.Exp(-lambda)
k := 1
p := t.rnd.Float64()
for p > L {
k++
p *= t.rnd.Float64()
}
return k - 1
}