Image Source.

Some things don’t change. Even back in 1846, the planet detection business in our solar system was a rough-and-tumble game. No sooner had Urbain Jean Joseph LeVerrier announced his prediction of the existence and position of Neptune, and had it dramatically verified by Galle and d’Arrest, than the British tried to jump all over the discovery and claim priority for Adams! Not to mention that tricky issue of names. LeVerrier tried various jostling maneuvers with the French Academy to try to get his planet named after himself, but his machinations were unsuccesful and Neptune stuck.

At least LeVerrier didn’t have to wrangle with Nineteenth Century player haters rushing to either (1) strip his newfound world of its planet status, or (2) consign it to the marginalia of the solar system. Neptune packs 8,065.34 times the mass of Pluto and 68,286.7 times the mass of Charon. It’s a planet with a capital P.

Here’s what I find amazing. The dramatic tension of the Prague IAU meeting apparently hinges on the future nomenclature for 2003 UB-313 and its ilk. Even the New York Times, the United States paper of record, makes the episode seem like one of the scientific Big Deals of the year. Oklo dot org manages to climb out of its summer visitors slump on the basis of two chatty posts on the Pluto debate. While all this is going on, an amazing multiple-planet system orbiting the nearby solar-type star HD 160691 receives a new planet and a dramatically improved characterization, and hardly anyone notices.

The downloadable systemic console and the systemic back-end contain several data-sets for HD 169061. The exact data set used by the Geneva group in their astro-ph paper from yesterday is listed in the system menu as HD160691_M04P06CH. (No sooner had I laboriously typed in the table from the .pdf file than Eugenio pointed out that one can simply copy-paste from the text file source on astro-ph. Doh!)

The data used by the Geneva group comes from three telescopes. It includes AAT data from McCarthy et al. 2004, as well as older data from Coralie and new, extremely high-precision data from HARPS. The console therefore loads with three offset sliders.

The periodogram of the data shows that the power spectrum peaks at P=645 days. Activate a planet on the console, and give it this period. Then optimize on mean anomaly and mass. Then send the three offset sliders, the mass, and the period out for a Levenberg-Marquardt minimization or “polish”. This gives the best circular 1-planet fit to the data set:

A periodogram of the residuals shows a high peak at 2971 days. Activate a second planet with this period and employ a similar minimization procedure to get the best 2-planet fit using circular orbits:

Now, the highest residual peak in the periodogram is near 300 days. This corresponds to the newly announced planet:

A three-planet fit consisting of circular orbits does a reasonable job of running through the data:

With three circular orbits fitted out, the residuals periodogram shows a spike at a period just less than 10 days. This is the Neptune-mass planet announced in 2004 by the Swiss team, using the early high-precision results from HARPS:

The planet produces a high-frequency modulation on the overall radial velocity curve. In this snapshot, the orbits are all still circular:

Introducing eccentricities and longitudes of periastron as free parameters allows for a quick improvement of the fit. One quickly recovers the fit published in Table 1 of Pepe et al.’s paper (note that the console reports the chi-square statistic, whereas Pepe et al. report the square root of the chi-square statistic, as is often the custom in reporting planetary fits. Here’s the resulting orbital plot:

The console screenshot looks like this:

Whenever you get a fit that you like, be sure to save it:

And upload it to the systemic back-end:

Stefano has made sure that the back-end keeps things scrupulously honest. When I uploaded my 4-planet fit, the back-end reported a miserable chi-square statistic:

The reason for this degradation of the fit becomes clear when you click the integrate button on the console. The planets in this fit exert gravitational perturbations on each other which are strong enough to screw up the fit over the timescale that the system has been observed. The Keplerian radial velocity curve and the N-body radial velocity curve are deviating significantly near the end of the radial velocity baseline:

In order to get an accurate characterization of the system, you need to use the console with the integrate button activated. Pepe et al. pointed out this issue, and quite a bit of their paper is devoted to exploring self consistent fits. We’ll take these up in an upcoming post. The best way to get a self-consistent fit is to reproduce the best Keplerian fit, and then optimize with integration turned on:

This ups the computation time dramatically. The fan really starts to roar on my dual-G5 as soon as I hit the polish button…