SIGGRAPH中海洋的研究学习

演示demo：

https://v.qq.com/x/page/g08391vnr97.html

从海岛奇兵的海水一路改进过来，但总感觉还是不够好看。想来想去还是重新写一个新版海水。总体思路不再是优先考虑性能，而是先做效果，只要手机上还能支持，就先试试看。

打算先做Gerstner Wave。

数学部分知识如下：（来自https://zhuanlan.zhihu.com/p/31670275）

实际实现的时候还是挺麻烦的。首先要自己创建一个网格，因为要做效果，这个网格的顶点数要多一点，我用的是程序动态生成，可以调整精细度。生成网格代码就不再赘述（因为又臭又长）。

波形公式有了，但具体用几个波进行叠加，怎么叠加却没有明确的说法。我查了好多资料，参考了UWA上面的一个项目https://lab.uwa4d.com/lab/5b55ee58d7f10a201fd760a9，最终决定分层随机叠加的方式。通过参数进行调节。以下代码让波长逐渐增加，同时生成了相位和角度。角度是用来控制波的方向。

     public void GenerateWaveData(int componentsPerOctave, ref float[] wavelengths, ref float[] anglesDeg, ref float[] phases)
        {
            int totalComponents = NUM_OCTAVES * componentsPerOctave;
     
            if (wavelengths == null || wavelengths.Length != totalComponents) wavelengths = new float[totalComponents];
            if (anglesDeg == null || anglesDeg.Length != totalComponents) anglesDeg = new float[totalComponents];
            if (phases == null || phases.Length != totalComponents) phases = new float[totalComponents];
     
            float minWavelength = Mathf.Pow(2f, SMALLEST_WL_POW_2);
            float invComponentsPerOctave = 1f / componentsPerOctave;
     
            for (int octave = 0; octave < NUM_OCTAVES; octave++)
            {
                for (int i = 0; i < componentsPerOctave; i++)
                {
                    int index = octave * componentsPerOctave + i;
                    float minWavelengthi = minWavelength + invComponentsPerOctave * minWavelength * i;
                    float maxWavelengthi = Mathf.Min(minWavelengthi + invComponentsPerOctave * minWavelength, 2f * minWavelength);
                    wavelengths[index] = Mathf.Lerp(minWavelengthi, maxWavelengthi, Random.value);
     
                    float rnd;
     
                    rnd = (i + Random.value) * invComponentsPerOctave;
                    anglesDeg[index] = (2f * rnd - 1f) * _waveDirectionVariance;
     
                    rnd = (i + Random.value) * invComponentsPerOctave;
                    phases[index] = 2f * Mathf.PI * rnd;
                }
     
                minWavelength *= 2f;
            }
        }

光有相位和波长还不够，还需要振幅。根据一篇论文里的说法，海洋是可以根据相位和波长算出合理振幅的。论文地址如下：

https://hal.archives-ouvertes.fr/file/index/docid/307938/filename/frechot_realistic_simulation_of_ocean_surface_using_wave_spectra.pdf

我自己也没具体看，而是直接拿了结果：

     public float GetAmplitude(float wavelength, float componentsPerOctave)
     {
     
            float wl_pow2 = Mathf.Log(wavelength) / Mathf.Log(2f);
            wl_pow2 = Mathf.Clamp(wl_pow2, SMALLEST_WL_POW_2, SMALLEST_WL_POW_2 + NUM_OCTAVES - 1f);
     
            int index = (int)(wl_pow2 - SMALLEST_WL_POW_2);
            float wl_lo = Mathf.Pow(2f, Mathf.Floor(wl_pow2));
            float k_lo = 2f * Mathf.PI / wl_lo;
            float omega_lo = k_lo * ComputeWaveSpeed(wl_lo);
            float wl_hi = 2f * wl_lo;
            float k_hi = 2f * Mathf.PI / wl_hi;
            float omega_hi = k_hi * ComputeWaveSpeed(wl_hi);
     
            float domega = (omega_lo - omega_hi) / componentsPerOctave;
     
            float a_2 = 2f * Mathf.Pow(10f, _powerLog[index]) * domega;
            var a = Mathf.Sqrt(a_2);
            return a;
    }

对于Gerstner Wave的处理，uwa那个项目还有一种非常神奇的做法，一般都是在顶点着色器里对多个波形叠加，通过增加顶点数来提高精度，而它直接用cb先在片段着色器里画出波形并且存到贴图中，然后再对贴图进行采样得到位置。毋庸置疑这种做法得到的波是非常平滑自然的，特别美妙。具体步骤如下：

1.创建好海面网格。可以是普通平面或者是回字形平面。后者更适合优化和平一点的视角。

2.用程序计算生成Gerstner Wave的一系列参数，传递给材质进行渲染。

3.渲染流程是通过commandbuff去做的。Gerstner Wave前面说了，是通过片段着色器去渲染，所以直接画一个四边形就行。

shader直接从UWA项目中抄录如下：

    //四边形uv是0-1，调整到[-0.5,0.5]之间。texelSize是生成的贴图的大小，i_res是缩放过的系数。假设是放大8倍的四边形，那么i_res就是32/size,也就是32*[-0.5,0.5]， 就是[-16,16]，而回字形海面刚好是4x4的格子，对应正确。再从中心进行偏移，就成功从uv转到世界坐标了（这里其实是回字形特有的算法，不必深究，只要知道是从uv得到世界坐标就好）。
    float2 LD_UVToWorld(in float2 i_uv, in float2 i_centerPos, in float i_res, in float i_texelSize)
    {
        return i_texelSize * i_res * (i_uv - 0.5) + i_centerPos;
    }
     
     
    float2 LD_0_UVToWorld(in float2 i_uv)
    {
        return LD_UVToWorld(i_uv, _LD_Pos_Scale_0.xy, _LD_Params_0.y, _LD_Params_0.x);
    }
     
    v2f vert( appdata_t v )
    {
        v2f o;
        o.vertex = float4(v.vertex.x, -v.vertex.y, 0., .5);
        float2 worldXZ = LD_0_UVToWorld(v.uv);
        o.worldPos_wt.xy = worldXZ;
        o.uv = v.uv;
        return o;
    }
     
    //GridSize是每个像素代表的长度，由外部传入
    float MinWavelengthForCurrentOrthoCamera()
    {
        return _GridSize * _TexelsPerWave;
    }
     
    //波速可以通过公式获得，具体参考下面链接地址
     
    float ComputeWaveSpeed(float wavelength, float g)
    {
        // wave speed of deep sea ocean waves: https://en.wikipedia.org/wiki/Wind_wave
        // https://en.wikipedia.org/wiki/Dispersion_(water_waves)#Wave_propagation_and_dispersion
        //float g = 9.81; float k = 2. * 3.141593 / wavelength; float cp = sqrt(g / k); return cp;
        const float one_over_2pi = 0.15915494;
        return sqrt(wavelength*g*one_over_2pi);
    }
     
    half4 frag (v2f i) : SV_Target
    {
        const half minWavelength = MinWavelengthForCurrentOrthoCamera();
        half3 result = (half3)0.;
     
        for (uint vi = 0; vi < BATCH_SIZE / 4; vi++)
        {
            [unroll]
            for (uint ei = 0; ei < 4; ei++)
            {
                if (_Wavelengths[vi][ei] == 0.)
                {
                    return half4(result, 0.);
                }
     
     
                half wt = 1;
     
    //按照求解公式，我们可以找到对应项，D是方向，点乘P，也就是位置，C和NowTime对应tφ，k就是频率，2π/波长，最后就是振幅A和Q，Q这里是_Chop，从外部传入，到这里，公式已经计算完毕，就可以得到最终的波形图了，存在rgbafloat的贴图中
                half C = ComputeWaveSpeed(_Wavelengths[vi][ei], _Gravity * _GravityScales[vi][ei]);
                half2 D = half2(cos(_Angles[vi][ei]), sin(_Angles[vi][ei]));
     
                half k = TWOPI / _Wavelengths[vi][ei];
                half x = dot(D, i.worldPos_wt.xy);
     
                half3 result_i = wt * _Amplitudes[vi][ei];
                result_i.y *= cos(k*(x + C * NowTime) + _Phases[vi][ei]);
                result_i.xz *= -_Chop * _ChopScales[vi][ei] * D * sin(k*(x + C * NowTime) + _Phases[vi][ei]);
                result += result_i;
            }
            }
     
            return half4(i.worldPos_wt.z * result, 0.);
    }

拿到这张波形图之后，就可以对我们前面生成的网格进行扰动了。在这之前，UWA项目里面对回字形的两层LOD进行了混合叠加处理，用来使过度更加自然。代码如下：

    //这里采样两次，也是回字形造成的
    void SampleDisplacements(in sampler2D i_dispSampler, in float2 i_uv, in float i_wt, inout float3 io_worldPos)
    {
        const half3 disp = tex2Dlod(i_dispSampler, float4(i_uv, 0., 0.)).xyz;
        io_worldPos += i_wt * disp;
    }
     
    half4 frag (v2f i) : SV_Target
    {
        const float2 worldPosXZ = LD_0_UVToWorld(i.uv);
        // sample the shape 1 texture at this world pos
        const float2 uv_1 = LD_1_WorldToUV(worldPosXZ);
        float3 result = 0.;
            SampleDisplacements(_LD_Sampler_AnimatedWaves_0, i.uv, 1.0, result);
        // waves to combine down from the next lod up the chain
        SampleDisplacements(_LD_Sampler_AnimatedWaves_1, uv_1, 1.0, result);
        return half4(result, 1.);
    }

好了，现在终于可以进入顶点着色器看怎么进行波形扰动了。

     
        void OnWillRenderObject()
        {
            Camera.current.depthTextureMode |= DepthTextureMode.Depth;
     
            // per instance data
     
            if (_mpb == null)
            {
                _mpb = new MaterialPropertyBlock();
            }
            _rend.GetPropertyBlock(_mpb);
     
            float meshScaleLerp = 0f;
            float farNormalsWeight = 1f;
            _mpb.SetVector("_InstanceData", new Vector4(meshScaleLerp, farNormalsWeight, _lodIndex));
     
            //每个小格子的长度
            float squareSize = Mathf.Pow(2f, Mathf.Round(Mathf.Log(transform.lossyScale.x) / Mathf.Log(2f))) / _baseVertDensity;
            float mul = 1.875f; // fudge 1
            float pow = 1.4f; // fudge 2
            float normalScrollSpeed0 = Mathf.Pow(Mathf.Log(1f + 2f * squareSize) * mul, pow);
            float normalScrollSpeed1 = Mathf.Pow(Mathf.Log(1f + 4f * squareSize) * mul, pow);
            _mpb.SetVector("_GeomData", new Vector3(squareSize, normalScrollSpeed0, normalScrollSpeed1));
     
            // assign lod data to ocean shader
            var ldaws = Ocean.Instance._lodDataAnimWaves;
     
            ldaws.BindResultData(_lodIndex, 0, _mpb);
     
            if (_lodIndex + 1 < Ocean.Instance.CurrentLodCount)
            {
                ldaws.BindResultData(_lodIndex + 1, 1, _mpb);
            }
     
     
            _mpb.SetTexture(_reflectionTexId, Texture2D.blackTexture);
     
            _rend.SetPropertyBlock(_mpb);
        }
     
     
     
        public void SetInstanceData(int lodIndex, int totalLodCount, float baseVertDensity)
        {
            _lodIndex = lodIndex; _totalLodCount = totalLodCount; _baseVertDensity = baseVertDensity;
        }

着色器代码如下：

    //这里有一个小技巧，就是以最小格子为单位进行移动，因为在顶点数有限的情况下，如果顶点移动不是跳跃式，那么中间插值会导致轻微闪烁。这里采用的是2倍最小格子，原因是三角形分布式2x2对称的，保持稳定性。后面是因为回字形缩放，要让边缘部分逐渐放大到两倍，和下一个lod完美对齐
     
    float ComputeLodAlpha(float3 i_worldPos, float i_meshScaleAlpha)
    {
     
        float2 offsetFromCenter = float2(abs(i_worldPos.x - _OceanCenterPosWorld.x), abs(i_worldPos.z - _OceanCenterPosWorld.z));
        float taxicab_norm = max(offsetFromCenter.x, offsetFromCenter.y);
     
        float lodAlpha = taxicab_norm / _LD_Pos_Scale_0.z - 1.0;
            
        const float BLACK_POINT = 0.15, WHITE_POINT = 0.85;
        lodAlpha = max((lodAlpha - BLACK_POINT) / (WHITE_POINT - BLACK_POINT), 0.);
     
        lodAlpha = min(lodAlpha, 1.);
        return lodAlpha;
    }
     
    void SnapAndTransitionVertLayout(float i_meshScaleAlpha, inout float3 io_worldPos, out float o_lodAlpha)
    {
        const float SQUARE_SIZE_2 = 2.0*_GeomData.x, SQUARE_SIZE_4 = 4.0*_GeomData.x;
        io_worldPos.xz -= frac(unity_ObjectToWorld._m03_m23 / SQUARE_SIZE_2) * SQUARE_SIZE_2;
        o_lodAlpha = ComputeLodAlpha(io_worldPos, i_meshScaleAlpha);
        float2 m = frac(io_worldPos.xz / SQUARE_SIZE_4); // this always returns positive
        float2 offset = m - 0.5;               
        const float minRadius = 0.26;
        if (abs(offset.x) < minRadius) io_worldPos.x += offset.x * o_lodAlpha * SQUARE_SIZE_4;
        if (abs(offset.y) < minRadius) io_worldPos.z += offset.y * o_lodAlpha * SQUARE_SIZE_4;
    }
     
    //采样偏移之后，还需要计算法线，通过xz两个方向，分别进行采样，相减再叉乘，就会得到法线，可以画图求解
    void SampleDisplacementsNormals(in sampler2D i_dispSampler, in float2 i_uv, in float i_wt, in float i_invRes, in float i_texelSize, inout float3 io_worldPos, inout half2 io_nxz)
    {
        const float4 uv = float4(i_uv, 0., 0.);
        const half3 disp = tex2Dlod(i_dispSampler, uv).xyz;
        io_worldPos += i_wt * disp;
        float3 n;
            {
            float3 dd = float3(i_invRes, 0.0, i_texelSize);
            half3 disp_x = dd.zyy + tex2Dlod(i_dispSampler, uv + dd.xyyy).xyz;
            half3 disp_z = dd.yyz + tex2Dlod(i_dispSampler, uv + dd.yxyy).xyz;
            n = normalize(cross(disp_z - disp, disp_x - disp));
        }
        io_nxz += i_wt * n.xz;
    }
     
    v2f vert( appdata_t v )
    {
        v2f o;
        o.worldPos = mul(unity_ObjectToWorld, v.vertex);
        float lodAlpha;
        SnapAndTransitionVertLayout(_InstanceData.x, o.worldPos, lodAlpha);
        o.lodAlpha_worldXZUndisplaced_oceanDepth.x = lodAlpha;
        o.lodAlpha_worldXZUndisplaced_oceanDepth.yz = o.worldPos.xz;
     
        o.n_shadow = half4(0., 0., 0., 0.);
        o.foam_screenPos.x = 0.;
        o.lodAlpha_worldXZUndisplaced_oceanDepth.w = 0.;
     
        //根据权重，可以对两个lod分别采样混合
        float wt_0 = (1. - lodAlpha) * _LD_Params_0.z;
        float wt_1 = (1. - wt_0) * _LD_Params_1.z;
        // sample displacement textures, add results to current world pos / normal / foam
        const float2 worldXZBefore = o.worldPos.xz;
        if (wt_0 > 0.001)
        {
            const float2 uv_0 = LD_0_WorldToUV(worldXZBefore);
            SampleDisplacementsNormals(_LD_Sampler_AnimatedWaves_0, uv_0, wt_0, _LD_Params_0.w, _LD_Params_0.x, o.worldPos, o.n_shadow.xy);
            }
     
        if (wt_1 > 0.001)
        {
            const float2 uv_1 = LD_1_WorldToUV(worldXZBefore);
            SampleDisplacementsNormals(_LD_Sampler_AnimatedWaves_1, uv_1, wt_1, _LD_Params_1.w, _LD_Params_1.x, o.worldPos, o.n_shadow.xy);
     
        }
     
        // convert height above -1000m to depth below surface
        o.lodAlpha_worldXZUndisplaced_oceanDepth.w = DEPTH_BASELINE - o.lodAlpha_worldXZUndisplaced_oceanDepth.w;
     
        // foam can saturate
        o.foam_screenPos.x = saturate(o.foam_screenPos.x);
        // view-projection
        o.vertex = mul(UNITY_MATRIX_VP, float4(o.worldPos, 1.));
     
        UNITY_TRANSFER_FOG(o, o.vertex);            
        return o;
    }

做完以上步骤后，波形效果就出来了。

做完波形扰动后，就是要考虑开始着色，首先还是法线图一张，结合本身的法线进行基本的颜色显示。

    //法线贴图采样，uv通过两个魔数进行滚动。为了保证连续性，直接对下一个lod也进行一样的采样，但nstretch要翻倍，因为lod翻倍了，采样完毕后，把法线的值返回
    half2 SampleNormalMaps(float2 worldXZUndisplaced, float lodAlpha)
    {
        const float2 v0 = float2(0.94, 0.34), v1 = float2(-0.85, -0.53);
        const float geomSquareSize = _GeomData.x;
        float nstretch = _NormalsScale * geomSquareSize; // normals scaled with geometry
        const float spdmulL = _GeomData.y;
        half2 norm =
        UnpackNormal(tex2D(_Normals, (v0*NowTime*spdmulL + worldXZUndisplaced) / nstretch)).xy +
        UnpackNormal(tex2D(_Normals, (v1*NowTime*spdmulL + worldXZUndisplaced) / nstretch)).xy;
     
        // blend in next higher scale of normals to obtain continuity
        const float farNormalsWeight = _InstanceData.y;
        const half nblend = lodAlpha * farNormalsWeight;
        if (nblend > 0.001)
        {
            // next lod level
            nstretch *= 2.;
            const float spdmulH = _GeomData.z;
            norm = lerp(norm,
                UnpackNormal(tex2D(_Normals, (v0*NowTime*spdmulH + worldXZUndisplaced) / nstretch)).xy +
                UnpackNormal(tex2D(_Normals, (v1*NowTime*spdmulH + worldXZUndisplaced) / nstretch)).xy,
                nblend);
        }
     
        // approximate combine of normals. would be better if normals applied in local frame.
        return _NormalsStrength * norm;
    }
     
    //拿到法线后，和原始法线进行叠加混合
    float pixelZ = LinearEyeDepth(i.vertex.z);
    half3 screenPos = i.foam_screenPos.yzw;
    half2 uvDepth = screenPos.xy / screenPos.z;
    float sceneZ01 = tex2D(_CameraDepthTexture, uvDepth).x;
    float sceneZ = LinearEyeDepth(sceneZ01);
     
    float3 lightDir = WorldSpaceLightDir(i.worldPos);
    // Soft shadow, hard shadow
    fixed2 shadow = (fixed2)1.0;
     
    // Normal - geom + normal mapping
    half3 n_geom = normalize(half3(i.n_shadow.x, 1., i.n_shadow.y));
    if (underwater) n_geom = -n_geom;
    half3 n_pixel = n_geom;
     
     
    n_pixel.xz += (underwater ? -1. : 1.) * SampleNormalMaps(i.lodAlpha_worldXZUndisplaced_oceanDepth.yz, i.lodAlpha_worldXZUndisplaced_oceanDepth.x);
    n_pixel = normalize(n_pixel);
     
     
     
    half3 OceanEmission(in const half3 i_view, in const half3 i_n_pixel, in const float3 i_lightDir,in const half4 i_grabPos, in const float i_pixelZ, in const half2 i_uvDepth, in const float i_sceneZ, in const float i_sceneZ01,in const half3 i_bubbleCol, in sampler2D i_normals, in sampler2D i_cameraDepths, in const bool i_underwater, in const half3 i_scatterCol)
    {
        half3 col = i_scatterCol;
     
        // underwater bubbles reflect in light
        col += i_bubbleCol;
     
        return col;
    }
     
     
    先把天空盒的光线算进去，根据视线和法线，可以算出反射光线，采样天空盒，在用菲尼尔处理一下，菲尼尔用的是schlick 近似公式https://en.wikipedia.org/wiki/Schlick%27s_approximation
     
    void ApplyReflectionSky(half3 view, half3 n_pixel, half3 lightDir, half shadow, half4 i_screenPos, inout half3 col)
    {
        // Reflection
        half3 refl = reflect(-view, n_pixel);
        half3 skyColour;
        skyColour = texCUBE(_Skybox, refl).rgb;
     
        // Fresnel
        const float IOR_AIR = 1.0;
        const float IOR_WATER = 1.33;
        // reflectance at facing angle
        float R_0 = (IOR_AIR - IOR_WATER) / (IOR_AIR + IOR_WATER); R_0 *= R_0;
        // schlick's approximation
        float R_theta = R_0 + (1.0 - R_0) * pow(1.0 - max(dot(n_pixel, view), 0.), _FresnelPower);
        col = lerp(col, skyColour, R_theta);
    }

完成以上着色部分后，海水看上去是这样：

接下来，我们要看准方向加一个平行光，让海面亮起来

    高光用的是传统Phong模型就可以达到效果。
    skyColour += pow(max(0., dot(refl, lightDir)), _DirectionalLightFallOff) * _DirectionalLightBoost * _LightColor0 * shadow;

然而还是非常丑，主要还是因为光照太过简单，而海必须要考虑的就是散射。本来SSS也是一个大命题，可以看好几本书，所幸的是对于海来说，用近似次表面散射也可以得到好的效果，基础原理就是越看向太阳，就越亮

    #if _SUBSURFACESCATTERING_ON
        {
     
    // light
    // use the constant term (0th order) of SH stuff - this is the average. it seems to give the right kind of colour
            col *= half3(unity_SHAr.w, unity_SHAg.w, unity_SHAb.w);
     
    // Approximate subsurface scattering - add light when surface faces viewer. Use geometry normal - don't need high freqs.
            half towardsSun = pow(max(0., dot(i_lightDir, -i_view)), _SubSurfaceSunFallOff);
            col += (_SubSurfaceBase + _SubSurfaceSun * towardsSun) * _SubSurfaceColour.rgb * _LightColor0 * shadow;
        }
    #endif // _SUBSURFACESCATTERING_ON

一下子就好看不少，更加通透，有一种散射的小感觉了。不过这样还远远不够，我们放入一个地形，就会发现地形和海衔接的部分还完全没有考虑。

波形需要处理，随着地形的阻挡，波应该要逐步减弱，这个可以通过深度计算去处理。

在合适的位置放一个摄像机，往地形拍摄，将深度写入图中。

    v2f vert( appdata_t v )
    {
        v2f o;
        o.vertex = UnityObjectToClipPos( v.vertex );
     
        float altitude = mul(unity_ObjectToWorld, v.vertex).y;
        o.depth = altitude - (_OceanCenterPosWorld.y - depthMax);
     
        return o;
    }
     
    float frag (v2f i) : SV_Target
    {
        return i.depth;
    }

生成这张图后，我们要回到波形生成的地方，根据深度值，重新调整波的振幅。

    拿到深度，并且把depth还原成离水面的距离，如果depth很小，那么波长就要变小
    const half depth = depthMax - tex2D(_LD_Sampler_SeaFloorDepth_0, i.uv).x;
    half wt = 1;
    half depth_wt = saturate(depth / (0.5 * _Wavelengths[vi][ei]));
    wt *= .1 + .9 * depth_wt;

这样子处理之后，岸边的波浪就小下去了。

边缘硬切很难看，首先要处理透明问题，透明的基本原则是深度越浅越透明，深度越深越不透明。

    const half2 uvBackground = i_grabPos.xy / i_grabPos.w;
    //根据法线方向折射处理
    half2 uvBackgroundRefract = uvBackground + _RefractionStrength * i_n_pixel.xz;
    half3 sceneColour;
    half3 alpha = 0.;
    float depthFogDistance;
     
     
    //从深度贴图获得深度，并和顶点的深度作比较，如果顶点深度大于背景深度，那么就把距离算出来
    否则说明水在物体下面
    const half2 uvDepthRefract = i_uvDepth + _RefractionStrength * i_n_pixel.xz;
    const float sceneZRefract = LinearEyeDepth(tex2D(i_cameraDepths, uvDepthRefract).x);
     
    // Compute depth fog alpha based on refracted position if it landed on an underwater surface, or on unrefracted depth otherwise
    if (sceneZRefract > i_pixelZ)
    {
        depthFogDistance = sceneZRefract - i_pixelZ;
    }
    else
    {
        depthFogDistance = i_sceneZ - i_pixelZ;
        uvBackgroundRefract = uvBackground;
    }
     
    sceneColour = tex2D(_BackgroundTexture, uvBackgroundRefract).rgb;
     
    //对透明度根据距离进行处理
    alpha = 1. - exp(-_DepthFogDensity.xyz * depthFogDistance);
     
    // blend from water colour to the scene colour
    col = lerp(sceneColour, col, alpha);

透明和扰动都有了，但是边缘部分的切边还是很明显。一般这种时候就需要泡沫来帮忙了。 以前我曾经用两层泡沫图叠加的方式去做，效果一般般。而且根据深度去产生泡沫也并不正确，在海浪的波峰也是有可能产生泡沫的，泡沫产生的原因主要是因为运动撕裂程度大。在海洋统计学里可以用雅克比行列式（完全看不懂）去做，这里也可以模仿。

    //这是大猫知乎上关于雅克比行列式求解过程
    for (int i = 0; i < resolution; i++)
    {
        for (int j = 0; j < resolution; j++)
        {
            int index = i * resolution + j;
            Vector2 dDdx = Vector2.zero;
            Vector2 dDdy = Vector2.zero;
                    //ddx就是将改点的偏移减去x轴一个像素的偏移,ddy对应y轴
                    if (i != resolution - 1)
                {
                dDdx = 0.5f * (hds[index] - hds[index + resolution]);
            }
            if (j != resolution - 1)
            {
                dDdy = 0.5f * (hds[index] - hds[index + 1]);
            }
                    //这是行列式的值，后面应该是调整的值
            float jacobian = (1 + dDdx.x) * (1 + dDdy.y) - dDdx.y * dDdy.x;
            Vector2 noise = new Vector2(Mathf.Abs(normals[index].x), Mathf.Abs(normals[index].z)) * 0.3f;
            float turb = Mathf.Max(1f - jacobian + noise.magnitude, 0f);
            float xx = 1f + 3f * Mathf.SmoothStep(1.2f, 1.8f, turb);
            xx = Mathf.Min(turb, 1.0f);
            xx = Mathf.SmoothStep(0f, 1f, turb);
            colors[index] = new Color(xx, xx, xx, xx);
        }
    }

    half frag(v2f i) : SV_Target
    {
        float4 uv = float4(i.uv_uv_lastframe.xy, 0., 0.);
        float4 uv_lastframe = float4(i.uv_uv_lastframe.zw, 0., 0.);
        // #if _FLOW_ON
        half4 velocity = half4(tex2Dlod(_LD_Sampler_Flow_1, uv).xy, 0., 0.);
        half foam = tex2Dlod(_LD_Sampler_Foam_0, uv_lastframe
                    - ((_SimDeltaTime * _LD_Params_0.w) * velocity)
                    ).x;
     
        half2 r = abs(uv_lastframe.xy - 0.5);
        if (max(r.x, r.y) > 0.5 - _LD_Params_0.w)
        {
            // no border wrap mode for RTs in unity it seems, so make any off-texture reads 0 manually
            foam = 0.;
        }
     
        // fade
        foam *= max(0.0, 1.0 - _FoamFadeRate * _SimDeltaTime);
     
        // sample displacement texture and generate foam from it
        const float3 dd = float3(_LD_Params_1.w, 0.0, _LD_Params_1.x);
        half3 s = tex2Dlod(_LD_Sampler_AnimatedWaves_1, uv).xyz;
        half3 sx = tex2Dlod(_LD_Sampler_AnimatedWaves_1, uv + dd.xyyy).xyz;
        half3 sz = tex2Dlod(_LD_Sampler_AnimatedWaves_1, uv + dd.yxyy).xyz;
        float3 disp = s.xyz;
        float3 disp_x = dd.zyy + sx.xyz;
        float3 disp_z = dd.yyz + sz.xyz;
        // The determinant of the displacement Jacobian is a good measure for turbulence:
        // > 1: Stretch
        // < 1: Squash
        // < 0: Overlap
            //把两边偏移相减，这里是直接算行列式，没有+1的操作，算出foam后，还要根据深度去加强foam
        float4 du = float4(disp_x.xz, disp_z.xz) - disp.xzxz;
        float det = (du.x * du.w - du.y * du.z) / (_LD_Params_1.x * _LD_Params_1.x);
        foam += 5. * _SimDeltaTime * _WaveFoamStrength * saturate(_WaveFoamCoverage - det);
     
        // add foam in shallow water. use the displaced position to ensure we add foam where world objects are.
        float4 uv_1_displaced = float4(LD_1_WorldToUV(i.worldXZ + disp.xz), 0., 1.);
        float signedOceanDepth = depthMax - tex2Dlod(_LD_Sampler_SeaFloorDepth_1, uv_1_displaced).x + disp.y;
        foam += _ShorelineFoamStrength * _SimDeltaTime * saturate(1. - signedOceanDepth / _ShorelineFoamMaxDepth);
     
        return foam;
    }

拿到生成的foam贴图后就可以开始渲染。按照ppt里的说法，泡沫分成两层，顶部是白色泡沫，下面是褪色的海浪。具体的数学公式我没有查到，非常遗憾，只能有一个大概解释。

    void SampleFoam(in sampler2D i_oceanFoamSampler, float2 i_uv, in float i_wt, inout half io_foam)
    {
        io_foam += i_wt * tex2Dlod(i_oceanFoamSampler, float4(i_uv, 0., 0.)).x;
    }
     
    half WhiteFoamTexture(half i_foam, float2 i_worldXZUndisplaced)
    {
    //这里负责白色泡沫
        half ft = lerp(
            tex2D(_FoamTexture, (1.25*i_worldXZUndisplaced + NowTime / 10.) / _FoamScale).r,
            tex2D(_FoamTexture, (3.00*i_worldXZUndisplaced - NowTime / 10.) / _FoamScale).r,
            0.5);
     
        // black point fade
        i_foam = saturate(1. - i_foam);
        return smoothstep(i_foam, i_foam + _WaveFoamFeather, ft);
    }
     
    void ComputeFoam(half i_foam, float2 i_worldXZUndisplaced, float2 i_worldXZ, half3 i_n, float i_pixelZ, float i_sceneZ, half3 i_view, float3 i_lightDir, half i_shadow, out half3 o_bubbleCol, out half4 o_whiteFoamCol)
    {
        half foamAmount = i_foam;
     
        //海岸线衰减
        foamAmount *= saturate((i_sceneZ - i_pixelZ) / _ShorelineFoamMinDepth);
     
        // Additive underwater foam - use same foam texture but add mip bias to blur for free
            //这里进行了偏移，类似于模糊处理
        float2 foamUVBubbles = (lerp(i_worldXZUndisplaced, i_worldXZ, 0.05) + 0.5 * NowTime * _WindDirXZ) / _FoamScale + 0.125 * i_n.xz;
        half bubbleFoamTexValue = tex2Dlod(_FoamTexture, float4(.74 * foamUVBubbles - _FoamBubbleParallax * i_view.xz / i_view.y, 0., 5.)).r;
        o_bubbleCol = (half3)bubbleFoamTexValue * _FoamBubbleColor.rgb * saturate(i_foam * _WaveFoamBubblesCoverage) * AmbientLight();
     
        // White foam on top, with black-point fading
        half whiteFoam = WhiteFoamTexture(foamAmount, i_worldXZUndisplaced);
     
     
        o_whiteFoamCol.rgb = _FoamWhiteColor.rgb * (AmbientLight() +    _WaveFoamLightScale * _LightColor0 * i_shadow);
     
     
        o_whiteFoamCol.a = _FoamWhiteColor.a * whiteFoam;
    }

这样处理后，浪花效果还不错。

但是边缘效果依然丑陋，主要是海岸线和海洋中心差别还是挺大的。想了下还是希望走类似于下图这样的波浪。

这里我做了简化处理，没有做多层，只做了简单的一层。

然后开始继续优化效果，首先是散射的问题，在视角增高后，海水看起来很暗，其实是因为没有正确散射引起的，视角增高的时候，会增强散射效果。

col += pow(saturate(0.5 + 2.0 * waveHeight / _SubSurfaceHeightMax), _SubSurfaceHeightPower) * _SubSurfaceCrestColour.rgb;

这还不够，浅水的地方海的散射会更强。

    void SampleSeaFloorHeightAboveBaseline(in sampler2D i_oceanDepthSampler, float2 i_uv, in float i_wt, inout half io_oceanDepth)
    {
        io_oceanDepth += i_wt * (tex2Dlod(i_oceanDepthSampler, float4(i_uv, 0., 0.)).x);
    }
     
    #if _SUBSURFACESHALLOWCOLOUR_ON
        float shallowness = pow(1. - saturate(depth / _SubSurfaceDepthMax), _SubSurfaceDepthPower);
        half3 shallowCol = _SubSurfaceShallowCol;
        col = lerp(col, shallowCol, shallowness);
    #endif

最后是焦散

    void ApplyCaustics(in const half3 i_view, in const half3 i_lightDir, in const float i_sceneZ, in sampler2D i_normals, inout half3 io_sceneColour)
    {
        // could sample from the screen space shadow texture to attenuate this..
        // underwater caustics - dedicated to P
        float3 camForward = mul((float3x3)unity_CameraToWorld, float3(0., 0., 1.));
        float3 scenePos = _WorldSpaceCameraPos - i_view * i_sceneZ / dot(camForward, -i_view);
        const float2 scenePosUV = LD_1_WorldToUV(scenePos.xz);
        half3 disp = 0.;
        // this gives height at displaced position, not exactly at query position.. but it helps. i cant pass this from vert shader
        // because i dont know it at scene pos.
                                        SampleDisplacements(_LD_Sampler_AnimatedWaves_1, scenePosUV, 1.0, disp);
        half waterHeight = _OceanCenterPosWorld.y + disp.y;
        half sceneDepth = waterHeight - scenePos.y;
        half bias = abs(sceneDepth - _CausticsFocalDepth) / _CausticsDepthOfField;
        // project along light dir, but multiply by a fudge factor reduce the angle bit - compensates for fact that in real life
        // caustics come from many directions and don't exhibit such a strong directonality
        float2 surfacePosXZ = scenePos.xz + i_lightDir.xz * sceneDepth / (4.*i_lightDir.y);
        half2 causticN = _CausticsDistortionStrength * UnpackNormal(tex2D(i_normals, surfacePosXZ / _CausticsDistortionScale)).xy;
        half4 cuv1 = half4((surfacePosXZ / _CausticsTextureScale + 1.3 *causticN + half2(0.044*NowTime + 17.16, -0.169*NowTime)), 0., bias);
        half4 cuv2 = half4((1.37*surfacePosXZ / _CausticsTextureScale + 1.77*causticN + half2(0.248*NowTime, 0.117*NowTime)), 0., bias);
     
        half causticsStrength = _CausticsStrength;
     
        io_sceneColour *= 1. + causticsStrength *
            (0.5*tex2Dbias(_CausticsTexture, cuv1).x + 0.5*tex2Dbias(_CausticsTexture, cuv2).x - _CausticsTextureAverage);
    }

全部效果叠加有点闪烁，自己简化了代码，没有严格按照文档的做法，所以我自己修改了边界条件，修复了这个问题。其次，由于没有缩放考虑，摄像机拉高的时候海面有很多噪点，我通过线性减少扰动和浅滩散射来处理。暂时就处理到这里。



总结一下，完整的演讲中的海水远远比我这个复杂，而且即便是实现其中的这么一小部分，我也有大量的细节没有理解清楚，或者没有找到对应的公式。再自己复原效果的过程中，大量简化了一些实现，勉强达到了可以看的效果，不过由于为了让每个参数明显，海面看上去稍显油腻或者说卡通了一点。在手机上跑几乎是不可能了，也难以简化到那个程度。等我再补补数学，再来继续搞这个海水吧。

最后插件地址是：

https://www.assetstore.unity3d.com/cn/?stay#!/content/140114
---------------------

牛逼，虽然看不懂

66+

感谢楼主分享！