|
3 | 3 | //
|
4 | 4 | // Authors: Stefan Gustavson ([email protected])
|
5 | 5 | // and Ian McEwan ([email protected])
|
6 |
| -// Version 2021-11-04, published under the MIT license (see below) |
| 6 | +// Version 2021-12-02, published under the MIT license (see below) |
7 | 7 | //
|
8 | 8 | // Copyright (c) 2021 Stefan Gustavson and Ian McEwan.
|
9 | 9 | //
|
|
56 | 56 | // and "out vec3 dg" returns the second order derivatives as
|
57 | 57 | // dg = (dn2/dx2, dn2/dy2, dn2/dxy)
|
58 | 58 | //
|
59 |
| -// Setting both periods to 0.0 or a negative value will skip the periodic |
60 |
| -// wrap and make the function execute 10% faster. If you want periodicity |
61 |
| -// for one dimensions but not the other, set the period to 0.0 only for |
62 |
| -// one of the dimensions. In that case, there is no significant speedup. |
| 59 | +// Setting either period to 0.0 or a negative value will skip the wrapping |
| 60 | +// along that dimension. Setting both periods to 0.0 makes the function |
| 61 | +// execute about 10% faster. |
63 | 62 | //
|
64 | 63 | // Not using the return value for some or all of the first or second order
|
65 | 64 | // derivatives will make the compiler eliminate the code for computing them.
|
|
70 | 69 | // of psrddnoise(), setting alpha == 0.0 gives no speedup.
|
71 | 70 | //
|
72 | 71 | float psrddnoise(vec2 x, vec2 period, float alpha, out vec2 gradient,
|
73 |
| - out vec3 dg) { |
74 |
| - |
75 |
| - // Transform to simplex space (axis-aligned hexagonal grid) |
76 |
| - vec2 uv = vec2(x.x + x.y*0.5, x.y); |
77 |
| - |
78 |
| - // Determine which simplex we're in, with i0 being the "base" |
79 |
| - vec2 i0 = floor(uv); |
80 |
| - vec2 f0 = fract(uv); |
81 |
| - // o1 is the offset in simplex space to the second corner |
82 |
| - float cmp = step(f0.y, f0.x); |
83 |
| - vec2 o1 = vec2(cmp, 1.0-cmp); |
84 |
| - |
85 |
| - // Enumerate the remaining simplex corners |
86 |
| - vec2 i1 = i0 + o1; |
87 |
| - vec2 i2 = i0 + vec2(1.0, 1.0); |
88 |
| - |
89 |
| - // Transform corners back to texture space |
90 |
| - vec2 v0 = vec2(i0.x - i0.y * 0.5, i0.y); |
91 |
| - vec2 v1 = vec2(v0.x + o1.x - o1.y * 0.5, v0.y + o1.y); |
92 |
| - vec2 v2 = vec2(v0.x + 0.5, v0.y + 1.0); |
93 |
| - |
94 |
| - // Compute vectors from v to each of the simplex corners |
95 |
| - vec2 x0 = x - v0; |
96 |
| - vec2 x1 = x - v1; |
97 |
| - vec2 x2 = x - v2; |
98 |
| - |
99 |
| - vec3 iu, iv; |
100 |
| - vec3 xw, yw; |
101 |
| - |
102 |
| - // Wrap to periods, if desired |
103 |
| - if(any(greaterThan(period, vec2(0.0)))) { |
104 |
| - xw = vec3(v0.x, v1.x, v2.x); |
105 |
| - yw = vec3(v0.y, v1.y, v2.y); |
106 |
| - if(period.x > 0.0) |
107 |
| - xw = mod(vec3(v0.x, v1.x, v2.x), period.x); |
108 |
| - if(period.y > 0.0) |
109 |
| - yw = mod(vec3(v0.y, v1.y, v2.y), period.y); |
110 |
| - // Transform back to simplex space and fix rounding errors |
111 |
| - iu = floor(xw + 0.5*yw + 0.5); |
112 |
| - iv = floor(yw + 0.5); |
113 |
| - } else { // Shortcut if neither x nor y periods are specified |
114 |
| - iu = vec3(i0.x, i1.x, i2.x); |
115 |
| - iv = vec3(i0.y, i1.y, i2.y); |
116 |
| - } |
117 |
| - |
118 |
| - // Compute one pseudo-random hash value for each corner |
119 |
| - vec3 hash = mod(iu, 289.0); |
120 |
| - hash = mod((hash*51.0 + 2.0)*hash + iv, 289.0); |
121 |
| - hash = mod((hash*34.0 + 10.0)*hash, 289.0); |
122 |
| - |
123 |
| - // Pick a pseudo-random angle and add the desired rotation |
124 |
| - vec3 psi = hash * 0.07482 + alpha; |
125 |
| - vec3 gx = cos(psi); |
126 |
| - vec3 gy = sin(psi); |
127 |
| - |
128 |
| - // Reorganize for dot products below |
129 |
| - vec2 g0 = vec2(gx.x,gy.x); |
130 |
| - vec2 g1 = vec2(gx.y,gy.y); |
131 |
| - vec2 g2 = vec2(gx.z,gy.z); |
132 |
| - |
133 |
| - // Radial decay with distance from each simplex corner |
134 |
| - vec3 w = 0.8 - vec3(dot(x0, x0), dot(x1, x1), dot(x2, x2)); |
135 |
| - w = max(w, 0.0); |
136 |
| - vec3 w2 = w * w; |
137 |
| - vec3 w4 = w2 * w2; |
138 |
| - |
139 |
| - // The value of the linear ramp from each of the corners |
140 |
| - vec3 gdotx = vec3(dot(g0, x0), dot(g1, x1), dot(g2, x2)); |
141 |
| - |
142 |
| - // Multiply by the radial decay and sum up the noise value |
143 |
| - float n = dot(w4, gdotx); |
144 |
| - |
145 |
| - // Compute the first order partial derivatives |
146 |
| - vec3 w3 = w2 * w; |
147 |
| - vec3 dw = -8.0 * w3 * gdotx; |
148 |
| - vec2 dn0 = w4.x * g0 + dw.x * x0; |
149 |
| - vec2 dn1 = w4.y * g1 + dw.y * x1; |
150 |
| - vec2 dn2 = w4.z * g2 + dw.z * x2; |
151 |
| - gradient = 10.9 * (dn0 + dn1 + dn2); |
152 |
| - |
153 |
| - // Compute the second order partial derivatives |
154 |
| - vec3 dg0, dg1, dg2; |
155 |
| - vec3 dw2 = 48.0 * w2 * gdotx; |
156 |
| - // d2n/dx2 and d2n/dy2 |
157 |
| - dg0.xy = dw2.x * x0 * x0 - 8.0 * w3.x * (2.0 * g0 * x0 + gdotx.x); |
158 |
| - dg1.xy = dw2.y * x1 * x1 - 8.0 * w3.y * (2.0 * g1 * x1 + gdotx.y); |
159 |
| - dg2.xy = dw2.z * x2 * x2 - 8.0 * w3.z * (2.0 * g2 * x2 + gdotx.z); |
160 |
| - // d2n/dxy |
161 |
| - dg0.z = dw2.x * x0.x * x0.y - 8.0 * w3.x * dot(g0, x0.yx); |
162 |
| - dg1.z = dw2.y * x1.x * x1.y - 8.0 * w3.y * dot(g1, x1.yx); |
163 |
| - dg2.z = dw2.z * x2.x * x2.y - 8.0 * w3.z * dot(g2, x2.yx); |
164 |
| - dg = 10.9 * (dg0 + dg1 + dg2); |
165 |
| - |
166 |
| - // Scale the return value to fit nicely into the range [-1,1] |
167 |
| - return 10.9 * n; |
| 72 | + out vec3 dg) { |
| 73 | + |
| 74 | + // Transform to simplex space (axis-aligned hexagonal grid) |
| 75 | + vec2 uv = vec2(x.x + x.y*0.5, x.y); |
| 76 | + |
| 77 | + // Determine which simplex we're in, with i0 being the "base" |
| 78 | + vec2 i0 = floor(uv); |
| 79 | + vec2 f0 = fract(uv); |
| 80 | + // o1 is the offset in simplex space to the second corner |
| 81 | + float cmp = step(f0.y, f0.x); |
| 82 | + vec2 o1 = vec2(cmp, 1.0-cmp); |
| 83 | + |
| 84 | + // Enumerate the remaining simplex corners |
| 85 | + vec2 i1 = i0 + o1; |
| 86 | + vec2 i2 = i0 + vec2(1.0, 1.0); |
| 87 | + |
| 88 | + // Transform corners back to texture space |
| 89 | + vec2 v0 = vec2(i0.x - i0.y * 0.5, i0.y); |
| 90 | + vec2 v1 = vec2(v0.x + o1.x - o1.y * 0.5, v0.y + o1.y); |
| 91 | + vec2 v2 = vec2(v0.x + 0.5, v0.y + 1.0); |
| 92 | + |
| 93 | + // Compute vectors from v to each of the simplex corners |
| 94 | + vec2 x0 = x - v0; |
| 95 | + vec2 x1 = x - v1; |
| 96 | + vec2 x2 = x - v2; |
| 97 | + |
| 98 | + vec3 iu, iv; |
| 99 | + vec3 xw, yw; |
| 100 | + |
| 101 | + // Wrap to periods, if desired |
| 102 | + if(any(greaterThan(period, vec2(0.0)))) { |
| 103 | + xw = vec3(v0.x, v1.x, v2.x); |
| 104 | + yw = vec3(v0.y, v1.y, v2.y); |
| 105 | + if(period.x > 0.0) |
| 106 | + xw = mod(vec3(v0.x, v1.x, v2.x), period.x); |
| 107 | + if(period.y > 0.0) |
| 108 | + yw = mod(vec3(v0.y, v1.y, v2.y), period.y); |
| 109 | + // Transform back to simplex space and fix rounding errors |
| 110 | + iu = floor(xw + 0.5*yw + 0.5); |
| 111 | + iv = floor(yw + 0.5); |
| 112 | + } else { // Shortcut if neither x nor y periods are specified |
| 113 | + iu = vec3(i0.x, i1.x, i2.x); |
| 114 | + iv = vec3(i0.y, i1.y, i2.y); |
| 115 | + } |
| 116 | + |
| 117 | + // Compute one pseudo-random hash value for each corner |
| 118 | + vec3 hash = mod(iu, 289.0); |
| 119 | + hash = mod((hash*51.0 + 2.0)*hash + iv, 289.0); |
| 120 | + hash = mod((hash*34.0 + 10.0)*hash, 289.0); |
| 121 | + |
| 122 | + // Pick a pseudo-random angle and add the desired rotation |
| 123 | + vec3 psi = hash * 0.07482 + alpha; |
| 124 | + vec3 gx = cos(psi); |
| 125 | + vec3 gy = sin(psi); |
| 126 | + |
| 127 | + // Reorganize for dot products below |
| 128 | + vec2 g0 = vec2(gx.x,gy.x); |
| 129 | + vec2 g1 = vec2(gx.y,gy.y); |
| 130 | + vec2 g2 = vec2(gx.z,gy.z); |
| 131 | + |
| 132 | + // Radial decay with distance from each simplex corner |
| 133 | + vec3 w = 0.8 - vec3(dot(x0, x0), dot(x1, x1), dot(x2, x2)); |
| 134 | + w = max(w, 0.0); |
| 135 | + vec3 w2 = w * w; |
| 136 | + vec3 w4 = w2 * w2; |
| 137 | + |
| 138 | + // The value of the linear ramp from each of the corners |
| 139 | + vec3 gdotx = vec3(dot(g0, x0), dot(g1, x1), dot(g2, x2)); |
| 140 | + |
| 141 | + // Multiply by the radial decay and sum up the noise value |
| 142 | + float n = dot(w4, gdotx); |
| 143 | + |
| 144 | + // Compute the first order partial derivatives |
| 145 | + vec3 w3 = w2 * w; |
| 146 | + vec3 dw = -8.0 * w3 * gdotx; |
| 147 | + vec2 dn0 = w4.x * g0 + dw.x * x0; |
| 148 | + vec2 dn1 = w4.y * g1 + dw.y * x1; |
| 149 | + vec2 dn2 = w4.z * g2 + dw.z * x2; |
| 150 | + gradient = 10.9 * (dn0 + dn1 + dn2); |
| 151 | + |
| 152 | + // Compute the second order partial derivatives |
| 153 | + vec3 dg0, dg1, dg2; |
| 154 | + vec3 dw2 = 48.0 * w2 * gdotx; |
| 155 | + // d2n/dx2 and d2n/dy2 |
| 156 | + dg0.xy = dw2.x * x0 * x0 - 8.0 * w3.x * (2.0 * g0 * x0 + gdotx.x); |
| 157 | + dg1.xy = dw2.y * x1 * x1 - 8.0 * w3.y * (2.0 * g1 * x1 + gdotx.y); |
| 158 | + dg2.xy = dw2.z * x2 * x2 - 8.0 * w3.z * (2.0 * g2 * x2 + gdotx.z); |
| 159 | + // d2n/dxy |
| 160 | + dg0.z = dw2.x * x0.x * x0.y - 8.0 * w3.x * dot(g0, x0.yx); |
| 161 | + dg1.z = dw2.y * x1.x * x1.y - 8.0 * w3.y * dot(g1, x1.yx); |
| 162 | + dg2.z = dw2.z * x2.x * x2.y - 8.0 * w3.z * dot(g2, x2.yx); |
| 163 | + dg = 10.9 * (dg0 + dg1 + dg2); |
| 164 | + |
| 165 | + // Scale the return value to fit nicely into the range [-1,1] |
| 166 | + return 10.9 * n; |
168 | 167 | }
|
0 commit comments