Use Peter-Pike's shader code

ScriptableRenderPipeline/HDRenderPipeline/HDRP/Lighting/Volumetrics/Resources/VolumetricLighting.compute

     float3 radianceNoPhase;
 };
 
-// RGB SphericalHarmonics
-struct LightProbe
-{
-    float3 coeffs[9]; // 3 bands
-};
-
-LightProbe GetPreconvolvedAmbientProbe()
-{
-    LightProbe probe;
-
-    for (int i = 0; i < 9; i++)
-    {
-        probe.coeffs[i] = _AmbientProbeCoeffs[i].rgb;
-    }
-
-    return probe;
-}
-
-// Warning: since it's a SH, the results can be negative!
-float3 EvaluateLightProbe(LightProbe probe, float3 d)
-{
-    float3 radiance = 0;
-
-    // Coefficients are premultiplied and reordered in the MAD form:
-    // HornerForm[c_0 + c_1 y + c_2 z + c_3 x + c_4 x y + c_5 y z + c_6 (3 z^2 - 1) + c_7 x z + c_8 (x^2 - y^2)]
-    // = z (3 c_6 z + c_7 x + c_2) + y (-c_8 y + c_5 z + c_4 x + c_1) + x (c_8 x + c_3) + (c_0 - c_6)
-    return probe.coeffs[0] + d.x * (probe.coeffs[1] + d.x * probe.coeffs[2])
-                           + d.y * (probe.coeffs[3] + d.x * probe.coeffs[4] + d.z * probe.coeffs[5] - d.y * probe.coeffs[2])
-                           + d.z * (probe.coeffs[6] + d.x * probe.coeffs[7] + d.z * probe.coeffs[8]);
-}
-
 // Computes the light integral (in-scattered radiance) within the voxel.
 // Multiplication by the scattering coefficient and the phase function is performed outside.
 VoxelLighting EvaluateVoxelLighting(LightLoopContext context, uint featureFlags, PositionInputs posInput, float3 centerWS,
 
     // The contribution of the ambient probe does not depend on the position,
     // only on the direction and the length of the interval.
-    float3 probeInScatteredRadiance = EvaluateLightProbe(GetPreconvolvedAmbientProbe(), ray.centerDirWS);
+    // SampleSH9() evaluates the 3-band SH in a given direction.
+    float3 probeInScatteredRadiance = SampleSH9(_AmbientProbeCoeffs, ray.centerDirWS);
 
     float3 totalRadiance = 0;
     float  opticalDepth  = 0;

ScriptableRenderPipeline/HDRenderPipeline/HDRP/Lighting/Volumetrics/VolumetricLighting.cs

     public struct ZonalHarmonicsL2
     {
         public float[] coeffs; 
-    };
+    }
 
     public static ZonalHarmonicsL2 GetHenyeyGreensteinPhaseFunction(float asymmetry)
     {
         return sh;
     }
 
-    // Undoes coefficient normalization to obtain the canonical values of SH.
-    public static SphericalHarmonicsL2 DenormalizeSH(SphericalHarmonicsL2 sh)
+    // Undoes coefficient rescaling due to the convolution with the clamped cosine kernel
+    // to obtain the canonical values of SH.
+    public static SphericalHarmonicsL2 UndoCosineRescaling(SphericalHarmonicsL2 sh)
     {    
         float sqrtPi = Mathf.Sqrt(Mathf.PI);
 
         const float c4 = 0.13656855382400988382f; // 1/16 * sqrt(15/Pi)
 
         // Compute the inverse of SphericalHarmonicsL2::kNormalizationConstants.
-        // See SetSHEMapConstants() in "Stupid Spherical Harmonics Tricks". Note that we do not multiply by 3 here.
+        // See SetSHEMapConstants() in "Stupid Spherical Harmonics Tricks".
         float[] invNormConsts = { 1.0f / c0, -1.0f / c1, 1.0f / c1, -1.0f / c1, 1.0f / c2, -1.0f / c2, 1.0f / c3, -1.0f / c2, 1.0f / c4 };
 
         for (int c = 0; c < 3; c++)
     // Premultiplies the SH with the polynomial coefficients of SH basis functions,
     // which avoids using any constants during SH evaluation.
     // The resulting evaluation takes the form:
-    // c_0 + c_1 y + c_2 z + c_3 x + c_4 x y + c_5 y z + c_6 (3 z^2 - 1) + c_7 x z + c_8 (x^2 - y^2)
-    public static SphericalHarmonicsL2 PremultiplySH(SphericalHarmonicsL2 sh)
+    // (c_0 - c_6) + c_1 y + c_2 z + c_3 x + c_4 x y + c_5 y z + c_6 (3 z^2) + c_7 x z + c_8 (x^2 - y^2)
+    public static SphericalHarmonicsL2 PremultiplyCoefficients(SphericalHarmonicsL2 sh)
     {
         const float k0 = 0.28209479177387814347f; // {0, 0} : 1/2 * sqrt(1/Pi)
         const float k1 = 0.48860251190291992159f; // {1, 0} : 1/2 * sqrt(3/Pi)
 
     void SetPreconvolvedAmbientLightProbe(CommandBuffer cmd, float asymmetry)
     {
-        SphericalHarmonicsL2 probeSH = DenormalizeSH(RenderSettings.ambientProbe);
+        SphericalHarmonicsL2 probeSH = UndoCosineRescaling(RenderSettings.ambientProbe);
-        SphericalHarmonicsL2 finalSH = PremultiplySH(Convolve(probeSH, phaseZH));
+        SphericalHarmonicsL2 finalSH = PremultiplyCoefficients(Convolve(probeSH, phaseZH));
-        // Reorder coefficients in the MAD form:
-        // HornerForm[c_0 + c_1 y + c_2 z + c_3 x + c_4 x y + c_5 y z + c_6 (3 z^2 - 1) + c_7 x z + c_8 (x^2 - y^2)]
-        // = z (3 c_6 z + c_7 x + c_2) + y (-c_8 y + c_5 z + c_4 x + c_1) + x (c_8 x + c_3) + (c_0 - c_6)
-        Vector4[] coeffs = new Vector4[9];
+        // Pack coefficients so that we can use Peter-Pike Sloan's shader code.
+        // See SetSHEMapConstants() in "Stupid Spherical Harmonics Tricks".
+        Vector4[] coeffs = new Vector4[7];
-        const int r = 0, g = 1, b = 2;
+        // Constant + linear
+        for (int c = 0; c < 3; c++)
+        {
+            coeffs[c].x = finalSH[c, 3];
+            coeffs[c].y = finalSH[c, 1];
+            coeffs[c].z = finalSH[c, 2];
+            coeffs[c].w = finalSH[c, 0] - finalSH[c, 6];
+        }
-        coeffs[0] = new Vector3(finalSH[r, 0], finalSH[g, 0], finalSH[b, 0])
-                  - new Vector3(finalSH[r, 6], finalSH[g, 6], finalSH[b, 6]);
-        coeffs[1] = new Vector3(finalSH[r, 3], finalSH[g, 3], finalSH[b, 3]);
-        coeffs[2] = new Vector3(finalSH[r, 8], finalSH[g, 8], finalSH[b, 8]);
-        coeffs[3] = new Vector3(finalSH[r, 1], finalSH[g, 1], finalSH[b, 1]);
-        coeffs[4] = new Vector3(finalSH[r, 4], finalSH[g, 4], finalSH[b, 4]);
-        coeffs[5] = new Vector3(finalSH[r, 5], finalSH[g, 5], finalSH[b, 5]);
-        // Avoid reduplicating c_8.
-        coeffs[6] = new Vector3(finalSH[r, 2], finalSH[g, 2], finalSH[b, 2]);
-        coeffs[7] = new Vector3(finalSH[r, 7], finalSH[g, 7], finalSH[b, 7]);
-        coeffs[8] = new Vector3(finalSH[r, 6], finalSH[g, 6], finalSH[b, 6]) * 3.0f;
+        // Quadratic (4/5)
+        for (int c = 0; c < 3; c++)
+        {
+            coeffs[3 + c].x = finalSH[c, 4];
+            coeffs[3 + c].y = finalSH[c, 5];
+            coeffs[3 + c].z = finalSH[c, 6] * 3.0f;
+            coeffs[3 + c].w = finalSH[c, 7];
+        }
+
+        // Quadratic (5)
+        coeffs[6].x = finalSH[0, 8];
+        coeffs[6].y = finalSH[1, 8];
+        coeffs[6].z = finalSH[2, 8];
+        coeffs[6].w = 1.0f;
 
         cmd.SetGlobalVectorArray(HDShaderIDs._AmbientProbeCoeffs, coeffs);
     }

ScriptableRenderPipeline/HDRenderPipeline/HDRP/ShaderVariables.hlsl

     float2 _TaaFrameRotation; // {x = sin(_TaaFrameIndex * PI/2), y = cos(_TaaFrameIndex * PI/2), z = unused}
     uint   _TaaFrameIndex;    // [0, 7]
     // Volumetric lighting.
-    float4 _AmbientProbeCoeffs[9];      // 3 bands of SH, pre-convolved with the phase function. TODO: alpha is unused, pack better?
+    float4 _AmbientProbeCoeffs[7];      // 3 bands of SH, packed, rescaled and convolved with the phase function
     float  _GlobalFog_Asymmetry;
     float3 _GlobalFog_Scattering;
     float  _GlobalFog_Extinction;