[Workaround included]How to get rid of the view angle dependent gaussian artifacts generated by mcmc?

[WindyYam]

I'm currently doing nerfstudio with mcmc strategy on a simple scene. For what I've learned, the mcmc strategy tend to generate very large and spiky gaussians.

I'm ok with those, but I don't understand why they create view dependent artifact?

For example, I can see the very large gaussian when I zoom out the whole scene, and they are actually behind the room wall structure.

But if I rotate the camera in some angles, those gaussians defy their scene order, invade into the scene and completely block the view.

The other one is those spiky gaussians, if view from far away they all in correct position, but if they are viewed close up, partly occluded outside the camera view, they will point to the angle that totally subject to my camera look angle

If I export the scene and view it in supersplat it looks better, but they are the same as in the nerfstudio viewer when viewing from Unity gaussian plugin.

My guess is supersplat renders differently for those gaussians, maintaining their correct view properties.

Is there a way to overcome this? I really like the mcmc strategy as it produces overall better scene, but I can't use it on Unity due to the anonying artifacts which is a pity.

For now I can only turn on the use_scale_regularization, and manually add a max size regulator. It will look better, but not completely rule out the issue.

[WindyYam]

I'm not familiar with gsplat code, but the same artifact issue happen in also this Unity gaussian rendering repoaras-p/UnityGaussianSplatting#196, and I think I start to understand a little bit on that one, don't know if it is related to the one in gsplat
In Unity I aligned a real box geometry to a splat, they have same 3d center, so the splat cover the whole box geometry

But when close up and oblique the camera, on screen splat start to be like this. The position of the axes is the world space center

Apparently the shader try to render the splat on the on screen projection of the splat center, but on close up and oblique, the center in world space is not necessarily the center of where the geometry covers in screen space, and in some extreme angle this happens

[WindyYam]

After looking into gsplat cuda kernel I believe gsplat did the same projection of 3d mean into 2d mean and use it as on screen center.
I'm not a cuda guy at the moment, here is a workaround from my faithful chatGPT5:

inline __device__ void persp_proj(
    // inputs
    const vec3<T> mean3d,
    const mat3<T> cov3d,
    const T fx,
    const T fy,
    const T cx,
    const T cy,
    const uint32_t width,
    const uint32_t height,
    // outputs
    mat2<T> &cov2d,
    vec2<T> &mean2d
) {
    T x = mean3d[0], y = mean3d[1], z = mean3d[2];

    T tan_fovx = 0.5f * width / fx;
    T tan_fovy = 0.5f * height / fy;
    T lim_x_pos = (width - cx) / fx + 0.3f * tan_fovx;
    T lim_x_neg = cx / fx + 0.3f * tan_fovx;
    T lim_y_pos = (height - cy) / fy + 0.3f * tan_fovy;
    T lim_y_neg = cy / fy + 0.3f * tan_fovy;

    T rz = 1.f / z;
    T rz2 = rz * rz;
    T tx = z * min(lim_x_pos, max(-lim_x_neg, x * rz));
    T ty = z * min(lim_y_pos, max(-lim_y_neg, y * rz));

    // mat3x2 is 3 columns x 2 rows.
    mat3x2<T> J = mat3x2<T>(
        fx * rz,
        0.f, // 1st column
        0.f,
        fy * rz, // 2nd column
        -fx * tx * rz2,
        -fy * ty * rz2 // 3rd column
    );
    cov2d = J * cov3d * glm::transpose(J);
    mean2d = vec2<T>({fx * x * rz + cx, fy * y * rz + cy});

    // Mean-shift correction (second-order term under perspective).
    // Uses camera-space covariance: Vxz = cov3d[0][2], Vyz = cov3d[1][2], Vzz = cov3d[2][2].
    const T Vxz = cov3d[0][2];
    const T Vyz = cov3d[1][2];
    const T Vzz = cov3d[2][2];
    const T rz3 = rz2 * rz;
    vec2<T> meanShift = vec2<T>(
        fx * (-Vxz * rz2 + tx * Vzz * rz3),
        fy * (-Vyz * rz2 + ty * Vzz * rz3)
    );
    mean2d += (T)GSPLAT_MEAN_SHIFT_COEFF * meanShift;
}

We might need to add the backward pass term as well? Please correct this if I'm wrong

inline __device__ void persp_proj_vjp(
    // fwd inputs
    const vec3<T> mean3d,
    const mat3<T> cov3d,
    const T fx,
    const T fy,
    const T cx,
    const T cy,
    const uint32_t width,
    const uint32_t height,
    // grad outputs
    const mat2<T> v_cov2d,
    const vec2<T> v_mean2d,
    // grad inputs
    vec3<T> &v_mean3d,
    mat3<T> &v_cov3d
) {
    T x = mean3d[0], y = mean3d[1], z = mean3d[2];

    T tan_fovx = 0.5f * width / fx;
    T tan_fovy = 0.5f * height / fy;
    T lim_x_pos = (width - cx) / fx + 0.3f * tan_fovx;
    T lim_x_neg = cx / fx + 0.3f * tan_fovx;
    T lim_y_pos = (height - cy) / fy + 0.3f * tan_fovy;
    T lim_y_neg = cy / fy + 0.3f * tan_fovy;

    T rz = 1.f / z;
    T rz2 = rz * rz;
    T tx = z * min(lim_x_pos, max(-lim_x_neg, x * rz));
    T ty = z * min(lim_y_pos, max(-lim_y_neg, y * rz));

    // mat3x2 is 3 columns x 2 rows.
    mat3x2<T> J = mat3x2<T>(
        fx * rz,
        0.f, // 1st column
        0.f,
        fy * rz, // 2nd column
        -fx * tx * rz2,
        -fy * ty * rz2 // 3rd column
    );

    // cov = J * V * Jt; G = df/dcov = v_cov
    // -> df/dV = Jt * G * J
    // -> df/dJ = G * J * Vt + Gt * J * V
    v_cov3d += glm::transpose(J) * v_cov2d * J;

    // df/dx = fx * rz * df/dpixx
    // df/dy = fy * rz * df/dpixy
    // df/dz = - fx * mean.x * rz2 * df/dpixx - fy * mean.y * rz2 * df/dpixy
    v_mean3d += vec3<T>(
        fx * rz * v_mean2d[0],
        fy * rz * v_mean2d[1],
        -(fx * x * v_mean2d[0] + fy * y * v_mean2d[1]) * rz2
    );

    // df/dx = -fx * rz2 * df/dJ_02
    // df/dy = -fy * rz2 * df/dJ_12
    // df/dz = -fx * rz2 * df/dJ_00 - fy * rz2 * df/dJ_11
    //         + 2 * fx * tx * rz3 * df/dJ_02 + 2 * fy * ty * rz3
    T rz3 = rz2 * rz;
    mat3x2<T> v_J = v_cov2d * J * glm::transpose(cov3d) +
                    glm::transpose(v_cov2d) * J * cov3d;

    // fov clipping
    if (x * rz <= lim_x_pos && x * rz >= -lim_x_neg) {
        v_mean3d.x += -fx * rz2 * v_J[2][0];
    } else {
        v_mean3d.z += -fx * rz3 * v_J[2][0] * tx;
    }
    if (y * rz <= lim_y_pos && y * rz >= -lim_y_neg) {
        v_mean3d.y += -fy * rz2 * v_J[2][1];
    } else {
        v_mean3d.z += -fy * rz3 * v_J[2][1] * ty;
    }
    v_mean3d.z += -fx * rz2 * v_J[0][0] - fy * rz2 * v_J[1][1] +
                  2.f * fx * tx * rz3 * v_J[2][0] +
                  2.f * fy * ty * rz3 * v_J[2][1];
    
    // Mean-shift correction gradients
    const T Vxz = cov3d[0][2];
    const T Vyz = cov3d[1][2]; 
    const T Vzz = cov3d[2][2];
    const T rz4 = rz3 * rz;

    // Gradients w.r.t. mean3d from mean-shift term
    const T coeff = (T)GSPLAT_MEAN_SHIFT_COEFF;
    v_mean3d.x += coeff * fx * Vzz * rz3 * v_mean2d[0];
    v_mean3d.y += coeff * fy * Vzz * rz3 * v_mean2d[1];
    v_mean3d.z += coeff * (
        fx * (2.0f * Vxz * rz3 - 3.0f * tx * Vzz * rz4) * v_mean2d[0] +
        fy * (2.0f * Vyz * rz3 - 3.0f * ty * Vzz * rz4) * v_mean2d[1]
    );

    // Gradients w.r.t. cov3d from mean-shift term
    v_cov3d[0][2] += coeff * (-fx * rz2 * v_mean2d[0]);
    v_cov3d[1][2] += coeff * (-fy * rz2 * v_mean2d[1]);
    v_cov3d[2][0] += coeff * (-fx * rz2 * v_mean2d[0]); // symmetric
    v_cov3d[2][1] += coeff * (-fy * rz2 * v_mean2d[1]); // symmetric
    v_cov3d[2][2] += coeff * (
        fx * tx * rz3 * v_mean2d[0] + 
        fy * ty * rz3 * v_mean2d[1]
    );
}

This basically add an offset on screen to better represent the on screen mean2d. GSPLAT_MEAN_SHIFT_COEFF is offset scale factor. Whether or not it is mathmatical accurate? I don't believe, and I'm not a mathmatician myself. It does not completely rule out those artifacts, but it does help a bit on reducing the artifacts in mcmc scene
Before:
https://github.com/user-attachments/assets/90aea4bd-864a-4f81-8fbd-ebc92d1518f5
After (GSPLAT_MEAN_SHIFT_COEFF 2.0):
https://github.com/user-attachments/assets/2592f460-76d5-4c60-bf65-192bd6b6439b
GSPLAT_MEAN_SHIFT_COEFF 3.0:
https://github.com/user-attachments/assets/8c34927a-07ac-4b86-acaa-56006fd9532a