Shader "CrispyPin/Lenia"
{
	Properties
	{
		_LastFrame ("Texture", 2D) = "white" {}
		_GrowtCenter ("Growth fn center (mu)", Range(0, 1)) = 0.2
		_GrowthWidth ("Growth fn width (sigma / std deviation)", Range(0, 1)) = 0.07
		_Speed ("Speed factor", Range(0.001, 0.5)) = 0.1
	}
	SubShader
	{
		Tags { "RenderType"="Opaque" }
		LOD 100

		Pass
		{
			CGPROGRAM
			#pragma vertex vert
			#pragma fragment frag

			#include "UnityCG.cginc"

			struct appdata
			{
				float4 vertex : POSITION;
				float2 uv : TEXCOORD0;
			};

			struct v2f
			{
				float2 uv : TEXCOORD0;

				float4 vertex : SV_POSITION;
			};

			sampler2D _LastFrame;
			float _GrowtCenter;
			float _GrowthWidth;
			float _Speed;

			v2f vert (appdata v)
			{
				v2f o;
				o.vertex = UnityObjectToClipPos(v.vertex);
				o.uv = v.uv;
				return o;
			}

			float kernel(float p) {
				float r = p / 8.0;
				float k_sharpness = 28;
				float k_offset = 0.435;
				return exp(-((r - k_offset) * (r - k_offset)) * k_sharpness);
				// float r = p / _Radius;
				// // -- normal
				// return exp(-((r - _KOffset)*(r - _KOffset)) * _KSharpness);

				// -- example from the lenia paper
				// if (r >= 1) return 0;
				// const float alpha = 4;
				// return exp(alpha - alpha/(4.0 * r * (1.0 - r)));
			}

			inline float activation(float neighbors) {
				// return old_state * (neighbors > 2 && neighbors < 5) +
				// ((1 - old_state) * neighbors == 3);

				// const float sharpness = 1000.0;
				// float x = neighbors - _GrowtCenter;
				// return exp(-(x*x) * _GrowthWidth) * 2.0 - 1.0;
				const float mu = _GrowtCenter;
				const float sigma = _GrowthWidth;
				const float u = neighbors;
				return exp(-((u-mu) * (u-mu)) / (2 * sigma * sigma)) * 2.0 - 1.0;
			}



			inline half value(float2 center, float x, float y) {
				return tex2D(_LastFrame, center + float2(x, y)).r;
			}

			fixed4 frag (v2f i) : SV_Target
			{
				if(_ProjectionParams.z > 1) discard;

				const float resolution = 512.0;
				const float d = 1.0 / resolution;
				// Defines RADIUS Radius Kernel total_max
				#include "lenia_generated_kernel.cginc"

				float total = 0.0;
				[unroll(RADIUS)]
				for (int y = 0; y < Radius; y++) {
					[unroll(RADIUS)]
					for (int x = 1; x <= Radius; x++) {
						const float xx = (float)x * d;
						const float yy = (float)y * d;
						total += value(i.uv, xx, yy) * Kernel[y][x-1];
						total += value(i.uv, -yy, xx) * Kernel[y][x-1];
						total += value(i.uv, -xx, -yy) * Kernel[y][x-1];
						total += value(i.uv, yy, -xx) * Kernel[y][x-1];
					}
				}
				// */
				/*
				float total_max = 0;
				float total = 0;
				for (int x = -Radius; x <= Radius; x++) {
					for (int y = -Radius; y <= Radius; y++) {
						float dist = sqrt(x*x+y*y);
						float kval = kernel(dist);
						total_max += kval;
						total += value(i.uv, x*d, y*d) * kval;
					}
				}
				// */
				float old_state = value(i.uv, 0.0, 0.0) ;
				float count = total / total_max;

				float state = activation(count) * _Speed + old_state;
				state = clamp(state, 0, 1);

				// kernel visualization: lookup table (SLOW)
				// float k = 0;
				// {
				// 	float2 p = (i.uv - 0.5 ) * resolution;
				// 	p = floor(p);
				// 	if (p.x > 0 && p.y >= 0) {
				// 		k = Kernel[p.y][p.x-1];
				// 	} else if (p.x <= 0 && p.y > 0) {
				// 		k = Kernel[-p.x][p.y-1];
				// 	} else if (p.x < 0 && p.y <= 0) {
				// 		k = Kernel[-p.y][-p.x-1];
				// 	} else if (p.x >= 0 && p.y < 0) {
				// 		k = Kernel[p.x][-p.y-1];
				// 	}
				// }

				// kernel visualisation: real size
				// float2 p = (i.uv - 0.5) * resolution;
				// float k = kernel(length(p)) * (max(abs(p.x), abs(p.y)) <= _Radius);

				// kernel visualisation: fill square
				// float k = kernel(length(i.uv - 0.5) * _Radius * 2);

				// float a = activation(i.uv.x);
				// float4 col = float4(state, k, a, 1);

				float4 col = float4(state, i.uv.x * state, i.uv.y * state, 1);
				return col;
			}
			ENDCG
		}
	}
}