Calculating satellite trajectories through ephemeris is not elliptical

Cesium for Unity
1

I calculated the satellite’s orbit through broadcast ephemeris in U3D. The calculation result is consistent with the precision ephemeris. However, when drawing the orbit, the satellite’s orbit is not an ellipse. How to deal with it? Thanks！

double3 CalutePosition()
    {
        n0 = Math.Sqrt(GM) / Math.Pow(sqrt_a, 3);
        n = n0 + Delta_n;

        t = TimetoGPSsec(year, month, day, hour, minute, second);
        Delta_t = a0 + a1 * (t - toc) + a2 * (t - toc) * (t - toc);
        t = t - Delta_t;
        tk = t - toe;
        if (tk > 302400)
            tk -= 604800;
        else
        if (tk < -302400)
            tk += 604800;

        Mk = M0 + n * tk;

        double E0 = Mk;
        double Ek = Mk + e * Math.Sin(E0);
        while (Math.Abs(Ek - E0) > 10e-12)
        {
            E0 = Ek;
            Ek = Mk + e * Math.Sin(E0);
        }

        Vk = Math.Atan2(Math.Sqrt(1 - e * e) * Math.Sin(Ek), (Math.Cos(Ek) - e));
        fai_k = Vk + w;

        _6u = Cuc * Math.Cos(2 * fai_k) + Cus * Math.Sin(2 * fai_k);
        _6r = Crc * Math.Cos(2 * fai_k) + Crs * Math.Sin(2 * fai_k);
        _6i = Cic * Math.Cos(2 * fai_k) + Cis * Math.Sin(2 * fai_k);

        uk = fai_k + _6u;
        rk = Math.Pow(sqrt_a, 2) * (1 - e * Math.Cos(Ek)) + _6r;
        ik = i0 + _6i + IDOT * tk;

        Xk = rk * Math.Cos(uk);
        Yk = rk * Math.Sin(uk);

        OMEGA_k = OMEGA_0 + (OMEGA_DOT - we) * tk - we * toe;
        X_obs_k = Xk * Math.Cos(OMEGA_k) - Yk * Math.Cos(ik) * Math.Sin(OMEGA_k);
        Y_obs_k = Xk * Math.Sin(OMEGA_k) + Yk * Math.Cos(ik) * Math.Cos(OMEGA_k);
        Z_obs_k = Yk * Math.Sin(ik);
        return new double3(X_obs_k / 1000f, Y_obs_k / 1000f, Z_obs_k / 1000f);

    }

  globeAnchor.positionGlobeFixed = CalutePosition();
2

This is my orbit.
image