//	NEB_GPU_Functions.metal
//
//	© 2026 by Jeff Weeks
//	See TermsOfUse.txt

#include <metal_stdlib>
using namespace metal;

//	Include "NEB_GPU_Definitions.h" *after* "using namespace metal"
//	if you might need to use uint32_t in "NEB_GPU_Definitions.h"
#include "../../shared code/GGVR_GPU_Definitions.h"
#include "../../shared code/GGVR_GPU_Functions.h"
#include "NEB_GPU_Definitions.h"


//	For future compatibility with Euclidean-only platforms
//	like Unity, we'll map the vertices into Euclidean space.
//	If the player wants to perceive each vertex at the same
//	distance that a native inhabitant of the curved space
//	would see it, we'll use the azimuthal equidistant projection;
//	otherwise we'll use the (stretched) central projection.
constant uchar	gProjectionCode			[[ function_constant(0) ]];	//	 8-bit uint

constant uchar	gGeometryCode			[[ function_constant(1) ]];	//	 8-bit uint
constant float	gSpaceRadius			[[ function_constant(2) ]];	//	 32-bit float
constant float	gSpaceRadiusInv = 1.0 / gSpaceRadius;	//	avoid division at run time

//	The balls and the cue stick get special lighting.
constant bool	gBilliardBallShading	[[ function_constant(3) ]];	//	 8-bit boolean
constant bool	gCueStickShading		[[ function_constant(4) ]];	//	 8-bit boolean

//	The billiards table uses a procedural texture, which of course isn't mipmapped.
//	Left unmodified, the table would flicker wherever the distance between
//	pixels become comparable to or greater than the procedural texture's
//	feature size. To avoid such flickering, we'll reduce the contrast
//	in the table's more distant and/or more slanted regions.
//	No meshes other than the billiards table use gAntiflickering.
constant bool	gAntiflickering		[[ function_constant(5) ]];	//	 8-bit boolean

//	Rendering a pocket clears the depth buffer underneath it.
constant bool	gClearDepthBuffer	[[ function_constant(6) ]];	//	 8-bit boolean

//	The checkboard tiling of polyhedra (used to create the lattice)
//	requires instanced rendering, but nothing else does.
constant bool	gInstancedRendering	[[ function_constant(7) ]];	//	 8-bit boolean

//	Only the checkboard tiling of polyhedra uses fog.
constant bool	gFog				[[ function_constant(8) ]];	//	 8-bit boolean
constant float	gFogCutoffInMeters	[[ function_constant(9)	]];	//	32-bit float
constant half	gFogCutoffInvSq = half(1.0 / (gFogCutoffInMeters * gFogCutoffInMeters));

//	When rendering S³ via central projection we must
//	clip away the southern hemisphere, for two reasons:
//	first, because the central projection doesn't include it
//	and second, because it would force the player to look
//	walleyed to focus on the objects that it contains,
//	and we don't want to impose that discomfort (and risk)
//	on our players.
//
//		Technical note: I would like to write
//
//		  constant bool	gClipSouthernHemisphere
//						= gProjectionCode == 'c' // central projection
//					   && gGeometryCode == 's';  // spherical geometry
//
//		and use it to enable or disable the custom clip_distance,
//		but Metal doesn't handle that reliably. In some apps it works,
//		in others it doesn't. So instead let's go with the fully
//		reliable approach of defining a separate (but redundant)
//		function_constant for this task.
//
constant bool	gClipSouthernHemisphere	[[ function_constant(10) ]];	//	 8-bit boolean


struct VertexInput
{
	float4	pos [[ attribute(VertexAttributePosition)	]];
	float3	tex [[ attribute(VertexAttributeTexCoords)	]];
};

struct VertexOutput
{
	float4	position		[[ position			]];
	float3	texCoords		[[ user(texcoords)  ]];	//	(u,v,-) for regular texture
													//	  or (u,v,w) for cube map
	half	brightness		[[ user(brightness)	]];
	half	constrast		[[ user(constrast), function_constant(gAntiflickering) ]];
	float	clipDistance	[[ clip_distance, function_constant(gClipSouthernHemisphere)]] [1];

};
struct FragmentInput
{
	float3	texCoords	[[ user(texcoords)	]];	//	(u,v,-) for regular texture
												//	  or (u,v,w) for cube map
	half	brightness	[[ user(brightness)	]];
};
struct BilliardsTableFragmentInput
{
	float3	texCoords	[[ user(texcoords)	]];	//	(u,v,-)
	half	constrast	[[ user(constrast)  ]];
};

half	textureContrast(float4 positionInBridgeOfNoseSpace, float4 u);


vertex VertexOutput nonEuclideanBilliardsVertexFunction(
				 VertexInput		in				 [[ stage_in							]],
				 ushort				amp_id			 [[ amplification_id	]],	//	strange name!
		constant GGVRViewProjection	& viewProjection [[ buffer(BufferIndexViewProjection)	]],
		constant simd_float4x4		& modelMatrix	 [[ buffer(BufferIndexModelMatrix)		]],
		constant simd_half3			& upVector		 [[ buffer(BufferIndexUpVector),
													function_constant(gBilliardBallShading)	]],
													//	dimensionless, unit length,
  const device   simd_float4x4		*tilingMatrices  [[ buffer(BufferIndexTilingMatrices),
													function_constant(gInstancedRendering)	]],
				 ushort				iid				 [[ instance_id,
													function_constant(gInstancedRendering)	]]
)
{
	VertexOutput	out;
	float4			positionInTableSpace,
					transformedPositionInTableSpace,
					positionInBridgeOfNoseSpace,	//	in (typically curved) space; in meters
					positionInTangentSpace,	//	in meters
					positionInClipSpace;
	float			centralProjectionFactor;
	half			distanceSquared,
					fogFactor;

	//	When gInstancedRendering is disabled, it would be slightly
	//	more efficient to pass in a single pre-composed modelView matrix,
	//	but the performance benefit would be negligible, so I've left
	//	the model and view matrices separate to keep the Swift code
	//	simpler and clearer.
	positionInTableSpace = modelMatrix * in.pos;
	if (gInstancedRendering) {
		transformedPositionInTableSpace
			= tilingMatrices[iid]
			* positionInTableSpace;
	}
	else
	{
		transformedPositionInTableSpace = positionInTableSpace;
	}
	positionInBridgeOfNoseSpace = viewProjection.gameSpaceToBridgeOfNoseSpace
								   * transformedPositionInTableSpace;

	if (gClipSouthernHemisphere) {
		out.clipDistance[0] = positionInBridgeOfNoseSpace.w;
	}

	//	See comment accompanying the definition of gAntiflickering.
	if (gAntiflickering)
	{
		out.constrast = textureContrast(
							positionInBridgeOfNoseSpace,
							viewProjection.gameSpaceToBridgeOfNoseSpace[1]);
	}

	switch (gProjectionCode) {
	
	case 'c':

		if (positionInBridgeOfNoseSpace.w > 0.0) {
		
			//	H³, E³ or near hemisphere of S³
		
			centralProjectionFactor = gSpaceRadius / positionInBridgeOfNoseSpace.w;
			positionInTangentSpace = positionInBridgeOfNoseSpace * centralProjectionFactor;

		} else {
		
			//	Equator or far hemisphere of S³
		
			//	Even though this vertex will ultimately get clipped away,
			//	we must assign it a plausible position anyhow,
			//	because any triangle that contains both this vertex
			//	and a near-hemisphere vertex will need to get rasterized.
			//	So let's position this vertex in the correct direction,
			//	at some large distance.

			positionInTangentSpace = float4(
										1024.0 * normalize(positionInBridgeOfNoseSpace.xyz),
										gSpaceRadius);
		}

		//	The point now sits in the Euclidean tangent space.
		//	For mathematical perfection, we'd need to apply
		//	the correction factors here. In practice, though, they are
		//	always equal to 1.0 to within a fraction of a percent,
		//	so we may ignore them without compromising the quality
		//	of the rendered image. For more details, please see the comment
		//	accompanying the definition of GGVRAzimuthalProjection.central .
// Uncomment these three lines to enabled the mathematically perfect algorithm.
//		positionInTangentSpace[0] *= viewProjection.xFactor;
//		positionInTangentSpace[1] *= viewProjection.yzFactor;
//		positionInTangentSpace[2] *= viewProjection.yzFactor;
		
		break;
	
	case 'e':

		//	See comment accompanying the definition of GGVRAzimuthalProjection.equidistant.

		positionInTangentSpace = azimuthalEquidistantMap(
									positionInBridgeOfNoseSpace,
									gGeometryCode,
									gSpaceRadius,
									gSpaceRadiusInv);

		break;
	}

	distanceSquared = dot(half3(positionInTangentSpace.xyz),
						  half3(positionInTangentSpace.xyz));

	positionInClipSpace
		= viewProjection.tangentSpaceToClipSpace[amp_id]
		* positionInTangentSpace;
	
	if (gClearDepthBuffer)
	{

		positionInClipSpace.z = 0.0;

	} else if (gProjectionCode == 'e'
			&& gGeometryCode == 's')
	{

		//	When rendering in spherical geometry,
		//	the azimuthal equidistant projection ensures that
		//	the player never needs to look walleyed, but it has
		//	the side effect that a short geodesic near the south pole
		//	can get mapped to a broad arc that's concave towards
		//	the player in the Euclidean tangent space.
		//	The standard rendering algorithm would approximate
		//	that arc as a straight line that "takes a shortcut"
		//	from one endpoint to the other, instead of following
		//	the concave arc. This shortcut is usually harmless,
		//	because it doesn't affect what the arc looks like
		//	in the player's field of view, it only affects
		//	the depth coordinate (the arc is concave in the radial direction,
		//	not the sideways direction). However, the shortcut
		//	sometimes causes a more distant edge to pass in front
		//	of a slightly closer one.
		//
		//	To avoid that artefact, when using an equidistant projection
		//	in spherical geometry we use radial depth instead
		//	of the usual planar depth. In principle there's some risk
		//	that radial depth could create new artefacts -- same effect
		//	described above, but in reverse -- but in practice
		//	this doesn't happen in any noticeable way.
		//
		//		Technical note:  We compute the radial depth
		//		relative to the bridge of the player's nose,
		//		not relative to each eye separately. This
		//		ensures that the intersection of two meshes
		//		is always a well-defined curve in space.
		//
		const float	piInverse = 0.31830988618379067154;
		const float normalizationFactor = gSpaceRadiusInv * piInverse;	//	1/(π*ρ)
		float normalizedDistance = sqrt(distanceSquared) * normalizationFactor;

		//	Metal will divide the z coordinate by positionInClipSpace.w
		//	so here we multiply by positionInClipSpace.w to compensate.
		
		positionInClipSpace.z = (1.0 - normalizedDistance) * positionInClipSpace.w;
	}
	
	out.position = positionInClipSpace;
	out.texCoords = in.tex;
	
	//	To do non-Euclidean lighting properly, we'd have
	//	to work with specific light sources (point lights
	//	or maybe more extended lamps). The Euclidean concept
	//	of a directional light generalizes poorly,
	//	especially in the 3-sphere. And to get it all
	//	to look anywhere near as good as RealityKit's
	//	physically-based lighting would be a huge amount
	//	of work on my part, and would require a lot more
	//	GPU time as well.
	//
	//	So let's just put some simple shading on the billiard balls
	//	to mimic "ambient occlusion" and declare that to be
	//	good enough. Even though it will be far from photorealistic,
	//	the scene will look convincingly 3D.
	//
	//	Similarly, let the cue stick be brighter along its upper surfaces
	//	and darker along its lower surfaces. It should always be
	//	roughly horizontal, so with luck it'll look good.
	//
	if (gBilliardBallShading)
	{
		//	The texture coordinates provide a purely 3-dimensional
		//	outward-pointing unit normal at each vertex, which is
		//	exactly what we need to compute a good-looking gradient
		//	from full brightness at the top to zero brightness
		//	at the bottom (due to "ambient occlusion").
		half ambientBrightness = 0.5h + 0.5h * dot( half3(in.tex), upVector );

		out.brightness = ambientBrightness;
	}
	else if (gCueStickShading)
	{
		//	The second texture coordinate has value
		//
		//		 0	along one side of the cue stick
		//		1/4	along the top of the cue stick
		//		1/2	along the other side of the cue stick
		//		3/4 along the bottom of the cue stick
		//
		out.brightness = 0.5h + 0.5h * sin(half(in.tex[1]) * 2.0h * M_PI_H);
	}
	else if (gFog)
	{
		switch (gGeometryCode) {
		
		case 's':
			//	0.3 is an aesthetic choice
			fogFactor = 0.3 * acos(positionInBridgeOfNoseSpace.w * gSpaceRadiusInv);
			break;
		
		case 'e':
		case 'h':
			fogFactor = min(distanceSquared * gFogCutoffInvSq, 1.0h);
			break;
		}
		out.brightness = 1.0h - fogFactor;
	}
	else
	{
		out.brightness = 1.0h;
	}

	return out;
}

half textureContrast(
	float4 positionInBridgeOfNoseSpace,	//	in (typically curved) space; in meters
	float4 u)	//	see definition of u below
{
	//	Let the fragment shader know how strongly it should
	//	reduce contrast to avoid flicker, as explained
	//	in the comment accompanying the definition of gAntiflickering.

	//	The player's lines of sight passing through a given pixel
	//	define a very thin pyramid P. Where that pyramid intersects
	//	the billiard table's surface, we want to know how the size
	//	of the intersection compares to the procedural texture's
	//	feature size.
	//
	//	Different pixels may subtend different solid angles,
	//	but given the fancy post-procesing that visionOS applies
	//	to the image (to account for the lens shape, etc.) we
	//	really don't know how the subtended angle varies
	//	from pixel to pixel. So instead we make the simplifying
	//	assumption that every pixel defines a pyramid P
	//	of solid angle ε steradians. In other words, we'll
	//	assume that each pixel fill a constant fraction ε/(4π)
	//	of the player's visual sphere.
	//
	//	Let S be a sphere centered at the player's eye and
	//	containing the positionInBridgeOfNoseSpace. We may easily
	//	compute the "effective radius" r of S (not to be confused
	//	with its intrinsic radius as measured in the space itself).
	//
	float r = length(positionInBridgeOfNoseSpace.xyz);
	
	//	When the player looks at the table obliquely,
	//	the pyramid P encloses a larger portion of the table's surface
	//	than it would if the player were looking orthogonally
	//	at the table from an equal distance. More precisely,
	//	if the player looks at the table from some angle θ, where
	//
	//		θ =  0  when the player views the table orthogonally and
	//		θ = π/2 when the player views the table tangentially,
	//
	//	then the table area enclosed by the pyramid P gets
	//	stretched by a factor of 1/cos(θ).
	//
	//		Note: Even though the enclosed table area gets
	//		stretched in only one dirction, we'll nevertheless
	//		multiply r by 1/cos(θ) to allow for the direction in which
	//		the pixel spacing along the table will be greatest.
	//
	//	Let
	//		u be the table's upward-pointing unit normal vector and
	//		v be a unit vector pointing towards the player's head.
	//	We may compute
	//		cos(θ) = <u,v>
	//	where <,> is the inner product appropriate for the space.
	
	//	In table coordinates, u_table = (0.0, 1.0, 0.0, 0.0).
	//	To convert to bridge-of-nose coordinates we could write
	//
	//		u = viewProjection.gameSpaceToBridgeOfNoseSpace * u_table
	//
	//	because the table's modelMatrix is the identity.
	//	While the compiler might be smart enough to simplify
	//	that matrix multiplication, let's go ahead and simplify
	//	it ourselves as
	//
	//		float4 u = viewProjection.gameSpaceToBridgeOfNoseSpace[1];
	//
	//	The caller has already passed us this vector u,
	//	to avoid passing the whole gameSpaceToBridgeOfNoseSpace.
	
	//	Compute the vector v.
	
	//	Let pos_unit be the position vector normalized to unit length.
	float4 pos_unit = positionInBridgeOfNoseSpace * gSpaceRadiusInv;
	
	//	Split pos_unit into its xyz and w parts.
	float3 theXYZPart	= pos_unit.xyz;
	float  theXYZLength	= length(theXYZPart);
	float  theWPart		= pos_unit.w;
	float  theWLength	= abs(theWPart);
	float4 v;
	if (theXYZLength < 0.0001)
	{
		//	The point is extremely close to the player's head,
		//	so assign v some arbitrary direction.
		v = float4(1.0, 0.0, 0.0, 0.0);
	}
	else
	{
		//	Compute the xyz and w parts of v separately, then combine them.
		float3 v_xyz;
		float  v_w;
		
		if (gGeometryCode == 's')
		{
			if (theWLength < 0.0001)
			{
				//	The point is near the player's equatorial 2-sphere,
				//	so v points in the w direction.
				v_xyz = float3(0.0, 0.0, 0.0);
				v_w = 1.0;
			}
			else
			{
				//	Take pos_unit and swap the lengths of its xyz and w parts,
				//	while negating the xyz part.
				v_xyz = -(theWLength / theXYZLength) * theXYZPart;
				v_w   =  (theXYZLength / theWLength) * theWPart;
			}
		}
		if (gGeometryCode == 'e')
		{
			v_xyz = -(theWLength / theXYZLength) * theXYZPart;	//	theWLength = 1.0
			v_w = 0.0;
		}
		if (gGeometryCode == 'h')
		{
			//	Take pos_unit and swap the lengths of its xyz and w parts,
			//	while negating both the xyz part and the w part.
			v_xyz = -(theWLength / theXYZLength) * theXYZPart;
			v_w   = -(theXYZLength / theWLength) * theWPart;
		}
		
		v = float4(v_xyz, v_w);
	}
	
	//	Compute cos(θ) = <u,v>.
	float cosTheta;
	if (gGeometryCode == 's')
	{
		cosTheta = u.x * v.x
				 + u.y * v.y
				 + u.z * v.z
				 + u.w * v.w;
	}
	if (gGeometryCode == 'e')
	{
		cosTheta = u.x * v.x
				 + u.y * v.y
				 + u.z * v.z;
				// u.w = v.w = 0.0
	}
	if (gGeometryCode == 'h')
	{
		cosTheta = u.x * v.x
				 + u.y * v.y
				 + u.z * v.z
				 - u.w * v.w;
	}
	
	//	Compute the slope factor 1/cos(θ) and apply it.
	float absCosTheta = abs(cosTheta);
	float theSlopeFactor = ( absCosTheta > 0.01 ? 1.0 / absCosTheta : 100.0);
	float rAdjusted = r * theSlopeFactor;

	//	Let rThreshold be the adjusted distance at which the muting begins
	//	and rSaturation be the adjusted distance at which the muting reaches 100%.
	const float rThreshold = 0.5;					//	in meters
	const float rSaturation = 1.0;					//	in meters
	const float rSpan = rSaturation - rThreshold;	//	in meters
	const float rSpanInv = 1.0 / rSpan;				//	in meters⁻¹

	half theContrast = clamp(half((rSaturation - rAdjusted) * rSpanInv), 0.0h, 1.0h);
	
	return theContrast;
}

float simplexNoise(float2 uv);

fragment half4 billiardsTableFragmentFunction(
	BilliardsTableFragmentInput	in	[[ stage_in	]] )
{
	const half3 theTableColor    = half3(0.00, 1.00, 0.25);
	const half3 theContrastColor = half3(0.00, 0.00, 0.00);
	const half3 theAverageColor  = 0.5h * (theTableColor + theContrastColor);
	const half3 theDelta         = theTableColor - theAverageColor;

	half3 color;
	
	//	Performance note:  Typically either all the threads
	//	in a SIMD group have contrast == 0.0 or none of them do.
	//	If all do, they can return theAverageColor immediately,
	//	completely avoiding the relatively expensive call to simplexNoise().
	//	(But if only some of the threads in a SIMD group
	//	have contrast == 0.0, they must "follow along" as
	//	the other threads execute the call simplexNoise().)

	if (in.constrast == 0.0h)
	{
		color = theAverageColor;
	}
	else
	{
		float2 uv = in.texCoords.xy;

		//	Get a noise value in the range (-1.0 ... +1.0)
		//	based on the texture coordinates.
		half t = half( simplexNoise(uv) );
		
		//	Multiply t by the constrast, so that the more distant parts
		//	of the table, where flickering would be a problem,
		//	will have less contrast.
		t *= in.constrast;
		
		color = theAverageColor + t * theDelta;
	}
	
	return half4(color, 1.0h);
}

//	Some useful functions
float2 mod289(float2 x) { return x - floor(x * (1.0/289.0)) * 289.0; }
float3 mod289(float3 x) { return x - floor(x * (1.0/289.0)) * 289.0; }
float3 permute(float3 x) { return mod289(((x * 34.0) + 1.0) * x); }

float simplexNoise(float2 uv)	//	Seems to return values in the range (-1.0 ... +1.0).
								//	(comment added by JRW)
{
	//	References
	//
	//		https://en.wikipedia.org/wiki/Simplex_noise
	//		https://thebookofshaders.com/11/
	//		https://thebookofshaders.com/edit.php#11/2d-gnoise.frag
	//
	
	//	This code is based on the "GLSL 2D simplex noise function"
	//	which comes with the description
	//
	//		Description : GLSL 2D simplex noise function
	//		     Author : Ian McEwan, Ashima Arts
	//		 Maintainer : ijm
	//		    Lastmod : 20110822 (ijm)
	//		    License :
	//		Copyright (C) 2011 Ashima Arts. All rights reserved.
	//		Distributed under the MIT License. See LICENSE file.
	//		https://github.com/ashima/webgl-noise
	//
	//	which refers to the license
	//
	//		Copyright (C) 2011 by Ashima Arts (Simplex noise)
	//
	//		Permission is hereby granted, free of charge, to any person obtaining a copy
	//		of this software and associated documentation files (the "Software"), to deal
	//		in the Software without restriction, including without limitation the rights
	//		to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
	//		copies of the Software, and to permit persons to whom the Software is
	//		furnished to do so, subject to the following conditions:
	//
	//		The above copyright notice and this permission notice shall be included in
	//		all copies or substantial portions of the Software.
	//
	//		THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
	//		IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
	//		FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
	//		AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
	//		LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
	//		OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
	//		THE SOFTWARE.
	//
	
	//	That GLSL source code includes the lines
	//
	//		#ifdef GL_ES
	//		precision mediump float;
	//		#endif
	//
	//	A quick online search suggests that "mediump" in GLSL
	//	typically means "half" precision in Metal.
	//	Yet many of the constants are given to 15 digits
	//	(for example 0.211324865405187), which goes way beyond
	//	the 3 or 4 decimal digits that half precision supports
	//	and in fact well beyond the 7 digits that single-precision
	//	float supports. I'm guessing that the author left
	//	all those digits there to allow the code to be run
	//	with double-precision floats (although even there
	//	it might just be for the satisfaction clear code --
	//	as a practical matter we are, after all, just computing noise).
	//
	//	That said... the code fails with half precision
	//	but succeeds with single-precision floats.

	//	Precomputed values for skewed triangular grid
	const float4 c = float4(
				 0.211324865405187,		//	(3.0 - sqrt(3.0))/6.0
				 0.366025403784439,		//	0.5*(sqrt(3.0) - 1.0)
				-0.577350269189626,		//	-1.0 + 2.0 * c.x = sqrt(1.0/3.0)
				 0.024390243902439);	//	1.0 / 41.0

	//	First corner (x0)
	float2  i = floor(uv + dot(uv, c.yy));
	float2 x0 = uv - i + dot(i, c.xx);

	//	Other two corners (x1, x2)
	float2 i1 = (x0.x > x0.y ? float2(1.0, 0.0) : float2(0.0, 1.0));
	float2 x1 = x0.xy + c.xx - i1;
	float2 x2 = x0.xy + c.zz;

	//	Do some permutations to avoid
	//	truncation effects in permutation [JRW says:  WTF?!? ]
	i = mod289(i);
	float3 p = permute(
			permute( i.y + float3(0.0, i1.y, 1.0))
				+ i.x + float3(0.0, i1.x, 1.0 ));

	float3 m = max(
				0.0,
				0.5 - float3(dot(x0,x0), dot(x1,x1), dot(x2,x2)));
	m = m * m;
	m = m * m;

	//	Gradients:
	//	 41 pts uniformly over a line, mapped onto a diamond
	//	 The ring size 17*17 = 289 is close to a multiple
	//	  of 41 (41*7 = 287)

	float3 x = 2.0 * fract(p * c.www) - 1.0;
	float3 h = abs(x) - 0.5;
	float3 ox = floor(x + 0.5);
	float3 a0 = x - ox;

	//	Normalise gradients implicitly by scaling m
	//	Approximation of: m *= inversesqrt(a0*a0 + h*h);
	//
	//		Note by JRW:  It seems this approximation
	//		is required to keep the return values
	//		in the range (-1.0 ... +1.0). When I swapped in
	//		the exact division by sqrt(a0*a0 + h*h)
	//		the return values spread out to a range
	//		something like (-1.4 ... +1.4).
	//
	m *= 1.79284291400159 - 0.85373472095314 * (a0*a0 + h*h);

	//	Compute final noise value at P
	float3 g = float3(0.0);
	g.x  = a0.x * x0.x  +  h.x * x0.y;
	g.yz = a0.yz * float2(x1.x,x2.x) + h.yz * float2(x1.y,x2.y);

	return 130.0 * dot(m, g);
}

fragment half4 groundingShadowFragmentFunction(
	FragmentInput	in	[[ stage_in	]] )
{
	const half3 blackColor = half3(0.0h, 0.0h, 0.0h);

	float2 uv = in.texCoords.xy;

	half r = uv.x;	//	relative radial distance, in range 0.0 ≤ r ≤ 1.0
	half opacity = 1.0h - r * r;

	return half4(blackColor, opacity);
}

fragment half4 pocketFragmentFunction()
{
	return half4(0.0h, 0.0h, 0.0h, 1.0h);	//	opaque black
}

fragment half4 billiardBallFragmentFunction(
	FragmentInput		in				[[ stage_in							]],
	texturecube<half>	cubeMapTexture	[[ texture(TextureIndexPrimary)		]],
	sampler				textureSampler	[[ sampler(SamplerIndexPrimary)		]],
	constant half		& saturation	[[ buffer(BufferIndexColoration)	]]	)
{
	half4	tmpTexelColor;
	half	tmpSaturation;
	half3	tmpModifiedColor,
			tmpShadedColor;

	tmpTexelColor = cubeMapTexture.sample(textureSampler, in.texCoords.xyz);

	tmpSaturation = half(saturation);
	if ( tmpSaturation < 1.0h )
	{
		//	Convert the texel color to greyscale.
		
		//	The luminance ("Y" value) coefficients
		//	for the linear sRGB color space are
		//
		//		(0.2126, 0.7152, 0.0722)
		//
		//	But visionOS uses linear Display P3, so we need to multiply
		//	those coefficients by the P3-to-sRGB conversion matrix
		//	as computed in GeometryGamesColor.swift, namely
		//
		//	gP3toSRGB = simd_double3x3(
		//		SIMD3<Double>( 1.22494017628056, -0.04205695470969, -0.01963755459033),
		//		SIMD3<Double>(-0.22494017628056,  1.04205695470969, -0.07863604555063),
		//		SIMD3<Double>( 0.00000000000000,  0.00000000000000,  1.09827360014097))
		//
		//	to get the luminance ("Y" value) coefficients for Display P3, which are
		//
		//		(0.2289, 0.6918, 0.0793)
		//
		half3	tmpLuminanceCoefs = half3(0.2289h, 0.6918h, 0.0793h);
		
		//	As well as converting to greyscale, let's fade
		//	towards white (or some other fixed shade of grey).
		half	a = 0.5h,
				b = 1.0h - a,
				fadeValue = 1.0;	//	= white
		
		//	Let's write the luminance computation and the fade
		//	as a single expression, to encourage the compiler
		//	to pre-multiply the various constants at compile time.
		
		half greyValue
		//	= a * dot(tmpLuminanceCoefs, tmpTexelColor.rgb)
			= dot(a * tmpLuminanceCoefs, tmpTexelColor.rgb)
			+ b * fadeValue;	//	or
								//		b * fadeValue * tmpTexelColor.a
								//	when alpha ≠ 1.0, to respect pre-multiplied alpha.

		tmpModifiedColor
			= (1.0h - tmpSaturation) * half3(greyValue, greyValue, greyValue)
			+     tmpSaturation      * tmpTexelColor.rgb;
	}
	else	//	saturation = 1.0
	{
		tmpModifiedColor = tmpTexelColor.rgb;
	}

	tmpShadedColor = in.brightness * tmpModifiedColor;

	return half4(tmpShadedColor, tmpTexelColor.a);
}

fragment half4 cueStickFragmentFunction(
	FragmentInput	in				[[ stage_in							]],
	texture2d<half>	texture			[[ texture(TextureIndexPrimary)		]],
	sampler			textureSampler	[[ sampler(SamplerIndexPrimary) 	]],
	constant half	& opacity		[[ buffer(BufferIndexColoration)	]])
{
	half4	tmpTexelColor,
			tmpOpaqueColor,
			tmpFragmentColor;

	tmpTexelColor = texture.sample(textureSampler, in.texCoords.xy);

	tmpOpaqueColor = half4(in.brightness, in.brightness, in.brightness, 1.0h)
				   * tmpTexelColor;
	
	tmpFragmentColor = opacity * tmpOpaqueColor;	//	pre-multiplied alpha

	return tmpFragmentColor;
}

fragment half4 handFragmentFunction(
	FragmentInput	in		[[ stage_in							]],
	constant half4	& color	[[ buffer(BufferIndexColoration)	]])
{
	return color;	//	pre-multiplied alpha
}

fragment half4 checkerboardTilingFragmentFunction(
	FragmentInput	in		[[ stage_in	]],
	constant half3	& color	[[ buffer(BufferIndexColoration) ]])
{
	half3 foggedColor;
	
	foggedColor = in.brightness * color;
	
	return half4(foggedColor, 1.0h);
}
