Waves Moving Vertices----gerstner wave

https://catlikecoding.com/unity/tutorials/flow/waves/

animate vertices
create gerstner waves
control wave direction
combine multiple waves.

1 sine waves
animating textures can create the illusion of a moving surface, but the mesh surface itself remains motionless. this is fine for small ripples, but can not represent larger waves. on large bodies of water–like an ocean of big lake—the wind can create big waves that can persist for a long time. to represent these wind waves, we will make new shader that displaces mesh vertices vertically, using a sine wave function.

1.1 adjusting vertices
create a new surface shader named Waves.
we will leave the fragment surface function unchanged. instead, add another function vert to adjust the vertex data. this function has a single vertex parameter, both for input and ouput. we will use unity’s default vertex data structure, appdata_full.

Shader "Custom/Waves" {
	Properties {
		_Color ("Color", Color) = (1,1,1,1)
		_MainTex ("Albedo (RGB)", 2D) = "white" {}
		_Glossiness ("Smoothness", Range(0,1)) = 0.5
		_Metallic ("Metallic", Range(0,1)) = 0.0
	}
	SubShader {
		Tags { "RenderType"="Opaque" }
		LOD 200

		CGPROGRAM
		#pragma surface surf Standard fullforwardshadows
		#pragma target 3.0

		sampler2D _MainTex;

		struct Input {
			float2 uv_MainTex;
		};

		half _Glossiness;
		half _Metallic;
		fixed4 _Color;

		void vert(inout appdata_full vertexData) {}

		void surf (Input IN, inout SurfaceOutputStandard o) {
			fixed4 c = tex2D (_MainTex, IN.uv_MainTex) * _Color;
			o.Albedo = c.rgb;
			o.Metallic = _Metallic;
			o.Smoothness = _Glossiness;
			o.Alpha = c.a;
		}
		ENDCG
	}
	FallBack "Diffuse"
}

to indicate that the surface shader should use the vertex function, add vertex:vert to the surface pragma directive.

#pragma surface surf Standard fullforwardshadows vertex:vert

1.2 adjusting Y
ignore the z dimension for now, the position of each vertex can be defined as P=(x,y), where P is its final position, x is the original X coordinate, and y is the original y coordinate, both in object space.
to create a wave, we have to adjust the y component of P. the simplest way to make a wave is to use a sine wave based on x, so y=sinx. the final point is then p = (x,sinx).

void vert(inout appdata_full vertexData) {
			float3 p = vertexData.vertex.xyz;

			p.y = sin(p.x);

			vertexData.vertex.xyz = p;
		}

the result is a sine wave along the x dimension, which is constant along the z dimension. the quads of the plane are of unit size, so the entire plane covers a 10x10 are centered on its local origin. so we end up seeing 10/2pi=1.59 periods of a sine wave.

1.3 amplitude
the default amplitude of a sine wave is 1, but we do not need to limit ourselves to that.
let us add a property to our shader so we can use Py=asinx instead, where a is the amplitude.

Properties {
		…
		_Amplitude ("Amplitude", Float) = 1
	}
	SubShader {
		…

		half _Glossiness;
		half _Metallic;
		fixed4 _Color;
		float _Amplitude;

		void vert(inout appdata_full vertexData) {
			float3 p = vertexData.vertex.xyz;

			p.y = _Amplitude * sin(p.x);

			vertexData.vertex.xyz = p;
		}

		…
	}

1.4 wavelength
in the case of sinx, the length of a full sine wave is 2pi约等于6.28. this is the wavelength and let us make it configurable too.
to easily control the wavelength, we first have to multiply x by 2pi then divide by the desired wavelength.
so we end up with sin(2pix/λ), where λ (lamda) is the wavelength.

2pi divided by λ is known as the wave number k = 2pi/λ.
we could use this as the shader property, so we do not need to perform a division in the shader.
that is a useful optimization, but in this tutorial we will stick with the more user-friendly wavelength.

inside the shader, we will explicitly use the wave number, so we end up with Py = asin(kx).
这a是振幅
k是波的个数，其实就w角频率
k的计算公式： k = 2pi/λ
这里λ就是波长

Shader "Custom/Waves" 
{
	Properties {
		…
		_Wavelength ("Wavelength", Float) = 10
	}
	SubShader {
		…
		float _Amplitude, _Wavelength;

		void vert(inout appdata_full vertexData) {
			float3 p = vertexData.vertex.xyz;

			float k = 2 * UNITY_PI / _Wavelength;
			p.y = _Amplitude * sin(k * p.x);

			vertexData.vertex.xyz = p;
		}

		…
	}

1.5 speed
the wave needs to move, so we have to define a speed.
it is most convenient to use the phase speed c, which defines how fast the entire wave moves in units per second.
this is done by using the time offset kct.
to make the wave more in the positive direction, we have to subtract this from kx, so we end up with Py = sin(kx-kct) = sin(k(x-ct)).

Properties {
		_Speed ("Speed", Float) = 1
	}
	SubShader {
		…
		float _Amplitude, _Wavelength, _Speed;

		void vert(inout appdata_full vertexData) {
			float3 p = vertexData.vertex.xyz;

			float k = 2 * UNITY_PI / _Wavelength;
			p.y = _Amplitude * sin(k * (p.x - _Speed * _Time.y));

			vertexData.vertex.xyz = p;
		}

		…
	}

1.6 normal vectors
our surface is curved and moving, but the lighting is still that of a motionless flat plane.
that is because we have not changed the vertex normals yet.
instead of directly calculating the normal vector, let’s first look at the surface tangent vetor in the X dimension, T.
for a flat surface T = (x’,0) = (1,0), which corresponds to the original plane’s tangent.
but for our wave we have to use T = P’ = (x’, asin(k(x-ct))’)
we end up with T = (1, kacosf)

this makes sense, because chaning the wavelength also changes the slope of the wave.

to get the final tangent vector in the shader, we have to normalize T.

float k = 2 * UNITY_PI / _Wavelength;
float f = k * (p.x - _Speed * _Time.y);
p.y = _Amplitude * sin(f);
float3 tangent = normalize(float3(1, k * _Amplitude * cos(f), 0));

the normal vector is the cross product of both tangent vectors.
as our wave is constant in the Z dimension, the binormal is always the unit vector and can be ignored, so we end up with
N = (-kacosf, 1).
we can just grab the normalized tangent components after normalizing them.

float3 tangent = normalize(float3(1, k * _Amplitude * cos(f), 0));
float3 normal = float3(-tangent.y, tangent.x, 0);
vertexData.vertex.xyz = p;
vertexData.normal = normal;

2 gerstner waves
sine waves are simple, but they do not match the shape of real water waves. big wind waves are realistically modeled by the Stokes wave function, but it is rather complex.
instead, gernster waves are often used for realtime animation of water surfaces.

Gerstner waves are named after František Josef Gerstner, who discovered them. They’re also known as trochoidal waves, named after their shape, or periodic surface grativity waves, which describes their physical nature.

2.1 moving back and forth
the fundamental observation is that while waves move across a water surface, the water itself does not move along with it.
in the case of a sine wave, each surface point goes up and down, but does not move horizontally.

but actual water is not just the surface. there is more water underneath.
when the surface water moves down, where does the water below it go?
when the surface moves up, what fills the space below it?
it turns out that the surface points not only move up and down, they more forward and backward too.
half the time move along the with the wave, but the other half they move in the opposite direction. the same is true for the water below the surface, but the deeper u go the less movement there is.

specifically, each surface points moves in a circle, orbiting a fixed anchor point. as the crest of a wave approaches, the point moves toward it. after the crest passes, it slides back, and then the next crest comes along.
the results is that water bunches up in crests and and spreads out in troughs 低谷, and the same will happen to our vertices.

in reality surface points do drift and do not describe perfect circles, but Gerstner waves do not model this.
that is fine, because we will use the original vertex positions as the anchor points.

we can turn our sine wave into a circle by using p = (acosf, asinf), but that would collapse the entire plane into a single circle. intead, we have to anchor each point on its original X coordinate, so we need p = (x+acosf, asinf).

float k = 2 * UNITY_PI / _Wavelength;
float f = k * (p.x - _Speed * _Time.y);
p.x += _Amplitude * cos(f);
p.y = _Amplitude * sin(f);

2.2 normal vectors
because we changed the surface function, its derivate has also changed.
the x component of T used to be x’=1, but now it is a bit more complicated.
the derivative of the cosine is the negative sine, so we end up with T = (1-kasinf, kacosf).

float3 tangent = normalize(float3(1 - k * _Amplitude * sin(f),k * _Amplitude * cos(f),0));
float3 normal = float3(-tangent.y, tangent.x, 0);

2.3 preventing loops

while the resulting waves might look fine, this is not always the case.
for example, reducing the wavelength to 20 while keeping the amplitude at 10 produces weird results.

because the amlitude is so large relative to the wavelength, the orbits of the surface points overshoot and form loops above the surface.
if this was real water, then the waves would break and fall apart, but we can not represent that with Gerstner waves.

there is a relation between the wavelength and the wave amplitude.
we can use a = e^kb / k,
where b has to do with surface pressure.
the stronger the pressure, the flatter the waves. in case of zero pressure, we end up with a = 1/k, which produces 0度 crests,
the sharpest possible before looping.
we can just use a = s/k instead, where s is a measure of steepness, between 0 and 1, which is easier to work with.
then we have p = (x+s/kcosf, s/ksinf), which simplies our tangent of T = (1-ssinf, scosf).

Shader "Custom/Waves" {
	Properties {
		…
		//_Amplitude ("Amplitude", Float) = 1
		_Steepness ("Steepness", Range(0, 1)) = 0.5
		_Wavelength ("Wavelength", Float) = 10
		_Speed ("Speed", Float) = 1
	}
	SubShader {
		…
		float _Steepness, _Wavelength, _Speed;

		void vert(inout appdata_full vertexData) {
			float3 p = vertexData.vertex.xyz;

			float k = 2 * UNITY_PI / _Wavelength;
			float f = k * (p.x - _Speed * _Time.y);
			float a = _Steepness / k;
			p.x += a * cos(f);
			p.y = a * sin(f);

			float3 tangent = normalize(float3(
				1 - _Steepness * sin(f),
				_Steepness * cos(f),
				0
			));
			…
		}

		…
	}

2.4 phase speed
in reality, waves do not have an arbitrary phase speed.
it is related to the wave number, c = sqrt(g/k) = sqrt(gλ/2pi),
where g is the pull of gravity, roughly 9.8 on earth.

float k = 2 * UNITY_PI / _Wavelength;
float c = sqrt(9.8 / k);
float f = k * (p.x - c * _Time.y);

now we can eliminate the speed property.

3 wave direction
up to this point our waves only move in the X dimension. we are now going to remove this restriction.
this makes our calculations a bit more complex, as both x and z are needed to construct the final wave and its tangent vectors.

3.1 direction vector
to indicate the travel direction of our wave we will introduce the direction vector D=(dx, dz). this is purely an indication of direction, so it is a vector of unit length, ||D||=1.

now how much x contributes to the wave function is modulated by the x component of D.
so we get f=k(Dx.x-ct).
but z now also plays a role, in the same way, which leads to
f=k(Dx.x+Dz.z-ct).
in other words, we are using the dot product of D and the original X and Z coordinates.
so we end up with f = k(D.(x,z) -ct).

add a direction property to our shader and incoporate it into our function.
it should be a unit-length vector, but to make it easier to work with we will normalize it in the shader.
note that all vector properties are 4D, so just ignore the Z and W components.

Shader "Custom/Waves" 
{
	Properties 
	{
		…
		_Wavelength ("Wavelength", Float) = 10
		_Direction ("Direction (2D)", Vector) = (1,0,0,0)
	}
	SubShader {
		…
		float _Steepness, _Wavelength;
		float2 _Direction;

		void vert(inout appdata_full vertexData) 
		{
			float3 p = vertexData.vertex.xyz;

			float k = 2 * UNITY_PI / _Wavelength;
			float c = sqrt(9.8 / k);
			float2 d = normalize(_Direction);
			float f = k * (dot(d, p.xz) - c * _Time.y);
			…
		}
		…
	}

we also have to adjust the horizontal offsets of Px and Pz so they align with the wave direction.
so instead of just adding the offset to x, we have to add it to z as well, in both cases modulated by the appropriate
component of D. so the final calculation becomes

p.x += d.x * (a * cos(f));
p.y = a * sin(f);
p.z += d.y * (a * cos(f));

3.2 normal vectors
once again, we have to adjust the calculation of our tangent, but not just for the X dimension.
we now also have to calculate the tangent in the Z dimension, the binormal vector B.

4.3 looping animations
now that we have two waves, u can observe that one with a longer wavelength indeed travels faster than a shorter one.
but the relationship between phase speed and wavelength is nonlinear, because

this is releant when u want to create a looping animation with multiple waves.
in the case of two waves, u have to find two wavelengths that produce phase speeds with the relationship ac1 = bc2, where a and b are integers.
u could do this by using even powers of two for wavelengths.