1 /**************************************************************************/
2 /* Copyright 2009 Tim Day */
3 /* */
4 /* This file is part of Fracplanet */
5 /* */
6 /* Fracplanet is free software: you can redistribute it and/or modify */
7 /* it under the terms of the GNU General Public License as published by */
8 /* the Free Software Foundation, either version 3 of the License, or */
9 /* (at your option) any later version. */
10 /* */
11 /* Fracplanet is distributed in the hope that it will be useful, */
12 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
13 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
14 /* GNU General Public License for more details. */
15 /* */
16 /* You should have received a copy of the GNU General Public License */
17 /* along with Fracplanet. If not, see <http://www.gnu.org/licenses/>. */
18 /**************************************************************************/
19
20 /*! \file
21 \brief Implementation of Noise class
22 */
23
24 #include "noise.h"
25
Noise(uint seed)26 Noise::Noise(uint seed)
27 {
28 // We use our own random number generator so saved pictures will be the same when reloaded!
29 // (Besides, one isn't actually conveniently available)
30 Random01 r_01(seed);
31
32 // Create an array of random gradient vectors uniformly on the unit sphere
33
34 unsigned int i;
35
36 for (i=0;i<N;i++)
37 _g[i]=RandomXYZSphereNormal(r_01);
38
39 // Create a pseudorandom permutation of [1..B]
40
41 for (i=0;i<=N;i++)
42 _p[i]=i;
43
44 for (i=N;i>0;i-=2)
45 {
46 const int j=(static_cast<int>(r_01()*N));
47 int k=_p[i];
48 _p[i]=_p[j];
49 _p[j]=k;
50 }
51
52 // Extend g and p arrays to allow for faster indexing
53
54 for (i=0;i<N+2;i++)
55 {
56 _p[N+i]=_p[i];
57 _g[N+i]=_g[i];
58 }
59 }
60
~Noise()61 Noise::~Noise()
62 {}
63
value(const XYZ & q,float rx,float ry,float rz)64 inline float value(const XYZ& q,float rx,float ry,float rz)
65 {
66 return rx*q.x+ry*q.y+rz*q.z;
67 }
68
surve(float t)69 inline float surve(float t)
70 {
71 return t*t*(3.0-2.0*t);
72 }
73
lerp(float t,float a,float b)74 inline float lerp(float t,float a,float b)
75 {
76 return a+t*(b-a);
77 }
78
operator ()(const XYZ & p) const79 float Noise::operator()(const XYZ& p) const
80 {
81 // Crank up the frequency a bit otherwise don't see much variation in base case
82 const float tx=2.0*p.x+10000.0f;
83 const float ty=2.0*p.y+10000.0f;
84 const float tz=2.0*p.z+10000.0f;
85
86 const int itx=static_cast<int>(tx);
87 const int ity=static_cast<int>(ty);
88 const int itz=static_cast<int>(tz);
89
90 const int bx0=( itx &(N-1));
91 const int bx1=((bx0+1)&(N-1));
92 const int by0=( ity &(N-1));
93 const int by1=((by0+1)&(N-1));
94 const int bz0=( itz &(N-1));
95 const int bz1=((bz0+1)&(N-1));
96
97 const int i=_p[bx0];
98 const int b00=_p[i+by0];
99 const int b01=_p[i+by1];
100
101 const int j=_p[bx1];
102 const int b10=_p[j+by0];
103 const int b11=_p[j+by1];
104
105 const float rx0=tx-itx;
106 const float ry0=ty-ity;
107 const float rz0=tz-itz;
108
109 const float rx1=rx0-1.0;
110 const float ry1=ry0-1.0;
111 const float rz1=rz0-1.0;
112
113 const float sx=surve(rx0);
114
115 const float a0=lerp(sx,value(_g[b00+bz0],rx0,ry0,rz0),value(_g[b10+bz0],rx1,ry0,rz0));
116 const float b0=lerp(sx,value(_g[b01+bz0],rx0,ry1,rz0),value(_g[b11+bz0],rx1,ry1,rz0));
117 const float a1=lerp(sx,value(_g[b00+bz1],rx0,ry0,rz1),value(_g[b10+bz1],rx1,ry0,rz1));
118 const float b1=lerp(sx,value(_g[b01+bz1],rx0,ry1,rz1),value(_g[b11+bz1],rx1,ry1,rz1));
119
120 const float sy=surve(ry0);
121
122 const float c=lerp(sy,a0,b0);
123 const float d=lerp(sy,a1,b1);
124
125 const float sz=surve(rz0);
126
127 return 1.5*lerp(sz,c,d);
128 }
129
MultiscaleNoise(uint seed,uint terms,float decay)130 MultiscaleNoise::MultiscaleNoise(uint seed,uint terms,float decay)
131 :_terms(terms)
132 ,_noise(new boost::scoped_ptr<const Noise>[terms])
133 ,_amplitude(new float[terms])
134 {
135 float k=1.0f;
136 float kt=0.0f;
137 for (uint i=0;i<_terms;i++)
138 {
139 _noise[i].reset(new Noise(seed+i));
140 _amplitude[i]=k;
141 kt+=k;
142 k*=decay;
143 }
144 for (uint i=0;i<_terms;i++)
145 {
146 _amplitude[i]/=kt;
147 }
148 }
149
~MultiscaleNoise()150 MultiscaleNoise::~MultiscaleNoise()
151 {}
152
153 //! Return noise value at a point.
operator ()(const XYZ & p) const154 float MultiscaleNoise::operator()(const XYZ& p) const
155 {
156 float v=0.0;
157 for (uint i=0;i<_terms;i++)
158 {
159 v+=_amplitude[i]*(*_noise[i])(p*(1<<i));
160 }
161 return v;
162 }
163