ゲルストナー波の法線を計算する

前回に引き続き、ゲルストナー波をやります。

ゲルストナー波の法線を計算してみます。安直な方法としては、CalcGerstnerWaveOffset関数が波のあらゆる地点の頂点位置の移動量を求められるものであることを利用し、x方向とz方向に1センチずらしてそれぞれCalcGerstnerWaveOffsetを呼び、これを三角形に見立てて法線を求めます。

	float3 ofsByGerstnerWave = CalcGerstnerWaveOffset(vertexPosition);
	float3 dx = float3(0.01, 0, 0) + CalcGerstnerWaveOffset(vertexPosition + float3(0.01, 0, 0));
	float3 dz = float3(0, 0, 0.01) + CalcGerstnerWaveOffset(vertexPosition + float3(0, 0, 0.01));
	float3 N = normalize(cross(dz - ofsByGerstnerWave, dx - ofsByGerstnerWave));

view raw calculate_gerstner_wave_normals.hlsl hosted with ❤ by GitHub

これでも実用上問題無いレベルになるとは思いますが、GPU GemsのEffective Water Simulation from Physical Modelsによると、Xで偏微分するとBinormal、Yで偏微分するとTangent、それらの外積からNormalを求めることができます。親切にもNormal, Binormal, Tangentそれぞれの式がEquation 10～12として提示されていますので、そのままHLSL化してみます。ただし、Island11はYが上、論文ではZが上なので、YとZを入れ替えてあります。そのせいで、Normal = Binormal X Tangentとはならず、Normal = Tangent X Binormalになっています。

	float3 CalcGerstnerWaveNormal(float3 P)
	{
	float3 normal = float3(0, 1, 0);
	[unroll]
	for (int i = 0; i < numWaves; i++)
	{
	Wave wave = waves[i];
	float wi = 2 / wave.waveLength;
	float WA = wi * wave.amplitude;
	float phi = speed * wi;
	float rad = wi * dot(wave.dir, P.xz) + phi * g_time;
	float Qi = steepness / (wave.amplitude * wi * numWaves);
	normal.xz -= wave.dir * WA * cos(rad);
	normal.y -= Qi * WA * sin(rad);
	}
	return normalize(normal);
	}
	float3 CalcGerstnerWaveBinormal(float3 P)
	{
	float3 binormal = float3(1, 0, 0);
	[unroll]
	for (int i = 0; i < numWaves; i++)
	{
	Wave wave = waves[i];
	float wi = 2 / wave.waveLength;
	float WA = wi * wave.amplitude;
	float phi = speed * wi;
	float rad = wi * dot(wave.dir, P.xz) + phi * g_time;
	float Qi = steepness / (wave.amplitude * wi * numWaves);
	binormal.x -= Qi * wave.dir.x * wave.dir.x * WA * sin(rad);
	binormal.z -= Qi * wave.dir.x * wave.dir.y * WA * sin(rad);
	binormal.y += wave.dir.x * WA * cos(rad);
	}
	return normalize(binormal);
	}
	float3 CalcGerstnerWaveTangent(float3 P)
	{
	float3 tangent = float3(0, 0, 1);
	[unroll]
	for (int i = 0; i < numWaves; i++)
	{
	Wave wave = waves[i];
	float wi = 2 / wave.waveLength;
	float WA = wi * wave.amplitude;
	float phi = speed * wi;
	float rad = wi * dot(wave.dir, P.xz) + phi * g_time;
	float Qi = steepness / (wave.amplitude * wi * numWaves);
	tangent.x -= Qi * wave.dir.x * wave.dir.y * WA * sin(rad);
	tangent.z -= Qi * wave.dir.y * wave.dir.y * WA * sin(rad);
	tangent.y += wave.dir.y * WA * cos(rad);
	}
	return normalize(tangent);
	}

view raw calculate_gerstner_wave_orientation.hlsl hosted with ❤ by GitHub

あえてくどくどと関数を分けてありますが、論文で使われている式や変数名をなるべく残すためで、現場で使う場合は１ループで一気にfloat3x3の回転行列まで求めると良いと思います。

Island11のドメインシェーダの既存の法線はそのままワールド空間の法線として使われていましたが、それをゲルストナー波の接空間（tangent space）とみなしてワールド空間に配置します。

	float3 gerstnerWaveN = CalcGerstnerWaveNormal(vertexPosition + ofsByGerstnerWave);
	float3 gerstnerWaveB = CalcGerstnerWaveBinormal(vertexPosition + ofsByGerstnerWave);
	float3 gerstnerWaveT = CalcGerstnerWaveTangent(vertexPosition + ofsByGerstnerWave);
	float3x3 rotator;
	rotator[0] = gerstnerWaveB;
	rotator[1] = gerstnerWaveN;
	rotator[2] = gerstnerWaveT;
	output.normal = mul(water_normal.xyz, rotator);

view raw convert_tangent_space_to_world.hlsl hosted with ❤ by GitHub

次回は水面に物を浮かべてみます。