Essentially Im wanting to project:
A 3D transformation matrix of a guassian(Or in this case for simplicity a unit sphere) that includes scale, rotation and translation in world space TO
A 2D transformation matrix for a unit ellipse in screen space
Im attempting to render guassians efficiently within the limitations of a games API. Its a unity game, and all I really have access to is meshes and shaders.
How do I go about doing this?
Im also wondering if this might be a question for mathmatics
Currently after caressing ChatGPT for awhile I managed to get it to spit out some snippets I managed to piece together, it almost works. It works perfectly with orthographic projection or perspective at distance, but perspective has major error up close. Heres my code:
float2x2 Project3DGaussianTo2D(
float3x3 gaussianLocalTransform,
float3x3 cameraRotationWorldToCamera,
float3 cameraRelativeDirection
) {
// Step 1: 3D Covariance Matrix
float3x3 G = gaussianLocalTransform;
float3x3 C = mul(G, transpose(G));
// Step 2: Transform covariance to camera/view space
float3x3 R = cameraRotationWorldToCamera;
float3x3 C_view = mul(R, mul(C, transpose(R)));
// Step 3: Projection Jacobian J
float3 center = mul(R, cameraRelativeDirection);
float x = center.x;
float y = center.y;
float z = center.z;
float invZ = 1.0 / z;
float invZ2 = invZ * invZ;
// Corrected Jacobian using perspective projection ∂proj/∂pos
float2x3 J;
J[0] = float3(invZ, 0, -x * invZ2);
J[1] = float3(0, invZ, -y * invZ2);
// Step 4: Project 3D covariance to 2D
float scaleZ = 1.0 / (z * z);
float2x2 C_2D = mul(J, mul(C_view * scaleZ, transpose(J)));
// Step 5: Eigen-decomposition to get ellipse axes
float a = C_2D[0][0], b = C_2D[0][1];
float c = C_2D[1][0], d = C_2D[1][1];
float T = a + d;
float D = a * d - b * c;
float disc = sqrt(max(T*T - 4*D, 0));
float lambda1 = 0.5 * (T + disc);
float lambda2 = 0.5 * (T - disc);
// Handle eigenvectors
float2 majorAxis = normalize(float2(b, lambda1 - a));
if (dot(majorAxis, majorAxis) < 1e-5)
majorAxis = normalize(float2(lambda1 - d, c));
float2 minorAxis = float2(majorAxis.y, -majorAxis.x);
float2 scale = float2(sqrt(max(lambda1, 0)), sqrt(max(lambda2, 0)));
float2x2 basis = float2x2(majorAxis, minorAxis);
return mul(basis, float2x2(scale.x, 0, 0, scale.y));
}
I dont really know how it all works if Im honest. I just need this magic sauce to work.
And rendered in unity using a plane mesh with uv matching plane, scaled down to zero in Blender, and using shader:
v2f vert(appdata v) {
v2f o;
float2 uv = float2(saturate(v.uv.x), v.uv.y);
uv = uv - 0.5;
o.uv = uv*2;
float3x3 sphereMatRot = ...;
float3x3 sphereMatScale = ...;
float3 cameraTangentX = mul(UNITY_MATRIX_V[0].xyz, transpose((float3x3)unity_WorldToObject));
float3 cameraTangentY = mul(UNITY_MATRIX_V[1].xyz, transpose((float3x3)unity_WorldToObject));
float3 cameraNormal = mul(UNITY_MATRIX_V[2].xyz, transpose((float3x3)unity_WorldToObject));
float3x3 cameraRot = ExtractRotation(UNITY_MATRIX_V);
// ===== PROJECT SPHERE TO ECLIPSE
float2x2 ellipseTrans = Project3DGaussianTo2D(
mul(sphereMatRot, sphereMatScale),
cameraRot,
cameraNormal
);
uv = mul(ellipseTrans, uv.xy);
float3 offset = uv.x * cameraTangentX + uv.y * cameraTangentY;
o.pos = UnityObjectToClipPos(v.vertex + float4(offset, 0));
return o;
}
float4 frag(v2f i) : SV_Target {
float circle = saturate(1-length(i.uv));
if (circle < 0.02) discard;
return float4(0,0,1,1);
}
```
