noise2d.cpp
Engine/source/util/noise2d.cpp
Detailed Description
Public Functions
clamp(F32 f, F32 m)
1 2//----------------------------------------------------------------------------- 3// Copyright (c) 2012 GarageGames, LLC 4// 5// Permission is hereby granted, free of charge, to any person obtaining a copy 6// of this software and associated documentation files (the "Software"), to 7// deal in the Software without restriction, including without limitation the 8// rights to use, copy, modify, merge, publish, distribute, sublicense, and/or 9// sell copies of the Software, and to permit persons to whom the Software is 10// furnished to do so, subject to the following conditions: 11// 12// The above copyright notice and this permission notice shall be included in 13// all copies or substantial portions of the Software. 14// 15// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR 16// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 17// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE 18// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 19// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING 20// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS 21// IN THE SOFTWARE. 22//----------------------------------------------------------------------------- 23 24#include "util/noise2d.h" 25#include "core/util/tVector.h" 26 27//-------------------------------------- 28Noise2D::Noise2D() 29{ 30 dMemset(mPermutation, 0, sizeof(mPermutation)); 31 dMemset(mGradient, 0, sizeof(mGradient)); 32 mSeed = 0; 33} 34 35Noise2D::~Noise2D() 36{ 37} 38 39 40//-------------------------------------- 41void Noise2D::normalize(F32 v[2]) 42{ 43 F32 s; 44 45 s = mSqrt(v[0] * v[0] + v[1] * v[1]); 46 v[0] = v[0] / s; 47 v[1] = v[1] / s; 48} 49 50 51//-------------------------------------- 52void Noise2D::setSeed(U32 seed) 53{ 54 if (mSeed == seed) 55 return; 56 mSeed = seed; 57 mRandom.setSeed(mSeed); 58 59 S32 i, j, k; 60 61 for (i = 0 ; i < SIZE ; i++) { 62 mPermutation[i] = i; 63 64 for (j = 0 ; j < 2 ; j++) 65 mGradient[i][j] = mRandom.randF( -1.0f, 1.0f ); 66 normalize(mGradient[i]); 67 } 68 69 while (--i) { 70 k = mPermutation[i]; 71 j = mRandom.randI(0, SIZE-1); 72 mPermutation[i] = mPermutation[j]; 73 mPermutation[j] = k; 74 } 75 76 // extend the size of the arrays x2 to get rid of a bunch of MODs 77 // we'd have to do later in the code 78 for (i = 0 ; i < SIZE + 2 ; i++) { 79 mPermutation[SIZE + i] = mPermutation[i]; 80 for (j = 0 ; j < 2 ; j++) 81 mGradient[SIZE + i][j] = mGradient[i][j]; 82 } 83} 84 85 86//-------------------------------------- 87U32 Noise2D::getSeed() 88{ 89 return mSeed; 90} 91 92 93inline F32 Noise2D::lerp(F32 t, F32 a, F32 b) 94{ 95 return a + t * (b - a); 96} 97 98 99inline F32 Noise2D::curve(F32 t) 100{ 101 return t * t * (3.0f - 2.0f * t); 102} 103 104 105inline F32 clamp(F32 f, F32 m) 106{ 107 while (f > m) 108 f -= m; 109 while (f < 0.0f) 110 f += m; 111 return f; 112} 113 114 115//-------------------------------------- 116void Noise2D::fBm( Vector<F32> *dst, U32 size, U32 interval, F32 h, F32 octaves ) 117{ 118 interval = getMin(U32(128), getMax(U32(1), interval)); 119 F32 H = getMin(1.0f, getMax(0.0f, h)); 120 octaves = getMin(5.0f, getMax(1.0f, octaves)); 121 F32 lacunarity = 2.0f; 122 123 F32 exponent_array[32]; 124 125 U32 shift = getBinLog2( size ); 126 127 // precompute and store spectral weights 128 // seize required memory for exponent_array 129 F32 frequency = 1.0; 130 for (U32 i=0; i<=octaves; i++) 131 { 132 // compute weight for each frequency 133 exponent_array[i] = mPow( frequency, -H ); 134 frequency *= lacunarity; 135 } 136 137 // initialize dst 138 for (S32 k=0; k < (size*size); k++) 139 (*dst)[k] = 0.0f; 140 141 F32 scale = 1.0f / (F32)size * interval; 142 for (S32 o=0; o<octaves; o++) 143 { 144 F32 exp = exponent_array[o]; 145 for (S32 y=0; y<size; y++) 146 { 147 F32 fy = (F32)y * scale; 148 for (S32 x=0; x<size; x++) 149 { 150 F32 fx = (F32)x * scale; 151 F32 noise = getValue(fx, fy, interval); 152 (*dst)[x + (y << shift)] += noise * exp; 153 } 154 } 155 scale *= lacunarity; 156 interval = (U32)(interval * lacunarity); 157 } 158} 159 160 161//-------------------------------------- 162void Noise2D::rigidMultiFractal(Vector<F32> *dst, Vector<F32> *sig, U32 size, U32 interval, F32 h, F32 octaves) 163{ 164 interval = getMin(U32(128), getMax(U32(1), interval)); 165 F32 H = getMin(1.0f, getMax(0.0f, h)); 166 octaves = getMin(5.0f, getMax(1.0f, octaves)); 167 F32 lacunarity = 2.0f; 168 F32 offset = 1.0f; 169 F32 gain = 2.0f; 170 171 U32 shift = getBinLog2( size ); 172 173 F32 exponent_array[32]; 174 175 // precompute and store spectral weights 176 // seize required memory for exponent_array 177 F32 frequency = 1.0; 178 for (U32 i=0; i<=octaves; i++) 179 { 180 // compute weight for each frequency 181 exponent_array[i] = mPow( frequency, -H ); 182 frequency *= lacunarity; 183 } 184 185 F32 scale = 1.0f / (F32)size * interval; 186 187 //-------------------------------------- 188 // compute first octave 189 for (S32 y=0; y<size; y++) 190 { 191 F32 fy = (F32)y * scale; 192 for (S32 x=0; x<size; x++) 193 { 194 F32 fx = (F32)x * scale; 195 196 F32 signal = mFabs(getValue(fx,fy,interval)); // get absolute value of signal (this creates the ridges) 197 //signal = mSqrt(signal); 198 signal = offset - signal; // invert and translate (note that "offset" should be ~= 1.0) 199 signal *= signal + 0.1; // square the signal, to increase "sharpness" of ridges 200 201 // assign initial values 202 (*dst)[x + (y << shift)] = signal; 203 (*sig)[x + (y << shift)] = signal; 204 } 205 } 206 207 //-------------------------------------- 208 // compute remaining octaves 209 for (S32 o=1; o<octaves; o++) 210 { 211 // increase the frequency 212 scale *= lacunarity; 213 interval = (U32)(interval * lacunarity); 214 F32 exp = exponent_array[o]; 215 216 for (S32 y=0; y<size; y++) 217 { 218 F32 fy = (F32)y * scale; 219 for (S32 x=0; x<size; x++) 220 { 221 F32 fx = (F32)x * scale; 222 U32 index = x + (y << shift); 223 F32 result = (*dst)[index]; 224 F32 signal = (*sig)[index]; 225 226 // weight successive contributions by previous signal 227 F32 weight = mClampF(signal * gain, 0.0f, 1.0f); 228 229 signal = mFabs(getValue( fx, fy, interval )); 230 231 signal = offset - signal; 232 signal *= signal + 0.2; 233 // weight the contribution 234 signal *= weight; 235 result += signal * exp; 236 237 (*dst)[index] = result; 238 (*sig)[index] = signal; 239 } 240 } 241 } 242 243 for (S32 k=0; k < (size*size); k++) 244 (*dst)[k] = ((*dst)[k]-1.0f)/2.0f; 245} 246 247bool Noise2D::erodeHydraulic( Vector<F32> *src, Vector<F32> *dst, U32 iterations, U32 size ) 248{ 249 // early out if there is nothing to do 250 if (iterations == 0 ) 251 { 252 *dst = *src; 253 return true; 254 } 255 256 F32 fmin, fmax; 257 getMinMax( src, &fmin, &fmax, size); 258 259 U32 shift = getBinLog2( size ); 260 U32 mask = size - 1; 261 262 263// currently using SCRATCH_3 for debugging -- Rick 264 Vector<F32> scratch = *src; 265 U32 *o = (U32*)scratch.address(); 266 Vector<F32> a = *src; 267 Vector<F32> b = *src; 268 Vector<F32> c = *src; 269 270 for (S32 k=0; k < (size*size); k++) 271 c[k] = 0.0f; 272 273 for (S32 i=0; i<iterations; i++) 274 { 275 b = a; 276 277 for (S32 y=0; y<size; y++) 278 { 279 for (S32 x=0; x<size; x++) 280 { 281 U32 srcOffset = (x + (y << shift)); 282 F32 height = a[srcOffset]; 283 o[srcOffset] = srcOffset; 284 for (S32 y1=y-1; y1 <= y+1; y1++) 285 { 286 F32 maxDelta = 0.0f; 287 S32 ywrap = (y1 & mask); 288 for (S32 x1=x-1; x1 <= x+1; x1++) 289 { 290 if (x1 != x && y1 != y) 291 { 292 U32 adjOffset = ((x1 & mask) + (ywrap << shift)); 293 F32 &adjHeight = a[adjOffset]; 294 F32 delta = height - adjHeight; 295 if (x1 != x || y1 != y) 296 delta *= 1.414213562f; // compensate for diagonals 297 if (delta > maxDelta) 298 { 299 maxDelta = delta; 300 o[srcOffset] = adjOffset; 301 } 302 } 303 } 304 } 305 } 306 } 307 for (S32 j=0; j < (size*size); j++) 308 { 309 F32 &s = a[j]; 310 F32 &d = b[ o[j] ]; 311 F32 delta = s - d; 312 if (delta > 0.0f) 313 { 314 F32 alt = (s-fmin) / (fmax-fmin); 315 F32 amt = delta * (0.1f * (1.0f-alt)); 316 s -= amt; 317 d += amt; 318 } 319 } 320// debug only 321 for (S32 k=0; k < (size*size); k++) 322 c[k] += b[k] - a[k]; 323 324 Vector<F32> tmp = a; 325 a = b; 326 b = tmp; 327 } 328 *dst = b; 329 //*dst = *c; 330 331 return true; 332} 333 334 335 336bool Noise2D::erodeThermal(Vector<F32> *src, Vector<F32> *dst, F32 slope, F32 materialLoss, U32 iterations, U32 size, U32 squareSize, F32 maxHeight ) 337{ 338 // early out if there is nothing to do 339 if (iterations == 0 ) 340 { 341 *dst = *src; 342 return true; 343 } 344 345 346 F32 fmin, fmax; 347 getMinMax(src, &fmin, &fmax, size); 348 349 Vector<F32> a = *src; 350 // Heightfield *b = getScratch(1); 351 Vector<F32> r; 352 r.setSize( size * size ); 353 //dMemset( r.address(), 0, r.memSize() ); 354 355 F32 conservation = 1.0f - mClampF(materialLoss, 0.0f, 100.0f)/100.0f; 356 slope = mClampF(conservation, 0.0f, 89.0f); // clamp to 0-89 degrees 357 358 F32 talusConst = mTan(mDegToRad(slope)) * squareSize; // in world units 359 talusConst = talusConst * (fmax-fmin) / maxHeight; // scale to current height units 360 F32 p = 0.1f; 361 362 U32 mask = size - 1; 363 U32 shift = getBinLog2( size ); 364 365 for (U32 i=0; i<iterations; i++) 366 { 367 // clear out the rubble accumulation field 368 dMemset( r.address(), 0, r.memSize() ); 369 370 for (S32 y=0; y<size; y++) 371 { 372 for (S32 x=0; x<size; x++) 373 { 374 F32 *height = &a[ x + ( y << shift )]; 375 F32 *dstHeight = &r[ x + ( y << shift )]; 376 377 // for each height look at the immediate surrounding heights 378 // if any are higher than talusConst erode on me 379 for (S32 y1=y-1; y1 <= y+1; y1++) 380 { 381 S32 ywrap = (y1 & mask); 382 for (S32 x1=x-1; x1 <= x+1; x1++) 383 { 384 if (x1 != x && y1 != y) 385 { 386 S32 adjOffset = ((x1 & mask) + (ywrap << shift)); 387 F32 adjHeight = a[adjOffset]; 388 F32 delta = adjHeight - *height; 389 if (delta > talusConst) 390 { 391 F32 rubble = p * (delta - talusConst); 392 r[adjOffset] -= rubble; 393 *dstHeight += rubble * conservation; 394 } 395 } 396 } 397 } 398 } 399 } 400 for (S32 k=0; k < (size*size); k++) 401 a[k] += r[k]; 402 } 403 *dst = a; 404 return true; 405} 406 407void Noise2D::getMinMax( Vector<F32> *src, F32 *fmin, F32 *fmax, U32 size ) 408{ 409 if (!src) 410 return; 411 412 F32 *p = (*src).address(); 413 *fmin = *p; 414 *fmax = *p; 415 for (S32 i=0; i < (size*size); i++, p++) 416 { 417 if (*fmin > *p) *fmin = *p; 418 if (*fmax < *p) *fmax = *p; 419 } 420} 421 422//-------------------------------------- 423F32 Noise2D::turbulence(F32 x, F32 y, F32 freq) 424{ 425 F32 t, x2, y2; 426 427 for ( t = 0.0f ; freq >= 3.0f ; freq /= 2.0f) 428 { 429 x2 = freq * x; 430 y2 = freq * y; 431 t += mFabs(getValue(x2, y2, (S32)freq)) / freq; 432 } 433 return t; 434} 435 436 437//-------------------------------------- 438inline void Noise2D::setup(F32 t, S32 &b0, S32 &b1, F32 &r0, F32 &r1) 439{ 440 // find the bounding integers of u 441 b0 = S32(t) & SIZE_MASK; 442 b1 = (b0+1) & SIZE_MASK; 443 444 // seperate the fractional components 445 r0 = t - (S32)t; 446 r1 = r0 - 1.0f; 447} 448 449inline F32 Noise2D::dot(const F32 *q, F32 rx, F32 ry) 450{ 451 return (rx * q[0] + ry * q[1] ); 452} 453 454 455 456//-------------------------------------- 457F32 Noise2D::getValue(F32 x, F32 y, S32 interval) 458{ 459 S32 bx0, bx1, by0, by1; 460 F32 rx0, rx1, ry0, ry1; 461 462 // Imagine having a square of the type 463 // p0---p1 Where p0 = (bx0, by0) +----> U 464 // |(u,v)| p1 = (bx1, by0) | 465 // | | p2 = (bx0, by1) | Coordinate System 466 // p2---p3 p3 = (bx1, by1) V 467 // The u, v point in 2D texture space is bounded by this rectangle. 468 469 // Goal, determine the scalar at the points p0, p1, p2, p3. 470 // Then the scalar of the point (u, v) will be found by linear interpolation. 471 472 // First step: Get the 2D coordinates of the points p0, p1, p2, p3. 473 // We also need vectors pointing from each point in the square above and 474 // ending at the (u,v) coordinate located inside the square. 475 // The vector (rx0, ry0) goes from P0 to the (u,v) coordinate. 476 // The vector (rx1, ry0) goes from P1 to the (u,v) coordinate. 477 // The vector (rx0, ry1) goes from P2 to the (u,v) coordinate. 478 // The vector (rx1, ry1) goes from P3 to the (u,v) coordinate. 479 480 setup(x, bx0, bx1, rx0, rx1); 481 setup(y, by0, by1, ry0, ry1); 482 483 // Make sure the box corners fall within the interval 484 // so that the final output will wrap on itself 485 bx0 = bx0 % interval; 486 bx1 = bx1 % interval; 487 by0 = by0 % interval; 488 by1 = by1 % interval; 489 490 S32 i = mPermutation[ bx0 ]; 491 S32 j = mPermutation[ bx1 ]; 492 493 S32 b00 = mPermutation[ i + by0 ]; 494 S32 b10 = mPermutation[ j + by0 ]; 495 S32 b01 = mPermutation[ i + by1 ]; 496 S32 b11 = mPermutation[ j + by1 ]; 497 498 // Next, calculate the dropoff component about the point p0. 499 F32 sx = curve(rx0); 500 F32 sy = curve(ry0); 501 502 // Now, for each point in the square shown above, calculate the dot 503 // product of the gradiant vector and the vector going from each square 504 // corner point to the (u,v) point inside the square. 505 F32 u = dot(mGradient[ b00 ], rx0,ry0); 506 F32 v = dot(mGradient[ b10 ], rx1,ry0); 507 508 // Interpolation along the X axis. 509 F32 a = lerp(sx, u, v); 510 511 u = dot(mGradient[ b01 ], rx0,ry1); 512 v = dot(mGradient[ b11 ], rx1,ry1); 513 514 // Interpolation along the Y axis. 515 F32 b = lerp(sx, u, v); 516 517 // Final Interpolation 518 return lerp(sy, a, b); 519} 520 521 522