本文目的是以一種通俗簡單的方式介紹Ken Perlin的改進版柏林噪聲算法，講解代碼采用c#編寫，開源免費使用。如果你只是想看完整代碼，點擊這里。
柏林噪聲是一個非常強大算法，經常用于程序生成隨機內容，在游戲和其他像電影等多媒體領域廣泛應用。算法發(fā)明者Ken Perlin也因此算法獲得奧斯卡科技成果獎（靠算法拿奧斯卡也是沒誰了666）。本文將剖析他于2002年發(fā)表的改進版柏林噪聲算法。在游戲開發(fā)領域，柏林噪聲可以用于生成波形，起伏不平的材質或者紋理。例如，它能用于程序生成地形（例如使用柏林噪聲來生成我的世界（Minecraft）里的地形），火焰燃燒特效，水和云等等。柏林噪聲絕大部分應用在2維，3維層面上，但某種意義上也能拓展到4維。柏林噪聲在1維層面上可用于卷軸地形、模擬手繪線條等。
如果將柏林噪聲拓展到4維層面，以第4維，即w軸代表時間，就能利用柏林噪聲做動畫。例如，2D柏林噪聲可以通過插值生成地形，而3D柏林噪聲則可以模擬海平面上起伏的波浪。下面是柏林噪聲在不同維度的圖像以及在游戲中的應用場景。
正如圖所示，柏林噪聲算法可以用來模擬許多自然中的噪聲現(xiàn)象。接下來讓我們從數理上分析算法的實現(xiàn)原理。
基本原理
注意：事先聲明，本節(jié)內容大多源于this wonderful article by Matt Zucker，不過該篇文章內容也是建立在1980年所發(fā)明的柏林噪聲算法基礎上的。本文我將使用2002年發(fā)明的改進版柏林噪聲算法。因此，我的算法版本跟Zucker的版本會有些不同。
讓我們從最基本的柏林噪聲函數看起：public double perlin(double x, double y, double z);
函數接收x,y,z三個坐標分量作為輸入，并返回0.0~1.0的double值作為輸出。那我們應該怎么處理輸入值？首先，我們取3個輸入值x,y,z的小數點部分，就可以表示為單元空間里的一個點了。為了方便講解，我們將問題降維到2維空間來討論（原理是一樣的），下圖是該點在2維空間上的表示：
圖1：小藍點代表輸入值在單元正方形里的空間坐標，其他4個點則是單元正方形的各頂點
接著，我們給4個頂點（在3維空間則是8個頂點）各自生成一個偽隨機的梯度向量。梯度向量代表該頂點相對單元正方形內某點的影響是正向還是反向的（向量指向方向為正向，相反方向為反向）。而偽隨機是指，對于任意組相同的輸入，必定得到相同的輸出。因此，雖然每個頂點生成的梯度向量看似隨機，實際上并不是。這也保證了在梯度向量在生成函數不變的情況下，每個坐標的梯度向量都是確定不變的。
舉個例子來理解偽隨機，比如我們從圓周率π（3.14159...）的小數部分中隨機抽取某一位數字，結果看似隨機，但如果抽取小數點后1位，結果必定為1；抽取小數點后2位，結果必定為4。
圖2：各頂點上的梯度向量隨機選取結果
請注意，上圖所示的梯度向量并不是完全準確的。在本文所介紹的改進版柏林噪聲中，這些梯度向量并不是完全隨機的，而是由12條單位正方體（3維）的中心點到各條邊中點的向量組成：(1,1,0),(-1,1,0),(1,-1,0),(-1,-1,0), (1,0,1),(-1,0,1),(1,0,-1),(-1,0,-1), (0,1,1),(0,-1,1),(0,1,-1),(0,-1,-1)
采用這些特殊梯度向量的原因在Ken Perlin's SIGGRAPH 2002 paper: Improving Noise這篇文章里有具體講解。
注意：許多介紹柏林噪聲算法的文章都是根據最初版柏林噪聲算法來講解的，預定義的梯度表不是本文所說的這12個向量。如圖2所示的梯度向量就是最初版算法所隨機出來的梯度向量，不過這兩種算法的原理都是一樣的。
接著，我們需要求出另外4個向量（在3維空間則是8個），它們分別從各頂點指向輸入點（藍色點）。下面有個2維空間下的例子：
圖3:各個距離向量
接著，對每個頂點的梯度向量和距離向量做點積運算，我們就可以得出每個頂點的影響值：grad.x * dist.x + grad.y * dist.y + grad.z * dist.z
這正是算法所需要的值，點積運算為兩向量長度之積，再乘以兩向量夾角余弦：
dot(vec1,vec2) = cos(angle(vec1,vec2)) * vec1.length * vec2.length
換句話說，如果兩向量指向同一方向，點積結果為：
1 * vec1.length * vec2.length
如果兩向量指向相反方向，則點積結果為：
-1 * vec1.length * vec2.length
如果兩向量互相垂直，則點積結果為0。
0 * vec1.length * vec2.length
點積也可以理解為向量a在向量b上的投影，當距離向量在梯度向量上的投影為同方向，點積結果為正數；當距離向量在梯度向量上的投影為反方向，點積結果為負數。因此，通過兩向量點積，我們就知道該頂點的影響值是正還是負的。不難看出，頂點的梯度向量直接決定了這一點。下面通過一副彩色圖，直觀地看下各頂點的影響值：
圖4：2D柏林噪聲的影響值
下一步，我們需要對4個頂點的影響值做插值，求得加權平均值（在3維空間則是8個）。算法非常簡單（2維空間下的解法）：
// Below are 4 influence values in the arrangement:// [g1] | [g2]// -----------// [g3] | [g4]int g1, g2, g3, g4;int u, v; // These coordinates are the location of the input coordinate in its unit square. // For example a value of (0.5,0.5) is in the exact center of its unit square.int x1 = lerp(g1,g2,u);int x2 = lerp(g3,h4,u);int average = lerp(x1,x2,v);
至此，整個柏林噪聲算法還剩下最后一塊拼圖了：如果直接使用上述代碼，由于是采用lerp線性插值計算得出的值，雖然運行效率高，但噪聲效果不好，看起來會不自然。我們需要采用一種更為平滑，非線性的插值函數：fade函數，通常也被稱為ease curve(也作為緩動函數在游戲中廣泛使用)：
圖5：ease curve
ease curve的值會用來計算前面代碼里的u和v，這樣插值變化不再是單調的線性變化，而是這樣一個過程：初始變化慢，中間變化快，結尾變化又慢下來（也就是在當數值趨近于整數時，變化變慢）。這個用于改善柏林噪聲算法的fade函數可以表示為以下數學形式：
基本上，這就是整個柏林噪聲算法的原理了！搞清了算法的各個實現(xiàn)關鍵步驟后，現(xiàn)在讓我們著手把代碼實現(xiàn)出來。
代碼實現(xiàn)
在本節(jié)開始前我需要重申一遍，代碼實現(xiàn)是C#版本。相比Ken Perlin的Java版本實現(xiàn)做了小小的改動，主要是增加了代碼的整潔性和可讀性，支持噪聲重復（瓦片重復）特性。代碼完全開源，可免費使用（考慮到這畢竟不是我原創(chuàng)發(fā)明的算法 - Ken Perlin才是?。?/p>
準備工作
第一步，我們需要先聲明一個排列表（permutation table），或者直接縮寫為p[]數組就行了。數組長度為256，分別隨機、無重復地存放了0-255這些數值。為了避免緩存溢出，我們再重復填充一次數組的值，所以數組最終長度為512：
private static readonly int[] permutation = { 151,160,137,91,90,15, // Hash lookup table as defined by Ken Perlin. This is a randomly
131,13,201,95,96,53,194,233,7,225,140,36,103,30,69,142,8,99,37,240,21,10,23, // arranged array of all numbers from 0-255 inclusive.
190, 6,148,247,120,234,75,0,26,197,62,94,252,219,203,117,35,11,32,57,177,33, 88,237,149,56,87,174,20,125,136,171,168, 68,175,74,165,71,134,139,48,27,166, 77,146,158,231,83,111,229,122,60,211,133,230,220,105,92,41,55,46,245,40,244, 102,143,54, 65,25,63,161, 1,216,80,73,209,76,132,187,208, 89,18,169,200,196, 135,130,116,188,159,86,164,100,109,198,173,186, 3,64,52,217,226,250,124,123, 5,202,38,147,118,126,255,82,85,212,207,206,59,227,47,16,58,17,182,189,28,42, 223,183,170,213,119,248,152, 2,44,154,163, 70,221,153,101,155,167, 43,172,9, 129,22,39,253, 19,98,108,110,79,113,224,232,178,185, 112,104,218,246,97,228, 251,34,242,193,238,210,144,12,191,179,162,241, 81,51,145,235,249,14,239,107, 49,192,214, 31,181,199,106,157,184, 84,204,176,115,121,50,45,127, 4,150,254, 138,236,205,93,222,114,67,29,24,72,243,141,128,195,78,66,215,61,156,180};private static readonly int[] p; // Doubled permutation to avoid overflowstatic Perlin() {
p = new int[512]; for(int x=0;x<512;x++) {
p[x] = permutation[x%256];
}
}p[]數組會在算法后續(xù)的哈希計算中使用到，用于確定一組輸入最終挑選哪個梯度向量（從前面所列出的12個梯度向量中挑選）。后續(xù)代碼會詳細展示p[]數組的用法。
接著，我們開始編寫柏林噪聲函數：
public double perlin(double x, double y, double z) { if(repeat > 0) { // If we have any repeat on, change the coordinates to their "local" repetitions
x = x%repeat;
y = y%repeat;
z = z%repeat;
}
int xi = (int)x & 255; // Calculate the "unit cube" that the point asked will be located in
int yi = (int)y & 255; // The left bound is ( |_x_|,|_y_|,|_z_| ) and the right bound is that
int zi = (int)z & 255; // plus 1. Next we calculate the location (from 0.0 to 1.0) in that cube.
double xf = x-(int)x; double yf = y-(int)y; double zf = z-(int)z; // ...}上面的代碼很直觀。首先，對輸入坐標使用求余運算符%，求出[0,repeat)范圍內的余數。緊接著聲明xi, yi, zi三個變量。它們代表了輸入坐標落在了哪個單元正方形里。我們還要限制坐標在[0,255]這個范圍內，避免訪問數組p[]時，出現(xiàn)數組越界錯誤。這也產生了一個副作用：柏林噪聲每隔256個整數就會再次重復。但這不是太大的問題，因為算法不僅能處理整數，還能處理小數。最后，我們通過xf, yf, zf三個變量（也就是x,y,z的小數部分值），確定了輸入坐標在單元正方形里的空間位置（就是前面所示的小藍點）。
Fade函數
現(xiàn)在我們需要用代碼表示前面所提到的fade函數（圖5）。正如上文所提，函數的數學表示：
代碼實現(xiàn)如下：
public static double fade(double t) { // Fade function as defined by Ken Perlin. This eases coordinate values
// so that they will ease towards integral values. This ends up smoothing
// the final output.
return t * t * t * (t * (t * 6 - 15) + 10); // 6t^5 - 15t^4 + 10t^3}public double perlin(double x, double y, double z) { // ...
double u = fade(xf); double v = fade(yf); double w = fade(zf); // ...}代碼所計算得出的u / v / w變量將在后面的插值計算中使用到。
哈希函數
柏林噪聲哈希函數用于給每組輸入計算返回一個唯一、確定值。哈希函數在維基百科的定義如下：
哈希函數是一種從任何一種數據中創(chuàng)建小的數字“指紋”的方法，輸入數據有任何細微的不同，都會令輸出結果完全不一樣
下面代碼就是柏林噪聲算法所使用的哈希函數。它使用了早前我們聲明的p[]數組：
public double perlin(double x, double y, double z) {
// ...
int aaa, aba, aab, abb, baa, bba, bab, bbb;
aaa = p[p[p[ xi ]+ yi ]+ zi ];
aba = p[p[p[ xi ]+inc(yi)]+ zi ];
aab = p[p[p[ xi ]+ yi ]+inc(zi)];
abb = p[p[p[ xi ]+inc(yi)]+inc(zi)];
baa = p[p[p[inc(xi)]+ yi ]+ zi ];
bba = p[p[p[inc(xi)]+inc(yi)]+ zi ];
bab = p[p[p[inc(xi)]+ yi ]+inc(zi)];
bbb = p[p[p[inc(xi)]+inc(yi)]+inc(zi)];
// ...
}
public int inc(int num) {
num++;
if (repeat > 0) num %= repeat;
return num;
}代碼的哈希函數，對包圍著輸入坐標（小藍點）的周圍8個單元正方形的索引坐標進行了哈希計算。inc()函數用于將輸入值增加1，同時保證范圍在[0,repeat)內。如果不需要噪聲重復，inc()函數可以簡化成單純將輸入值增加1。由于哈希結果值是從p[]數組中得到的，所以哈希函數的返回值范圍限定在[0,255]內。
梯度函數
我時常認為Ken Perlin的最初版算法里的grad()函數寫法過于復雜，令人費解。我們只要明白grad()函數的作用在于計算隨機選取的梯度向量以及頂點位置向量的點積。Ken Perlin巧妙地使用了位翻轉(bit-flipping)技巧來實現(xiàn)：
public static double grad(int hash, double x, double y, double z) { int h = hash & 15; // Take the hashed value and take the first 4 bits of it (15 == 0b1111)
double u = h < 8 /* 0b1000 */ ? x : y; // If the most significant bit (MSB) of the hash is 0 then set u = x. Otherwise y.
double v; // In Ken Perlin's original implementation this was another conditional operator (?:). I
// expanded it for readability.
if(h < 4 /* 0b0100 */) // If the first and second significant bits are 0 set v = y
v = y; else if(h == 12 /* 0b1100 */ || h == 14 /* 0b1110*/) // If the first and second significant bits are 1 set v = x
v = x; else // If the first and second significant bits are not equal (0/1, 1/0) set v = z
v = z;
return ((h&1) == 0 ? u : -u)+((h&2) == 0 ? v : -v); // Use the last 2 bits to decide if u and v are positive or negative. Then return their addition.}下面代碼則是以另一種令人容易理解的方式完成了這個任務（而且在很多語言版本的運行效率都優(yōu)于前面一種）：
// Source: http://riven8192.blogspot.com/2010/08/calculate-perlinnoise-twice-as-fast.htmlpublic static double grad(int hash, double x, double y, double z){ switch(hash & 0xF)
{ case 0x0: return x + y; case 0x1: return -x + y; case 0x2: return x - y; case 0x3: return -x - y; case 0x4: return x + z; case 0x5: return -x + z; case 0x6: return x - z; case 0x7: return -x - z; case 0x8: return y + z; case 0x9: return -y + z; case 0xA: return y - z; case 0xB: return -本文原文出自Understanding Perlin Noise，翻譯不易，如果你已經看到這里了，順手點個贊唄！
http://www.cnblogs.com/leoin2012/p/7218033.html
延伸閱讀
- 人生感悟之項目篇 2023-04-15
- 懂哥網 2022-09-20
- 多敲代碼，不斷試錯 2022-09-08
- 加油吧！代碼人 2021-05-27
- 代碼魅力之所在 2021-05-19
- 每行代碼都是進步 2021-05-19
- 我與代碼打交道 ——王*松 2020-06-19
- 寫代碼是一個痛并快樂著的過程——李*文 2019-11-27
- 我對如何寫代碼有了一定的新的認識——胡*翔 2019-03-12
- 開啟奇妙代碼之旅 2018-12-05
人生感悟之項目篇