barycenter

Barycentric correction of X-ray event files, in pure Python.

Tests Code of Conduct License: MIT

A photon’s arrival time at the spacecraft is not a useful clock. The spacecraft moves around the Earth, and the Earth around the Sun, so recorded arrival times wander by up to about 500 seconds over a year. Barycentring removes that: it reports each photon as it would have arrived at the solar system’s centre of mass, which is very nearly inertial, so times can be compared between missions and across decades.

That is the job of HEASOFT barycorr, XMM-SAS barycen and CIAO axbary. barycenter does it with nothing installed but Python, agrees with each of those tools to better than 100 nanoseconds, and is not restricted to the ephemerides they happen to ship.

Install

pip install git+https://github.com/matteobachetti/barycenter.git

The optional speed extra installs numba, which makes the spacecraft-clock interpolation about six times faster. It changes no answer — the compiled and pure-numpy kernels are bit-identical — and only pays for its compilation time on very large files:

pip install "barycenter[speed] @ git+https://github.com/matteobachetti/barycenter.git"

Use it from the command line

You need the event file and the mission’s orbit file, and a position on the sky:

barycenter nu30702012003A06_cl.evt nu30702012003A.attorb \
    --ra 148.9584 --dec 69.6794 --ephem DE440

The position can come from a pulsar parameter file instead, which saves copying coordinates by hand and also picks up the ephemeris the model was fitted with:

barycenter ni1013010101_0mpu7_cl.evt ni1013010101.orb --parfile Crab.par

The default engine takes the static RAJ/DECJ out of the model. If the source’s proper motion or parallax matters, add --engine pint, which hands the whole timing model over and so applies them.

Several orbit files can be given at once, for an observation split across them, and --clockfile names a spacecraft clock file (NuSTAR’s is fetched from the CALDB automatically if you do not). barycenter --help lists the rest; the options worth knowing about are --ephem, --engine, --dt and --fill-orbit-gaps (fit an orbit across long gaps in the orbit file, such as the Fermi LAT’s South Atlantic Anomaly).

If the good time intervals come in a separate file, as for SVOM, merge the ones you want into the event file first with barycenter-apply-gti (run it without -e to list them):

barycenter-apply-gti events.fits gtis.fits -e GTICAL-STA,GTICAL-NSA -o gti_events.fits

Use it from Python

The whole command line is one function:

from barycenter import apply_barycenter_correction

outfile = apply_barycenter_correction(
    "nu30702012003A06_cl.evt",
    "nu30702012003A.attorb",
    outfile="bary.evt",
    ra=148.9584,
    dec=69.6794,
    ephem="DE440",
)

If you only want the correction itself — to apply to times you already hold, or to plot — ask for the function of mission elapsed time and call it:

from barycenter import get_barycentric_correction, read_orbit

orbit = "nu30702012003A.attorb"
correction = get_barycentric_correction(orbit, ra=148.9584, dec=69.6794, ephem="DE440")

met = read_orbit(orbit)["MET"].value  # the times the orbit file covers
barycentred = met + correction(met)  # correction is in seconds, about +309 here

The returned object is a cubic spline over the orbit file’s own time span, so ask it only for times inside that span: outside it the spline extrapolates, quietly and badly. The command line checks this for you and warns.

Missions

Mission Orbit file Clock correction Agreement with the mission’s own tool
NuSTAR nu<obsid>A.attorb nuCclock*.fits, fetched automatically +22 ns vs barycorr
RXTE orbit/FPorbit_* tdc.dat, bundled +25 ns vs barycorr
Swift auxil/sw<obsid>sao.fits swclockcor*.fits (the UTCF) +37 ns vs barycorr
XMM-Newton PPS *ORBTSR*.FTZ none needed −42 ns vs SAS barycen
Chandra primary/orbitf*_eph1.fits none needed +47 ns vs CIAO axbary
NICER ni<obsid>.orb none needed +53 ns vs barycorr, one-off check
Fermi LAT FT2 spacecraft file none needed 30 ns vs gtbary
Fermi GBM glg_poshist_all_*.fit none needed gtbary refuses GBM files; identical to the LAT route given the same positions
SVOM SVOM_SVO-ORB-CNV_*.fits, one per pass; GTIs in a separate file none needed no official tool; 1 ulp vs barycorr on SVOM-format data relabelled as NICER
XRISM auxil/xa<obsid>.orb none needed +41 ns vs HEASOFT barycen
IXPE see the docs none implemented orbit file read, never checked against a tool

The first five rows, Fermi LAT and XRISM have a reference file committed to the repository, so continuous integration re-checks them on every change. The NICER figure comes from a single comparison against a barycorr run on a Crab observation too large to commit. GBM and SVOM have no official tool at all, so they are checked against the official tool of a mission with the same kind of orbit file, given the same numbers; see docs/missions.md. IXPE has never been compared against anything — its orbit file is read correctly, which is a different claim.

XMM-Newton and Chandra are worth singling out: barycorr cannot do either mission, and axbary can only reach DE200 and DE405, so this is the only way to barycentre them against a modern ephemeris. On the Chandra test data DE405 is 79 µs away from DE440.

Per-mission instructions are in docs/missions.md, and adding a new mission is about five lines plus a reference file to check it against.

Accuracy

The target is 100 nanoseconds against each mission’s own tool, and it is met on all five missions for which an official reference could be produced. Most of the residual scatter in the table above is not ours: a reference file stores times as 64-bit floats, which near a mission elapsed time of 5e8 seconds cannot express a difference finer than about 119 ns, so the tests add that floor to their tolerance.

Two engines compute the correction. The default, --engine native, uses astropy, ERFA and a JPL ephemeris directly, keeping every term as a small quantity so that it works identically on any platform. --engine pint routes the same job through PINT as an independent second opinion; the two agree to sub-nanosecond once the differing Shapiro-delay convention is accounted for. On Apple Silicon and Windows the PINT path is quantised at about 1.1 µs, because those platforms have no extended-precision float; the native engine is unaffected.

Documentation

Contributing

Bug reports and patches are welcome; see CONTRIBUTING.md and the Code of Conduct.

This site is open source. Improve this page.