The code below is how I calculate KRt from P:
def rq_decomposition(M):
Q, R = np.linalg.qr(np.flipud(M).T)
R = np.flipud(R.T)
Q = Q.T
R = np.fliplr(R)
Q = np.fliplr(Q)
return R, Q

def P_to_KRt(P):

    M = P[0:3, 0:3]
    R, Q = rq_decomposition(M)

    # Normalize K such that K[2,2] = 1
    K = R / R[-1, -1]

    # Ensure positive diagonal elements of K
    T = np.diag(np.sign(np.diag(K)))
    K = np.dot(K, T)
    Q = np.dot(T, Q)

    if np.linalg.det(Q) < 0:
        Q = -Q

    # Compute the translation vector t
    t = np.dot(np.linalg.inv(K), P[:, 3])

    return K, Q, t

[K1_est, R1_est, t1_est] = P_to_KRt(P1_est)

Then to use P = K[R|t],to verify if the P = P1_est or not:

def recompose_projection_matrix(K, R, t):
    # Ensure t is a column vector
    t = t.reshape(-1, 1) if t.ndim == 1 else t

    # Concatenate R and t to form [R | t]
    Rt = np.hstack((R, t))

    P_recomposed = np.dot(K, Rt)

    return P_recomposed

Clearly the K,R and t above is from P1_est decoposition, but when I recompose P_recompose with the same K, R and t, P_recompose is not equal to P1_est, I cannot find where is wrong