http://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf
        x_ecl = cos(w)cos(Om)-sin(w)sin(Om)cos(I) * x_prime +
            (-sin(w)cos(Om) - cos(w)sin(Om)cos(I)) * y_prime
        y_ecl = cos(w)sin(Om)-sin(w)cos(Om)cos(I) * x_prime +
            (-sin(w)cos(Om) - cos(w)sin(Om)cos(I)) * y_prime
        z_ecl = (sin(w)sin(I)) * x_prime +
            (cos(w)sin(I)) * y_prime

def getECLCoordinatesAtTime(self, date):
         # localize the orbital parameters
        a = self.semimajor_axis
        e = self.eccentricity
        I = self.inclination
        Om = self.longitude_ascending
        #n = 0.230605479
        n = 0.230652907

        w = self.argument_perihelion

        M = self.mean_anomaly
        d = date.tdb - self.epoch.tdb

        Om = Om / 360.0 * constants.tau
        w = w / 360.0 * constants.tau
        I = I / 360.0 * constants.tau
        M = M / 360.0 * constants.tau
        n = n / 360.0 * constants.tau

        M += d * n

        # calculate the mean anomaly in rads
        E = convergeEccentricAnomaly(M, e, 30)

        # calculate the initial primes
        x_prime = a * (cos(E) - e)
        y_prime = a * (1 - e ** 2.0) ** (0.5) * sin(E)

        """
        http://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf
        x_ecl = cos(w)cos(Om)-sin(w)sin(Om)cos(I) * x_prime +
            (-sin(w)cos(Om) - cos(w)sin(Om)cos(I)) * y_prime
        y_ecl = cos(w)sin(Om)-sin(w)cos(Om)cos(I) * x_prime +
            (-sin(w)cos(Om) - cos(w)sin(Om)cos(I)) * y_prime
        z_ecl = (sin(w)sin(I)) * x_prime +
            (cos(w)sin(I)) * y_prime
        """

        # calculate the ecliptic coordinates
        x_ecl = ((cos(w) * cos(Om) - sin(w) * sin(Om) * cos(I)) * x_prime +
                (-1 * sin(w) * cos(Om) - cos(w) * sin(Om) * cos(I)) * y_prime)
        y_ecl = ((cos(w) * sin(Om) + sin(w) * cos(Om) * cos(I)) * x_prime +
                (-1 * sin(w) * sin(Om) + cos(w) * cos(Om) * cos(I)) * y_prime)
        z_ecl = ((sin(w) * sin(I)) * x_prime + (cos(w) * sin(I)) * y_prime)

        return ICRCoordinates(x_ecl, y_ecl, z_ecl)