Get the positions and velocities of the earth and sun from date t0 to t1.
    Returned as a tuple q, v
    q: Nx3 array of positions (x, y, z) in the J2000.0 coordinate frame.

def configuration_ES(t0: date, t1: Optional[date] = None, 
                     steps_per_day: int = None) -> Tuple[np.ndarray, np.ndarray]:
    """
    Get the positions and velocities of the earth and sun from date t0 to t1.
    Returned as a tuple q, v
    q: Nx3 array of positions (x, y, z) in the J2000.0 coordinate frame.
    """
    # Default steps_per_day = 1
    if steps_per_day is None:
        steps_per_day = 1
        
    # Time step dt is 1.0 over steps per day
    dt: float = 1.0 / float(steps_per_day)

    # Default t1 to one day after t0
    if t1 is not None:
        # Convert t to a julian day
        jd0: int = julian_day(t0)
        jd1: int = julian_day(t1)
    else:
        jd0: int = julian_day(t0)
        jd1: int = jd0 + dt
    # Pass the times as an array of julian days
    jd: np.ndarray = np.arange(jd0, jd1, dt)

    # Position and velocity of the sun as arrays of length 3
    sun_id: int = jpl_body_id['sun']
    pos_sun, vel_sun = jpl_kernel[0, sun_id].compute_and_differentiate(jd)

    # Position and velocity of the earth as arrays of length 3
    earth_id: int = jpl_body_id['earth']
    pos_earth, vel_earth = jpl_kernel[0, earth_id].compute_and_differentiate(jd)

    # Convert positions from km to meters (multiply by km2m)
    q = np.vstack([pos_sun, pos_earth]).T * km2m
    # Convert velocities from km / day to meters / sec (multiply by km2m, divide by day2sec)
    v = np.vstack([vel_sun, vel_earth]).T * (km2m / day2sec)

    # Return tuple of Tx6 arrays for position q and velocity v
    return q, v