<p>Probably the most reliable guide to <code>vsop2000-p11.dat</code> is the reference implementation <code>vsop2000.for</code>.
The input line in question is the first of 5399 k=0 records contributing to <span class="math-container">$x$</span>.
Those are followed by 5842 k=1 records, 610 k=2 records, and so on up to k=8.
Then <span class="math-container">$y$</span> (k=0..8) and <span class="math-container">$z$</span> (k=0..7) follow the same pattern.</p>

<p>If I read the Fortran code correctly, the fixed-length (not blank-delimited) record format is:</p>

<ul>
<li>Columns 1-7 are an index <span class="math-container">$i$</span></li>
<li>Columns 9-63 are 17 per-body coefficients <span class="math-container">$b_{kij}$</span></li>
<li>Columns 64-72 are three zeroes, ignored</li>
<li>Columns 73-126 are coefficients <span class="math-container">$s_{ki}, c_{ki}$</span> in astronomical units</li>
</ul>

<p>To compute a position for time <span class="math-container">$t$</span> in Julian years since J2000:</p>

<ol>
<li>For each body <span class="math-container">$j$</span>, compute <span class="math-container">$\theta_j(t) = \dot{\theta}_j t + \theta_j(0)$</span>, where <span class="math-container">$\dot{\theta}_j$</span> is the mean motion and <span class="math-container">$\theta_j(0)$</span> is the J2000 mean longitude from <code>vsop2000-exp.txt</code></li>
<li>For each record <span class="math-container">$x_{ki}$</span>, compute <span class="math-container">$\phi_{ki} = \sum\limits_{j=1}^{17} b_{kij} \theta_j(t)$</span></li>
<li><span class="math-container">$x = \sum\limits_{k=0}^8 t^k \sum\limits_{i=1}^{n_k} (s_{ki} \sin \phi_{ki} + c_{ki} \cos \phi_{ki})$</span></li>
<li>For <span class="math-container">$y$</span> and <span class="math-container">$z$</span>, repeat steps 2 and 3 using the respective sets of records.</li>
</ol>

<p>This should give you a geocentric Moon position in astronomical units.</p>

CC BY-SA 4.0
2020-04-07T19:22:16.167Z
A note to any hypothetical future readers: this works like a charm.
John Dumancic
CC BY-SA 4.0
2020-04-07T19:31:53.127Z
@Tesseract I modified vsop2000.for to give Moon positions for 2020-2023 in spherical coordinates accounting for precession. lat/lon differences from JPL HORIZONS (DE430) were on the order of 15". Moisson and Bretagnon 2001 compare their results to DE403.
Mike G
CC BY-SA 4.0
2020-11-21T20:35:17.260Z
In the sample line given by John Dumancic in his original question, we find e.g. “-.1053243070325313D-09.” What does the D stand for? I’d be expecting an E here, for E-09 as in one billionth, but a D…?
Pierre Paquette
CC BY-SA 4.0
2020-11-22T01:55:58.363Z
@PierrePaquette I think that's just double precision floating point, no other difference from E.
Mike G
CC BY-SA 4.0
2020-11-27T05:00:28.690Z
@PierrePaquette $n_k$ is the number of records in the `.dat` file for a given value of $k$.
Mike G
CC BY-SA 4.0
2020-04-07T01:09:48.830Z
2020-04-07T01:09:48.830Z
Mike G
<p>For a project I'm doing, I need to implement the VSOP-2000 ephemeris for the Moon. However, unlike VSOP-87, there is little to no documentation online on how to use the data. For example, here is the first line of the x-term for the Moon in VSOP-2000:</p>

<pre><code>  1   0  0  1  0   0  0  0  0   0  1  0  0   0  0  0  0   0  0  0  0     -.1053243070325313D-09     -.2560241101399878D-02
</code></pre>

<p>This is very different from VSOP-87, which I know how to use. Are there any resources online that could inform me on how to use this?</p>

CC BY-SA 4.0
2020-04-06T15:51:34.700Z
Are you sure you need to use the VSOP-2000 ephemeris and you can not use another like one of the [JPL development ephemerides](https://en.wikipedia.org/wiki/Jet_Propulsion_Laboratory_Development_Ephemeris), which are widely used and well supported and well documented? You *really must use* an ephemeris that has "little to no documentation online"?
uhoh
CC BY-SA 4.0
2020-04-06T17:55:25.217Z
@uhoh as long as I can implement it in Mathematica, and it gets the Moon's position to >10'' for a couple millennia forward and back, I don't really care what it is. I just heard VSOP 2000 could do this.
John Dumancic
CC BY-SA 4.0
2020-04-06T20:45:01.750Z
Okay I see. So SPICE (mentioned in [this comment](https://astronomy.stackexchange.com/questions/35669/how-accurate-should-sun-moon-ephemeris-be-for-calculating-eclipses#comment70396_35669)) uses the JPL development ephemerides (DE's) as does the Python package [Skyfield](https://rhodesmill.org/skyfield/). The DEs are tables of coefficients that are used with a special form of Chebyshev polynomial interpolation. See [this answer](https://space.stackexchange.com/a/35806/12102) I am sure Mathematica can do this but it may be quite a lot of work to do what is designed to do in code.
uhoh
CC BY-SA 4.0
2020-04-06T20:49:13.657Z
If you know any Python at all then Skyfield is by far the easiest way to do implement an eclipse finding algorithm. I don't know if it will be for sure accurate to 10 arcsec over millennia but the documentation for each of the various DEs will go into detail on the expected accuracy. They use all known data including the centimeter-level laser ranging of the Moon that's been going on regularly since the Apollo missions put corner cube reflectors on the Moon in the 1970's.
uhoh
CC BY-SA 4.0
2020-04-06T20:51:23.013Z
more about Chebyshev interpolation of DE's also in [this answer](https://space.stackexchange.com/a/37882/12102), possibly also helpful are answers to [Is there a way to extract the Chebyshev coefficients for a body from a SPICE kernel?](https://space.stackexchange.com/q/10490/12102)
uhoh
<the-moon><ephemeris>
Implementing the VSOP-2000 Ephemeris
John Dumancic
# Implementing the VSOP-2000 Ephemeris

(Asked by: username_0 on 2020-04-06)

For a project I'm doing, I need to implement the VSOP-2000 ephemeris for the Moon. However, unlike VSOP-87, there is little to no documentation online on how to use the data. For example, here is the first line of the x-term for the Moon in VSOP-2000:
  1   0  0  1  0   0  0  0  0   0  1  0  0   0  0  0  0   0  0  0  0     -.1053243070325313D-09     -.2560241101399878D-02

This is very different from VSOP-87, which I know how to use. Are there any resources online that could inform me on how to use this?

**Tags:** the-moon, ephemeris
> **Comments:**
> **username_1**: Are you sure you need to use the VSOP-2000 ephemeris and you can not use another like one of the [JPL development ephemerides](https://en.wikipedia.org/wiki/Jet_Propulsion_Laboratory_Development_Ephemeris), which are widely used and well supported and well documented? You *really must use* an ephemeris that has "little to no documentation online"?  
> **username_0**: @uhoh as long as I can implement it in Mathematica, and it gets the Moon's position to >10'' for a couple millennia forward and back, I don't really care what it is. I just heard VSOP 2000 could do this.  
> **username_1**: Okay I see. So SPICE (mentioned in [this comment](https://astronomy.stackexchange.com/questions/35669/how-accurate-should-sun-moon-ephemeris-be-for-calculating-eclipses#comment70396_35669)) uses the JPL development ephemerides (DE's) as does the Python package [Skyfield](https://rhodesmill.org/skyfield/). The DEs are tables of coefficients that are used with a special form of Chebyshev polynomial interpolation. See [this answer](https://space.stackexchange.com/a/35806/12102) I am sure Mathematica can do this but it may be quite a lot of work to do what is designed to do in code.  
> **username_1**: If you know any Python at all then Skyfield is by far the easiest way to do implement an eclipse finding algorithm. I don't know if it will be for sure accurate to 10 arcsec over millennia but the documentation for each of the various DEs will go into detail on the expected accuracy. They use all known data including the centimeter-level laser ranging of the Moon that's been going on regularly since the Apollo missions put corner cube reflectors on the Moon in the 1970's.  
> **username_1**: more about Chebyshev interpolation of DE's also in [this answer](https://space.stackexchange.com/a/37882/12102), possibly also helpful are answers to [Is there a way to extract the Chebyshev coefficients for a body from a SPICE kernel?](https://space.stackexchange.com/q/10490/12102)  

---

## 1 answers

### Answer by username_2 on 2020-04-07. Score: 5 | Accepted answer

Probably the most reliable guide to `vsop2000-p11.dat` is the reference implementation `vsop2000.for`.
The input line in question is the first of 5399 k=0 records contributing to $x$.
Those are followed by 5842 k=1 records, 610 k=2 records, and so on up to k=8.
Then $y$ (k=0..8) and $z$ (k=0..7) follow the same pattern.

If I read the Fortran code correctly, the fixed-length (not blank-delimited) record format is:

* Columns 1-7 are an index $i$
* Columns 9-63 are 17 per-body coefficients $b\_{kij}$
* Columns 64-72 are three zeroes, ignored
* Columns 73-126 are coefficients $s\_{ki}, c\_{ki}$ in astronomical units

To compute a position for time $t$ in Julian years since J2000:

1. For each body $j$, compute $\theta\_j(t) = \dot{\theta}\_j t + \theta\_j(0)$, where $\dot{\theta}\_j$ is the mean motion and $\theta\_j(0)$ is the J2000 mean longitude from `vsop2000-exp.txt`
2. For each record $x\_{ki}$, compute $\phi\_{ki} = \sum\limits\_{j=1}^{17} b\_{kij} \theta\_j(t)$
3. $x = \sum\limits\_{k=0}^8 t^k \sum\limits\_{i=1}^{n\_k} (s\_{ki} \sin \phi\_{ki} + c\_{ki} \cos \phi\_{ki})$
4. For $y$ and $z$, repeat steps 2 and 3 using the respective sets of records.

This should give you a geocentric Moon position in astronomical units.

> **Comments:**
> **username_0**: A note to any hypothetical future readers: this works like a charm.  
> **username_2**: @Tesseract I modified vsop2000.for to give Moon positions for 2020-2023 in spherical coordinates accounting for precession. lat/lon differences from JPL HORIZONS (DE430) were on the order of 15". Moisson and Bretagnon 2001 compare their results to DE403.  
> **username_3**: In the sample line given by John Dumancic in his original question, we find e.g. “-.1053243070325313D-09.” What does the D stand for? I’d be expecting an E here, for E-09 as in one billionth, but a D…?  
> **username_2**: @PierrePaquette I think that's just double precision floating point, no other difference from E.  
> **username_2**: @PierrePaquette $n_k$ is the number of records in the `.dat` file for a given value of $k$.  

---

