首先说明为什么这篇博客叫这个题目，我刚刚开始学习Perlin Noise是从知乎上的一篇文章入门的，作者的题目是不只是噪声，我觉得很有韵味，就借鉴过来。这是链接：https://zhuanlan.zhihu.com/p/22337544
一. 背景
Perlin Noise，译作柏林噪声，是指Ken Perlin发明的噪声算法。Ken Perlin早在1983年就提出了Perlin noise，当时他正在参与制作迪士尼的动画电影《电子世界争霸战》（TRON），但是他不满足于当时计算机产生的那种非常不自然的纹理效果，因此提出了Perlin噪声。随后，他在1984年的SIGGRAPH Course上做了名为Advanced Image Synthesis1的课程演讲，并在SIGGRAPH 1985上发表了他的论文[Perlin K. An image synthesizer[J]. ACM Siggraph Computer Graphics, 1985, 19(3): 287-296.](Perlin K. An image synthesizer[J]. ACM Siggraph Computer Graphics, 1985, 19(3): 287-296.)。由于Perlin噪声的算法简单，被迅速应用到各种商业软件中。后来Perlin继续研究程序纹理的生成，并和他的一名学生又在SIGGRAPH 1989上发表了一篇文章[Perlin K, Hoffert E M. Hypertexture[C]//ACM SIGGRAPH Computer Graphics. ACM, 1989, 23(3): 253-262](Perlin K, Hoffert E M. Hypertexture[C]//ACM SIGGRAPH Computer Graphics. ACM, 1989, 23(3): 253-262)，提出了超级纹理（hypertexture）。他们使用噪声+fbm+ray marching实现了各种有趣的效果。到1990年，已经有大量公司在他们的产品中使用了Perlin噪声。在1999年的GDCHardCore大会上，Ken Perlin做了名为Making Noise的演讲[Perlin K. Making noise[C]//Proc. of the Game Developer Conference. 1999.](Perlin K. Making noise[C]//Proc. of the Game Developer Conference. 1999.)，系统地介绍了Perlin噪声的发展、实现细节和应用。在2002年他对原有的柏林噪声算法做了改进，并发表了论文ImprovingNoise。
Prelin由于出色的工作，获得了奥斯卡科技成果奖（奥斯卡是电影业的诺贝尔奖，靠着一个算法单枪匹马获得这个奖励，Perlin确实厉害）。
二. 分类
看了很多博客，博主都没有把分形噪声和Simplex噪声归为柏林噪声。但是我觉得这样不对，因为分形噪声和Simplex都是柏林先后提出来的，所以我把这2个和原始的柏林噪声一起归为柏林噪声，在后面都做介绍。事实上，维基百科也是这样做的，将这3个噪声全部归为柏林噪声。https://zh.wikipedia.org/wiki/Perlin%E5%99%AA%E5%A3%B0
根据wiki，由程序产生噪声的方法大致可以分为两类：
| 类别 | 名称 |
|---|---|
| 基于晶格的方法（Lattice based） | 又可细分为两种： 第一种是梯度噪声（Gradient noise），包括Perlin噪声， Simplex噪声，Wavelet噪声等； 第二种是Value噪声（Value noise）。 |
| 基于点的方法（Point based） | Worley噪声 |
从这个表格中，我们可以看到柏林噪声的产生是基于晶格的方法，从后面的介绍中相信一定会对”晶格”这个名字有一个深刻的理解。
下面贴出几张图片，先对这柏林噪声中的3种形式有一个直观的印象。
随机噪声
原始柏林噪声
分形噪声
Simplex噪声
从这4张图片中我们可以大概看出不同噪声的一些端倪，随机噪声就像黑白电视一样，杂乱无章。原始的柏林噪声中包含有随机的成分，但是相互之间又有关联，像素点之间是平滑的进行变化，而不是像白噪声那么尖锐。分形噪声就像我们在计算机图形学中学的分形树，简单的细节不断的重复，最后生成一个复杂的图形。
下面对白噪声，原始的柏林噪声，分形噪声，Simplex噪声的原理一一解释。
三. 白噪声
噪声的基础是随机数，如果给屏幕上的每个像素点赋予一个0和1之间的随机数来表示像素点的亮度，就会得到一幅杂乱无章的图片，就像我们小时候家里的电视在下雨天突然没有信号后屏幕上呈现的图像。从程序员的角度来看白噪声，我们无非就是引入一些随机变量到程序，给定一个输入（比如种子），程序给出一个输出。在C中，rand.h中对rand函数进行了声明，我们只要在自己的程序中引入这个头文件就行。不过这里生成的随机数是伪随机数，我们可以采用时间来作为种子。C++中也有实现了随机数的库，和C类似，毕竟是C的延伸和拓展。Python的第三方模块random对随机数的实现更加强大，因为在很多机器学习方法中需要用的随机数。U3D中通过Ramdom.Range()l来实现随机数的生成。基本的程序框架如下:(以Python为例)
import random
random.seed()
random.random()
我们从上面的代码生成随机数的过程中可以看出，这样产生的噪声非常不自然，听起来不自然，看起来也不好看。这学期学了信号与系统，不妨从信号的角度来考虑这个问题，一幅随机噪声的图像的功率谱密度在整个频域内均匀分布， 所有频率具有相同能量密度。也就是说，白噪声是指在较宽的频率范围内，各等带宽的频带所含的噪声能量相等。白噪声就像全彩色白颜色的光频。
白噪声虽然很好模拟，但是自身具有不可修复的缺点。观察现实生活中的自然噪声，它们不会长成上面这个样子。例如木头纹理、山脉起伏，它们的形状大多是趋于分形状（fractal）的，即包含了不同程度的细节。比如地形，它有起伏很大的山脉，也有起伏稍小的山丘，也有细节非常多的石子等，这些不同程度的细节共同组成了一个自然的地形表面。这些随机的成分并不是完全独立的，他们之间有一定的关联。采用微积分里面的说法，就是像素点之间变化是连续的，并没有发生跳变。很显然，白噪声没有做到这一点。
四. 原始的柏林噪声
从上一部分的讨论中，我们可以看到白噪声的缺点，原始柏林噪声的提出正是针对这一缺点。在大二学的数值分析中，我们学了很多插值方法，用来平滑变化趋势，这种方法在上学期X老师教的数字图像处理中处处可见。原始柏林噪声的原理就是基于插值的方法，平滑变化趋势，这让噪声可以更好地模拟自然世界的万物。下面讨论一维，二维，三维柏林噪声的原理及其实现。
1. 一维柏林噪声
首先，在X轴上每个整数座标随机生成一个数（范围为-1~1），我们称这个数为Gradient，译为梯度或者斜率。然后我们对相邻两个整数之间使用梯度进行插值计算，使得相邻两点之间平滑过渡。平滑度取决于所选用的插值函数，老版的柏林噪声使用f(t)=3*t*t-2*t*t*t，改进后的柏林噪声使用f(t)=t*t*t*(t*(t*6-15)+10)。如图所示：
上图来自于KenPerlin的Simplex Noise论文，论文中提到了经典的柏林噪声定义。
在这里，需要注意一个问题，上面说的梯度和我们在微积分里面学的梯度并不是一回事。在多元微分学中，梯度的定义是函数对各个座标轴的偏导数，只要空间中的位置确定，函数在这一点的梯度就已经确定。梯度的数学意义是空间沿着梯度的方向变化最快。这里的梯度是随机生成的，至于如何生成，这个关键问题在后面单独讨论。（起初把这里的梯度和微积分里的梯度等同，结果迟迟没有理解柏林噪声的原理）。
2. 二维柏林噪声
对于二维来说我们可以获取点P(x, y)最近的四个整数点ABCD，ABCD四个点的座标分别为A(i, j)、B(i+1, j)、C(i, j+1)、D(i+1, j+1)，随后获取ABCD四点的二维梯度值G(A)、G(B)、G(C)、G(D)，并且算出ABCD到P点的向量AP、BP、CP以及DP。如下图所示：
红色箭头表示该点处的梯度值，绿色箭头表示该点到P点的向量。接着，将G(A)与AP进行点乘，计算出A点对于P点的梯度贡献值，然后分别算出其余三个点对P点的梯度贡献值，最后将(u, v)代入插值函数中算出P点的最终噪声值。
3. 三维柏林噪声
将以上讨论类推到三维的柏林噪声，则需要算出八个顶点的梯度贡献值，然后进行插值计算。在后来的改进算法中，取的是立方体12条边的中点。这样处理来保证各个方向之间的不关联性，带来的是计算量的提高，增加了接近一半的计算量。
在一篇博客中，博主给出了详细的计算过程，想了解详细的过程可以看http://www.twinklingstar.cn/2015/2581/classical-perlin-noise/。
4. 随机梯度的选取
对于随机梯度的选取，柏林采用了这样一种方法：预先生成256个伪随机的梯度值（以及二维和三维的梯度向量）保存在G1[256]、G2[256][2]]以及G3[256][3]中，然后对于一维柏林噪声来说，我们可以直接去取G1[]数组中的梯度值使用，对于二维或者三维，预先又生成了一个排列P[256]，将0~255的下标随机存放在P数组中，然后通过下面公式来取到随机的梯度值：
G2[P[x] + y] ——二维
G3[P[P[x] + y] +z] ——三维
以上就完成了随机取梯度值的算法。
这里还要提一下在Improved Noise论文中柏林对于三维柏林噪声的改进。在论文中柏林的第一个改进就是插值函数的改进，也就是我们上面提到过的f(t)=t*t*t*(t*(t*6-15)+10)，它使得插值出来的噪声值更平滑了，特别是在三维中。
左边是以前的插值函数计算的结果，右边是改进后的插值函数计算的结果，不难看出，右边的效果提高了不少。
对于上面的过程，如果感觉不很懂，建议直接看柏林的论文。http://mrl.nyu.edu/~perlin/paper445.pdf
5. 实现
Perlin噪声实现比较简单的，在1983年的计算机上实现的算法也不允许对计算量、内存有多大的要求。概括来说，Perlin噪声的实现需要三个步骤：
# a. 定义一个晶格结构，每个晶格的顶点有一个“伪随机”的梯度向量（其实就是个向量啦）。对于二维的Perlin噪声来说，晶格结构就是一个平面网格，三维的就是一个立方体网格。
# b. 输入一个点（二维的话就是二维座标，三维就是三维座标，n维的就是n个座标），我们找到和它相邻的那些晶格顶点（二维下有4个，三维下有8个，n维下有2n X 2n个），计算该点到各个晶格顶点的距离向量，再分别与顶点上的梯度向量做点乘，得到2n X 2n个点乘结果。
# c. 使用缓和曲线（ease curves）来计算它们的权重和。在原始的Perlin噪声实现中，缓和曲线是s(t)=3t^2−2t^3 s(t)=3t^2−2t^3，在2002年的论文中[Perlin K. Improving noise[C]//ACM Transactions on Graphics (TOG). ACM, 2002, 21(3): 681-682](Perlin K. Improving noise[C]//ACM Transactions on Graphics (TOG). ACM, 2002, 21(3): 681-682)，Perlin改进为s(t)=6t^5−15t^4+10t^3 s(t)=6t^5−15t^4+10t^3。这里简单解释一下，为什么不直接使用s(t)=t，即线性插值。直接使用的线性插值的话，它的一阶导在晶格顶点处（即t = 0或t = 1）不为0，会造成明显的不连续性。s(t)=3t^2−2t^3 s(t)=3t^2−2t^3在一阶导满足连续性，s(t)=6t^5−15t^4+10t^3 s(t)=6t^5−15t^4+10t^3在二阶导上仍然满足连续性。
以二维为例，下面的图可以展示第一步和第二步的过程。
6. 算法时间复杂度分析
Perlin噪声算法的时间复杂度为O(2^n)（n是维数），2^n其实就是每个晶格的顶点个数。可见对于输入点需要进行指数级别的点乘操作。因此，随着维数的增大，例如计算三维四维的Perlin噪声，时间会迅速增加。为了解决这个问题，Perlin在2002年提出了优化后的一种噪声——Simplex噪声，这种噪声的复杂度是O(n^2)，在高纬度上时间优化明显,这样的时间复杂度类似于based-compared sort，时间复杂度显著下降。
对于算法在时间复杂度上的提高，我们可以把柏林噪声的实现算法当成斐波那契数列的递归实现，把Simplex算法当成一般的冒泡排序。这2个实验我们都做过，运行时间应该都有直观的印象。
五. 分形噪声
分形噪声可以用来模拟自然界的自相似过程，包括海岸线，地形，海浪等。分形噪声的原理是利用Perlin噪声频率受限的特性，通过不断叠加更高频率的perlin噪声达到自相似的效果。（定义来自维基百科https://zh.wikipedia.org/wiki/Perlin%E5%99%AA%E5%A3%B0）
关于分形噪声的详细介绍可以参看这篇博客。https://blog.csdn.net/syfolen/article/details/75051066
1. 一维分形噪声
在一维的情况下，设噪声函数为noise(x)，则通过 noise(2x), noise(4x)等就可以构造更高频率的噪声。下图第一张图是一维Perlin噪声频率频率在不断倍增，振幅不断倍减的情况，第二张图是第一张图中这些函数求和的图像。
不同频率下的分形噪声
合并后的分形噪声
上面第二张图像是自相似的，将它放大会得到和原图相似的效果。从数学上描述分形噪声，则是：
分形噪声有时也被称作1/f噪声
2. 二维分形噪声
对于二维的情况，可以使用类似的方法进行分形和运算，产生的效果类似于云层，如下图所示。
# 上图依次为：二维Simplex噪声原函数， 二维Simplex噪声叠加2倍频， 二维Simplex噪声叠加到4倍频， 二维Simplex噪声叠加到128倍频。
对分形和函数进行一些小修改可以模拟各种复杂效果，具体详见如下链接https://zh.wikipedia.org/wiki/Perlin%E5%99%AA%E5%A3%B0
六. Simplex噪声
在2001年SIGGRAPH Course上，Ken Perlin进行了一次演讲，他介绍了对Perlin噪声的一个改进版噪声——Simplex噪声。给出Perlin原版的Course Note链接https://www.csee.umbc.edu/~olano/s2002c36/ch02.pdf。，Simplex噪声的计算复杂度为O(n^2)，优于Perlin噪声的O(2^n)。而且在效果上，Simplex噪声也克服了经典的Perlin噪声在某些视觉问题。但Simplex噪声很难理解，我看了好几天也理解了一个大概，因此这一部分写的忐忐忑忑，错误在所难免。
Simplex 噪声使用 6t^5 – 15t^4 + 10t^3 作为插值函数，杜绝了导数中的线性部分。另外，Simplex 噪声算法有效地降低了插值次数。
1. 实现
二维经典Perlin噪声将二维空间用正方形填充，用四个顶点进行3次插值，而Simplex噪声将二维空间用等边三角形填充，使用三个顶点进行插值。三维经典Perlin噪声将三维空间用立方体进行填充，使用8个顶点进行7次插值，三维Simplex噪声使用正四面体填充空间，用4个顶点进行插值。对于更高维（N维）的情况，经典Perlin噪声将空间用超立方体填充，顶点数目是 2^N，而Simplex噪声使用高维正三角形（称为Simplex）进行填充，顶点数目为 N+1。Simplex噪声插值次数随维数增长线性增长。
对于二维的情况，使用正三角形填充空间使得直接判断某点落在哪个正三角形中，计算该三角形的顶点位置变得复杂。在实现上通常通过座标变换将正三角形映射成直角三角形。使用该方法进行变换可以使用和经典Perlin噪声相同的方法对顶点进行求值.http://staffwww.itn.liu.se/~stegu/simplexnoise/simplexnoise.pdf
网上一篇博客对Simplex噪声的实现进行了详细的介绍，下面对其进行总结。先给出博客地址，看不懂下面的总结可以直接去看博客，博客讲的比较详细。https://blog.csdn.net/candycat1992/article/details/50346469
Simplex噪声和Perlin噪声在实现上的区别是他的晶格不是方形，而是单形。单形就是N维空间中最简单紧凑的多边形。二维空间单形就是等边三角形，三维空间就是四面体。单形的顶点数目远少于N方体，因此计算梯度噪声时可以大大减少需要计算的顶点权重数目。
以二维空间为例：
每个顶点的噪声贡献度为：(r^2−|dist|^2)^4×dot(dist,grad)dist。dist 是输入点到该顶点的距离向量，gradgrad是该点存储的伪随机梯度向量，r2r2的取值是0.5或0.6。取0.5时可以保证没有不连续的间断点，在连续性并不那么明显时可以取0.6得到更好的视觉效果。
上面的实现中我们始终忽略了一个问题，就是如何找到输入点所在的单形？下面对这个问题进行探讨。
2.寻找单体
把单形进行座标偏斜（skewing），把平铺空间的单形变成一个新的网格结构，这个网格结构是由超立方体组成的，而每个超立方体又由一定数量的单形构成。
之前提到的单形网格如上图中的红色网格所示，它们有一些等边三角形组成（注意到这些等边三角形是沿空间对角线排列的）。经过座标倾斜后，它们变成了后面的黑色网格，这些网格由正方形组成，每个正方形是由之前两个等边三角形变形而来的三角形组成。这个把N维空间下的单形网格变形成新网格的公式如下：
x′=x+(x+y+…)⋅K1
y′=y+(x+y+…)⋅K1
其中，K1=(sqrt(n + 1) – 1) / n.
接下来就很难懂了，我可能也讲不清楚，如果对这方面感兴趣可以去参看上面给出的博客链接。
3. Simplex噪声和经典Perlin噪声的应用
Perlin噪声可以用来模拟自然界中的噪声现象。由于它的连续性，如果将二维噪声中的一个轴作为时间轴，得到的就是一个连续变化的一维函数。同样的也可以得到连续变化的二维图像。该噪声可以用来模拟人体的随机运动，蚂蚁行进的线路等。另外，还可以通过计算分形和模拟云朵，火焰等非常复杂的纹理。
Perlin噪声对各个点的计算是相互独立的，因此非常适合使用图形处理器进行计算。OpenGL在GLSL中定义了一维至四维的噪声函数noise1、noise2、noise3、noise4，在该规范中噪声的性质与上述Perlin提出的性质相同。在Mesa 3D实现中，这一组函数是使用Simplex算法实现的。在硬件的实现中，噪声的生成可以达到实时效果。
使用Simplex噪声模拟木纹
七. Unity中Perlin Noise的实现
首先给出出处：http://flafla2.github.io/2014/08/09/perlinnoise.html
下面介绍Simplex的算法实现.
public double perlin(double x, double y, double z);
函数接收x,y,z三个座标分量作为输入，并返回0.0~1.0的double值作为输出。那我们应该怎么处理输入值？首先，我们取3个输入值x,y,z的小数点部分，就可以表示为单元空间里的一个点了。为了方便讲解，我们将问题降维到2维空间来讨论（原理是一样的），下图是该点在2维空间上的表示：
小蓝点代表输入值在单元正方形里的空间座标，其他4个点则是单元正方形的各顶点
接着，我们给4个顶点（在3维空间则是8个顶点）各自生成一个伪随机的梯度向量。梯度向量代表该顶点相对单元正方形内某点的影响是正向还是反向的（向量指向方向为正向，相反方向为反向）。而伪随机是指，对于任意组相同的输入，必定得到相同的输出。因此，虽然每个顶点生成的梯度向量看似随机，实际上并不是。这保证了在生成函数不变的情况下，每个座标的梯度向量都是确定不变的。
各顶点上的梯度向量随机选取结果
上图所示的梯度向量并不是完全准确的。在本文所介绍的改进版柏林噪声中，这些梯度向量并不是完全随机的，而是由单位正方体（3维）的中心点指向各条边中点的12个向量：(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)，具体原因参见：http://mrl.nyu.edu/~perlin/paper445.pdf
接着，我们需要求出另外4个距离向量（在3维空间则是8个），它们分别从各顶点指向输入点（蓝色点）。下面有个2维空间下的例子：
各个距离向量
对每个顶点的梯度向量和距离向量做点积运算，就可以得出每个顶点的影响值：grad.x * dist.x + grad.y * dist.y + grad.z * dist.z
点积可以理解为向量a在向量b上的投影，当距离向量在梯度向量上的投影为同方向，点积结果为正数；当距离向量在梯度向量上的投影为反方向，点积结果为负数。因此，通过两向量点积，我们就知道该顶点的影响值是正还是负的。不难看出，顶点的梯度向量直接决定了这一点。下面通过一副彩色图，直观地看下各顶点的影响值：
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。
ease curve
ease curve的值会用来计算前面代码里的u和v，这样插值变化不再是单调的线性变化，而是这样一个过程：初始变化慢，中间变化快，结尾变化又慢下来（也就是在当数值趋近于整数时，变化变慢）。这个用于改善柏林噪声算法的fade函数可以表示为以下数学形式：
基本上，这就是整个柏林噪声算法的原理了！
代码实现
1. 准备工作
第一步，我们需要先声明一个排列表（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 overflow
static Perlin() {
p = new int[512];
for(int x=0;x<512;x++) {
p[x] = permutation[x%256];
}
}
p[]数组会在算法后续的哈希计算中使用到，用于确定一组输入最终挑选哪个梯度向量（从前面所列出的12个梯度向量中挑选）。后续代码会详细展示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[]时，出现数组越界错误。这也产生了一个副作用：柏林噪声每隔256个整数就会再次重复。但这不是太大的问题，因为算法不仅能处理整数，还能处理小数。最后，我们通过xf, yf, zf三个变量（也就是x,y,z的小数部分值），确定了输入座标在单元正方体里的空间位置（就是前面所示的小蓝点）。
2. Fade函数
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);
// ...
}
3. 哈希函数
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内。
4. 梯度函数
grad()函数的作用在于计算随机选取的梯度向量以及顶点位置向量的点积。Ken Perlin巧妙地使用了位翻转(bit-flipping)技巧来实现：
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.
}
下面代码则是以另一种令人容易理解的方式完成了这个任务（而且在很多语言版本的运行效率都优于前面一种）：
// Source: http://riven8192.blogspot.com/2010/08/calculate-perlinnoise-twice-as-fast.html
public 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 -y - z;
case 0xC: return y + x;
case 0xD: return -y + z;
case 0xE: return y - x;
case 0xF: return -y - z;
default: return 0; // never happens
}
}
上面的两种实现并没有实质差别。他们都是从以下12个向量里随机挑选一个作为梯度向量：(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).随机挑选结果其实取决于前一步所计算得出的哈希值（grad()函数的第一个参数）。后面3个参数则代表由输入点指向顶点的距离向量（最终拿来与梯度向量进行点积）。
5. 插值整合
经过前面的几步计算，我们得出了8个顶点的影响值，并将它们进行平滑插值，得出了最终结果：
public double perlin(double x, double y, double z) {
// ...
double x1, x2, y1, y2;
x1 = lerp( grad (aaa, xf , yf , zf), // The gradient function calculates the dot product between a pseudorandom
grad (baa, xf-1, yf , zf), // gradient vector and the vector from the input coordinate to the 8
u); // surrounding points in its unit cube.
x2 = lerp( grad (aba, xf , yf-1, zf), // This is all then lerped together as a sort of weighted average based on the faded (u,v,w)
grad (bba, xf-1, yf-1, zf), // values we made earlier.
u);
y1 = lerp(x1, x2, v);
x1 = lerp( grad (aab, xf , yf , zf-1),
grad (bab, xf-1, yf , zf-1),
u);
x2 = lerp( grad (abb, xf , yf-1, zf-1),
grad (bbb, xf-1, yf-1, zf-1),
u);
y2 = lerp (x1, x2, v);
return (lerp (y1, y2, w)+1)/2; // For convenience we bind the result to 0 - 1 (theoretical min/max before is [-1, 1])
}
// Linear Interpolate
public static double lerp(double a, double b, double x) {
return a + x * (b - a);
}
6. 利用倍频实现更自然的噪声
最后让我们再思考下，除了前面所讲的计算，还有其他办法可以令噪声结果更加自然吗？虽然柏林噪声算法一定程度上模拟了自然噪声，但仍没有完全表现出自然噪声的不规律性。举个现实例子，现实地形会有大段连绵、高耸的山地，也会有丘陵和蚀坑，更小点的有大块岩石，甚至更小的鹅卵石块，这都属于地形的一部分。那如何让柏林噪声算法模拟出这样的自然噪声特性，解决方法也很简单：我们可以使用不同的频率（frequencies）和振幅（amplitudes）参数进行多几次柏林噪声计算，然后将结果叠加在一起。频率是指采样数据的间隔，振幅是指返回值的幅度范围。
public double OctavePerlin(double x, double y, double z, int octaves, double persistence) {
double total = 0;
double frequency = 1;
double amplitude = 1;
double maxValue = 0; // Used for normalizing result to 0.0 - 1.0
for(int i=0;i<octaves;i++) {
total += perlin(x * frequency, y * frequency, z * frequency) * amplitude;
maxValue += amplitude;
amplitude *= persistence;
frequency *= 2;
}
return total/maxValue;
}
