How to flatten 3d astronomical coordinates into 2d polar form in Python

Question

I am making a 2d simulation of the solar system. For that you need a starting point for all planets. Right now I am storing that data in this form (ignore the color, radius, mass and Max_Trail_Length columns as they are cosmetics for the model):

Name	Color (RGB)	Radius (km)	Mass (kg)	Distance (AU)	Velocity (m/s)	Direction (radians)	Max_Trail_Length
Sun	(255, 255, 0)	6963400	1.98847e+30	0	0	0	100
Mercury	(255, 0, 0)	2440	3.3011e+23	0.393	47362	pi / 2	200
Venus	(255, 165, 0)	6052	4.8675e+24	0.723	35021.4	pi / 2	600
Earth	(0, 0, 255)	6371	5.9722e+24	1	29784.8	pi / 2	800
Mars	(255, 0, 0)	3390	6.4171e+23	1.52	24130.8	pi / 2	2000
Jupiter	(255, 222, 173)	69911	1.8982e+27	5.203	13070	pi / 2	10000
Saturn	(210, 180, 140)	58232	5.6834e+26	9.737	9690	pi / 2	20000
Uranus	(0, 0, 128)	25362	8.681e+25	19.61	6810	pi / 2	50000
Neptune	(0, 0, 255)	24622	1.02413e+26	29.897	5430	pi / 2	200000

As you can see, the direction is the same, so all planets are aligned, and the distance and velocity data is just the average.

For the record, I do not plan to make this into a 3d model, even though it would mitigate this problem entirely.

I want to collect real data from a date (January 1. 1750) using astropy, however when i collect the position and velocity using pos, vel = coordinates.get_body_barycentric_posvel(name, time) (coordinates from astropy), I get 3d data, as you can expect. But I want to convert it into data similar as I at the moment have.

I can't seem to do it right, as flattening the z axis doesn't work (I though it might be ignored at first as the solar system is relatively flat).

Here is the code I have right now for collecting this data.

import csv
import os

import numpy as np
from astropy import coordinates as coord
from astropy import units as u
from astropy.time import Time


# FIXME: Temporary wrong solution
def transform_data(pos, vel):
    pos = coord.CartesianRepresentation(pos.x, pos.y, 0 * u.au)
    vel = coord.CartesianDifferential(vel.x, vel.y, 0 * u.au / u.s)
    vel = vel.norm().to(u.m / u.s).value
    dist = pos.norm().to(u.AU).value

    return dist, vel


def get_path():
    parent_dir = os.path.dirname(os.path.dirname(os.path.abspath(__file__)))
    data_dir = os.path.join(parent_dir, 'data')
    os.makedirs(data_dir, exist_ok=True)
    file_path = os.path.join(data_dir, 'solsystem_data_1750.csv')
    return file_path


def process_data(time, writer):
    for name in ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Saturn', 'Uranus', 'Neptune']:
        pos, vel = coord.get_body_barycentric_posvel(name, time)

        dist, vel = transform_data(pos, vel)

        angle = np.arctan2(pos.y, pos.x).value

        writer.writerow([name, f'{dist:.5f}', f'{vel:.5f}', f'{angle:.5f}'])


def main():
    coord.solar_system_ephemeris.set('jpl')

    time = Time('1970-01-01')  # TODO: Change to 1750-01-01, but 1750 is a dubious year

    file_path = get_path()

    with open(file_path, 'w', newline='') as file:
        writer = csv.writer(file)
        writer.writerow(['Name', 'Distance (AU)', 'Velocity (m/s)', 'Direction (radians)'])

        process_data(time, writer)


if __name__ == '__main__':
    main()

I want to know if I can either modify my transform_data function, or collect the data in some sort of other way.

I couldn't find any 2d data from astropy, and I haven't found a way to convert the coordinates in a good way yet. I researched some alternatives to astropy, like just collecting the data from NASA's Horizons System, however I there run into the same problem of flattening the data.