MNRAS 000 ,

1-9 (2019)

## Photon-weighted barycentric correction and its importance for precise radial velocities

Ren´ e Tronsgaard, 1 ? Lars A. Buchhave, 1 Jason T. Wright ; 2 Jason D. Eastman, 3 and Ryan T. Blackman 4

- 1 DTU Space, National Space Institute, Technical University of Denmark, Elektrovej 328, DK-2800 Kgs. Lyngby, Denmark 2 Center for Exoplanets and Habitable Worlds, Department of Astronomy &amp; Astrophysics, The Pennsylvania State University, 525 Davey Laboratory, University Park, PA 16802, USA
- 3 Center for Astrophysics | Harvard &amp; Smithsonian, 60 Garden St, Cambridge, MA 02138, USA
- 4 Department of Astronomy, Yale University, 52 Hillhouse Avenue, New Haven, CT 06511, USA

Accepted 2019 July 25. Received 2019 July 9; in original form 2019 July 9.

## ABSTRACT

When applying the barycentric correction to a precise radial velocity measurement, it is common practice to calculate its value only at the photon-weighted midpoint time of the observation instead of integrating over the entire exposure. However, since the barycentric correction does not change linearly with time, this leads to systematic errors in the derived radial velocities. The typical magnitude of this second-order effect is of order 10 cms GLYPH&lt;0&gt; 1 , but it depends on several parameters, e.g. the latitude of the observatory, the position of the target on the sky, and the exposure time. We show that there are realistic observing scenarios, where the errors can amount to more than 1 ms 1 GLYPH&lt;0&gt; . We therefore recommend that instruments operating in this regime always record and store the exposure meter flux curve (or a similar measure) to be used as photon-weights for the barycentric correction. In existing data, if the flux curve is no longer available, we argue that second-order errors in the barycentric correction can be mitigated by adding a correction term assuming constant flux.

Key words: instrumentation: spectrographs - techniques: radial velocities

## 1 INTRODUCTION

In order to search for and characterise potential Earthanalogues outside our own Solar System, astronomers will need precise radial velocity (PRV) instruments that can deliver a yet unseen precision of a few cms GLYPH&lt;0&gt; 1 over timescales of several years (e.g. Fischer et al. 2016). Among the numerous prerequisites for being able to make such measurements is the ability to accurately transform the measured Doppler shift into a stationary reference frame with respect to the barycentre of our Solar System.

Wright &amp; Eastman (2014) presented a detailed review of how to calculate and perform barycentric corrections to an accuracy better than 1 cm s GLYPH&lt;0&gt; 1 . The authors provided an IDL implementation ZBARYCORR , which has later been ported to Python (Kanodia &amp; Wright 2018). However, one important question that has not been fully addressed in literature is how to best apply the instantaneous barycentric correction to an exposure that is extended in time. During the integration time of an exposure, usually tens of minutes, there can

- ? E-mail: rtr@space.dtu.dk (RT)

be a significant change in the barycentric correction, mainly due to Earth's rotation, and the recorded spectrum will get slightly smeared out on the detector, effectively broadening the stellar absorption lines. Time-variable changes in the flux that reaches the spectrograph, caused by clouds, variable seeing, poor telescope guiding, or other issues, will make this smearing non-uniform, shift the weight of the broadened lines, and introduce spurious Doppler shifts. This type of error, which we will refer to in this paper as the firstorder effect, is usually corrected for by monitoring the photon flux through the spectrograph with an exposure meter. From these measurements, the photon-weighted midpoint time of the exposure can be determined, ensuring that the barycentric velocity is calculated at the mean time of arrival of the collected photons instead of the 'geometric' midpoint time of the exposure, equidistant from the shutter open and close times. However, this method neglects the fact that the barycentric correction does not change linearly. As we will show in the following, a more correct approach is to weight the corrected velocity with the exposure meter flux curve and in that way capture the higher-order variation of the barycentric motion of the observatory. Neglecting what we will refer to as the second-order effect can lead to typical

systematic errors of order 10 cms GLYPH&lt;0&gt; 1 for exposure times of 30-60 minutes. In some quite realistic observing scenarios it will cause errors greater than 1 ms GLYPH&lt;0&gt; 1 .

The problem was first pointed out in Fischer et al. (2016), referring to an early draft version of this manuscript, and it was briefly mentioned again in Blackman et al. (2017), along with a description of the chromatic dependence of the barycentric correction. This paper aims to provide a more detailed and formal description of how to compute the barycentric correction for an extended exposure, estimate the error that arises when the second-order effect is neglected, and discuss some recommendations for future observations and existing data.

In Section 2, we describe how to calculate the barycentric correction for an extended exposure with the best possible accuracy, and in Section 3 we estimate the errors that result from using the geometric and photon-weighted midpoint times. In Section 4, we simulate the error using the full barycentric correction algorithm to confirm our analytic estimate and to understand how the error depends on the shape of the exposure meter flux curve. In Section 5, we discuss implications for exoplanet mass measurements, how to optimise observing strategies, and how to deal with and possibly improve existing data.

## 2 PHOTON-WEIGHTED BARYCENTRIC CORRECTION

Using the methods of Wright &amp; Eastman (2014), we are easily able to calculate the barycentric correction z B with a precision better than 1 cms GLYPH&lt;0&gt; 1 for a photon arriving from a star to a telescope on Earth at a given point in time. As emphasised in that paper, z B is a redshift, and it must be applied multiplicatively to the measurement, z meas , i.e.

$$1 + z = ( 1 + z _ { \text{meas} } ) ( 1 + z _ { \text{B} } ), & & ( 1 ) & & \text{The} \,.$$

where z is the 'true' barycentric-motion-corrected redshift. Thus the corrected radial velocity is calculated as

$$v = c z = c \left [ \left ( 1 + \frac { v _ { \text{meas} } } { c } \right ) \left ( 1 + \frac { v _ { \text{B} } } { c } \right ) - 1 \right ], & & ( 2 ) & & V _ { 0 } = \\ & & \text{white} \quad.$$

where we express the barycentric correction as a velocity, v B = cz B .

When taking an exposure of length GLYPH&lt;1&gt; t in a spectrograph, one would ideally calculate the instantaneous barycentric correction for each individual photon at the time it arrives to the telescope, then average all the calculated values for the collected photons in each pixel on the detector. If f ' ' t is proportional to the flux of photons in one pixel, and v ' ' t denotes the barycentric corrected velocity as a function of time, the average corrected velocity of the exposure in that pixel can be expressed as

$$\langle v \rangle = \int _ { t _ { 0 } - \Delta t / 2 } ^ { t _ { 0 } + \Delta t / 2 } v ( t ) f ( t ) d t \Big / \int _ { t _ { 0 } - \Delta t / 2 } ^ { t _ { 0 } + \Delta t / 2 } f ( t ) d t \,, \quad & ( 3 ) \quad & v ( t ^ { \prime } )$$

where t 0 is the geometric midpoint time of the exposure.

Depending on the type of detector, f ' ' t can be obtained in various ways. When infrared arrays are used with'up-theramp' sampling, which means sampling the voltage of each pixel every few ( GLYPH&lt;24&gt; 10) seconds throughout the exposure, it produces a direct measure of the time history of the photon arrivals. This is how instruments like the Habitable-Zone

Planet Finder (Mahadevan et al. 2012) are operated. CCD detectors work differently, and it is not possible to directly monitor the flux received by the detector without reading out the image. Instead one can use the exposure meter flux as a proxy. Landoni et al. (2014) and Blackman et al. (2017) describe how to use a multi-channel exposure meter to measure f ' ' t for a particular wavelength associated with a pixel or a section of the detector.

In this paper, we will refer to f ' ' t as the exposure meter flux curve, although our findings also apply to instruments with infrared arrays or other ways of monitoring flux variations. We will limit ourselves to a single wavelength channel (or a single pixel) and instead focus on quantifying the errors that can result from calculating the barycentric correction at a defined midpoint time rather than using Equation (3). The difference between the geometric and photon-weighted midpoint times, and how the photon-weighted barycentric correction can differ from what is calculated at the midpoint times, is illustrated by an example in Fig. 1.

## 3 QUANTIFYING THE ERRORS

On short timescales (hours), the change in the barycentric correction v B is dominated by the diurnal rotation of the Earth. For a telescope at location ' lat ; lon ' observing a star at position ' GLYPH&lt;11&gt; ; GLYPH&lt;14&gt; ' , the projected velocity of the telescope towards the star is given by

$$\mathbb { m } { \underset { \text{with} } { \text{are} } } \quad - v _ { B } ( t ) \approx \frac { 2 \pi R _ { \oplus } } { 2 4 h } \cos ( \text{lat} ) \cos ( \delta ) \sin \left ( \psi ( t ) \right ) + \text{constant} \quad \quad ( 4 )$$

where ' t ' is the local hour angle, expressed in radians, written here in terms of the local sidereal time (LST):

$$\text{be} \quad \psi ( t ) = \frac { 2 \pi } { 2 4 h } L \text{ST} ( t, \text{lon} ) - \alpha.$$

The hour angle, and thus v B ' ' t , is a sinusoidal variation with a period of 23 93 : hours and amplitude

$$V _ { 0 } = \frac { 2 \pi R _ { \oplus } } { 2 4 h } \cos ( \text{lat} ) \cos ( \delta ),$$

which can be up to 463 ms 1 GLYPH&lt;0&gt; (equatorial velocity of the Earth). The time derivative of ' t ' , which we will soon need, can be approximated as 2 GLYPH&lt;25&gt; GLYPH&lt;157&gt; 24h . For now, we consider a star with v meas = 0 , allowing us to set v = v B

The photon-weighted midpoint time of an exposure, h i t , can be defined in the same way as how we defined the average corrected velocity in Equation (3):

$$\int _ { \substack { t \\ \mathbf t } } ^ { J \, \mathbf v \, \prime } \int _ { \substack { t \\ \mathbf t } } ^ { t _ { 0 } + \Delta t / 2 } t f ( t ) d t \Big | \int _ { t _ { 0 } - \Delta t / 2 } ^ { t _ { 0 } + \Delta t / 2 } f ( t ) d t \,.$$

We can now expand v around h i t with a second-order Taylor polynomial:

$$( 3 ) \quad & \quad t ^ { \searrow \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \cdot \\ & \quad v ( t ^ { \prime } ) \approx v ( \langle t \rangle ) - 2 \pi V _ { 0 } \cos ( \psi ( \langle t \rangle ) ) \frac { t ^ { \prime } - \langle t \rangle } { 2 4 h } \\ \intertext { e. } \underset { \begin{smallmatrix} - \quad \\ \end{smallmatrix} } { \int } \quad & \quad + 2 \pi ^ { 2 } V _ { 0 } \sin ( \psi ( \langle t \rangle ) ) \left [ \frac { t ^ { \prime } - \langle t \rangle } { 2 4 h } \right ] ^ { 2 }. \quad & \quad ( 7 ) \\ \quad \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \$$

If we define the first-order error as the difference between calculating the barycentric correction at the geometric and photon-weighted midpoints, respectively, it can be approximated as:

Figure 1. Left panel shows the instantaneous barycentric correction (black curve) for a star during one 60-minutes exposure at Mauna Kea. The fictional, rising target is chosen such that the exposure starts when the star is 30 GLYPH&lt;14&gt; above the horizon (air mass = 2 0 : ), and such that we are looking due east at the geometric midpoint time of the exposure. The lower panel shows what the exposure meter flux could look like during the exposure, e.g. if a thin cloud passes over the telescope or if the seeing changes. We use the flux curve to compute the photon-weighted midpoint time (red) and the photon-weighted barycentric correction (blue). The zero point of the time axis is set at the geometric midpoint time of the exposure (green). The dashed grey line is the tangent to the barycentric correction curve at the geometric midpoint time, emphasising the curvature of the barycentric correction function. The right panel shows an enlarged view around the centre of the exposure. In this particular case there is a difference of 1 25 m s : GLYPH&lt;0&gt; 1 between the photon-weighted barycentric correction (blue) and the barycentric correction at the photon-weighted midpoint time (red).

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-9c16b917-6082-4ba6-8c77-03fcc230c526.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-ee83c15e-9850-41f5-8acf-b4abc80e5225.jpg)

Figure 2. Left panel : Locations of 25 different observatories that are hosting or will soon be hosting one or more PRV instruments (see e.g. Plavchan et al. 2015; Fischer et al. 2016; Wright &amp; Robertson 2017). Most of these sites are located around latitude GLYPH&lt;6&gt; 30 GLYPH&lt;14&gt; , yielding cos lat ' ' GLYPH&lt;25&gt; 0 87 : . Of these instruments, Mauna Kea (Hawaii) at latitude + 20 GLYPH&lt;14&gt; is the PRV site closest to equator, with cos lat ' ' = 0 94 : . Right panel: Systematic error for the realistic worst-case observing scenario. A rising (or setting) target due East (or West) is observed with the exposure starting (or ending) at 30 GLYPH&lt;14&gt; altitude (air mass = 2 0 : ). Solid and dotted lines indicate the shape of the exposure meter flux curve. Markers indicate the exposure time.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-88c5fc8b-d763-4dbd-82ae-95ee818c6f12.jpg)

$$v ( t _ { 0 } ) - v ( \langle t \rangle ) \approx - 2. 0 \, m \, s ^ { - 1 } \cos ( \text{lat} ) \cos ( \delta ) \cos ( \psi ( \langle t \rangle ) ) \frac { t _ { 0 } - \langle t \rangle } { 1 \, \min } \ \ ( 8 ) \quad \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \$$

A systematic error of up to 2 ms 1 GLYPH&lt;0&gt; for every minute of difference between the geometric and photon-weighted midpoint times is a quite severe error in PRV context. It becomes largest when observing stars on the meridian, which is where most observers usually try to observe their targets. The community is well-aware of this effect, and it has become a standard procedure at almost every PRV instrument to correct for this by recording the photon-weighted midpoint time for each exposure.

By inserting Equation (7) into Equation (3), we can estimate the second-order error, which we define as the difference between the corrected velocity at the photon-weighted midpoint time, v 'h i' t , and the photon-weighted average of the corrected velocity, h v i :

$$v ( \langle t \rangle ) - \langle v \rangle \approx - 2 \pi ^ { 2 } V _ { 0 } \sin ( \psi ( \langle t \rangle ) ) \frac { \langle t ^ { 2 } \rangle - \langle t \rangle ^ { 2 } } { ( 2 4 h ) ^ { 2 } }. \quad \ \ ( 9 ) \quad \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \$$

The second-order error is caused by the curvature of v B , as shown in Fig. 1. In principle, we are moving a smoothing filter over the instantaneous barycentric correction function, using the exposure meter flux curve as the smoothing kernel. This means that the point 'h i t ; h v i' will always land on the inside of the curvature, and the second-order error will thus always have the same sign as the diurnal component of v B (positive east of the meridian and negative west of the meridian). We also note that the magnitude of the secondorder error will usually have a dependence on GLYPH&lt;1&gt; t 2 .

Whereas the first-order error reaches its maximum when the target passes over the meridian, the second-order error gets largest at hour angle GLYPH&lt;6&gt; 6 h and declination zero. This is where the celestial equator intersects with the horizon, and it is thus not observable. A more useful number is the maximum error possible above a certain altitude on the sky, and we therefore transform Equation (9) to horizontal (alt,az) coordinates. With a bit of spherical trigonometry 1 , one can show that cos ' GLYPH&lt;14&gt; ' sin ' ' = GLYPH&lt;0&gt; cos alt ' ' sin az ' ' :

If we consider a uniform exposure (i.e. f ' t ' = constant ), we can evaluate the integrals contained in h i t and h t 2 i :

$$\rangle & \rangle \sim - 1. 3 2 \ m s ^ { - 1 } \cos ( \text{lat} ) \cos ( \delta ) \sin ( \psi ) \left ( \frac { \Delta t } { \text{lh} } \right ) ^ { 2 } & \quad \text{every} \\ & \quad \approx 1. 3 2 \ m s ^ { - 1 } \cos ( \text{lat} ) \cos ( \text{alt} ) \sin ( \text{az} ) \left ( \frac { \Delta t } { \text{lh} } \right ) ^ { 2 } & \quad ( 1 0 ) \quad \text{serva} \\ - \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cd.$$

The error is largest at az = GLYPH&lt;6&gt; 90 GLYPH&lt;14&gt; , i.e. when facing straight east or west. We adapt a lower limit on the altitude of 30 GLYPH&lt;14&gt; (air mass = 2 0 : ) based on the common practice we have seen at various PRV instruments. Since we would like to observe the entire exposure above the altitude limit, the actual limit for the 'alt' parameter in Equation (10) depends on the latitude and exposure time. Therefore, in the following example for Mauna Kea (Hawaii), we set the altitude limit to 37 GLYPH&lt;14&gt; for a 60 min exposure and 33 5 : GLYPH&lt;14&gt; for a 30 min exposure.

As shown in Fig. 2, the majority of current and planned PRV instruments are located at observatories near latitude GLYPH&lt;6&gt; 30 GLYPH&lt;14&gt; , with Mauna Kea at latitude + 20 GLYPH&lt;14&gt; being the closest to equator. If we expose for 30 min at this latitude and employ worst-case altitude and azimuth angles, as described above, the second-order error becomes 0 25 m s : GLYPH&lt;0&gt; 1 ; if we double the exposure time to 60 min, the error quadruples to 1 0 : ms 1 GLYPH&lt;0&gt; . In other words, using the photon-weighted midpoint time for the barycentric correction can lead to significant systematic errors even at the 1 ms 1 GLYPH&lt;0&gt; RV precision level.

For non-uniform flux curves, the result in Equation (10) needs to be altered only by a scaling factor. As an example, it can be shown that a V-shaped flux curve, symmetrically dropping to zero at the geometric midpoint time, yields an effect of exactly 1.5 times the uniform flux result. In general, if photons are missing in the middle of the exposure (e.g. due to a cloud passing over the telescope), the second-order error increases. If photons are concentrated towards one end of the exposure, the error second-order error decreases.

We finally note that although we obtained these results for v meas = 0 in Equation (2), they can safely be generalised to any v meas GLYPH&lt;28&gt; c . The factor ' 1 + v meas GLYPH&lt;157&gt; c ' has a significant effect on the absolute value of the barycentric correction (Wright &amp; Eastman 2014), but its influence on the size of the second-order error is negligible even at the 1 cms GLYPH&lt;0&gt; 1 precision level.

## 4 SIMULATED OBSERVATIONS

We have used the barycorrpy package (Kanodia &amp; Wright 2018) to simulate, under various conditions, the systematic second-order error that would result from not photonweighting the barycentric correction.

For each set of simulations, we let the same observation happen at the same time on multiple sky positions. We generate an exposure meter flux curve with a given shape, exposure time, and sampling cadence. Then, for a range of sky positions, we calculate the barycentric correction both as the photon-weighted average, h v i , and at the photon-weighted midpoint time, v 'h i' t . The integrals are carried out as simple Riemann sums, evaluated at the geometric midpoint of every exposure meter step. Finally, the second-order error is calculated as v 'h i' GLYPH&lt;0&gt; t h v i .

In Fig. 3-5, we use the location of Mauna Kea Observatory ( 19 8222 : GLYPH&lt;14&gt; N, 155 4749 : GLYPH&lt;14&gt; W, 4205 m) and select the time as 10:20 UTC (local midnight) on an arbitrarily chosen night (09-Sep 2017). For each curve, we fix the declination and vary the hour angle. We do this for various exposure times (Fig. 3), declinations (Fig. 4), and exposure meter flux curves (Fig. 5). The solid curves show the simulated errors, and the dotted curves in the background show the analytic estimates from Equation (10).

In Fig. 3 and 5, when the exposure starts at 30 GLYPH&lt;14&gt; altitude, the declination and the minimum hour angle are chosen such that the rising target is situated at az = 90 GLYPH&lt;14&gt; at the geometric midpoint time - and vice versa for the maximum hour angle. In Fig. 4 we simulate various declinations, while still choosing the minimum and maximum hour angle such that the exposure starts and ends at 30 GLYPH&lt;14&gt; altitude or above.

Figure 3. Second-order error, v 'h t i' GLYPH&lt;0&gt; h v i , simulated as function of local hour angle for various exposure times. Each point on the curves corresponds to an exposure at that hour angle. The declination and minimum/maximum hour angle is chosen for each exposure time such that the the telescope is pointing straight east/west at the geometric midpoint of the exposure, while the exposure starts/ends at 30 GLYPH&lt;14&gt; altitude. The exposure meter flux is uniform.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-86a8ddac-a016-458c-b067-ff97ced3d8c7.jpg)

The uniform 60-minute exposure at GLYPH&lt;14&gt; = + 11 7 : GLYPH&lt;14&gt; is present in all three plots, and the size of its error ranges between GLYPH&lt;6&gt; 1 00 : ms 1 GLYPH&lt;0&gt; in agreement with the analytic estimate from Section 3. In Fig. 5, we simulate different shapes of the exposure meter flux curve, and we verify that the error increases by a factor of 1.5 when the photon flux drops and rises again during the exposure.

In Fig. 6 we provide contour maps of the sky for three observatories at different latitudes, showing for a 30-minute uniform exposure how the second-order error depends on the sky position. The 30-minute curve in Fig. 3 corresponds to a horizontal cut through the central diagram in Fig. 6. The area defined by the dashed red line contains all the observations that both start and stop at an altitude higher than 30 GLYPH&lt;14&gt; . The asymmetric shape of the area gives a visual understanding why the maximum effect is not reached at declination zero.

All the simulations in Figures 3-6 are made with 1 Hz sampling cadence. When we generate the same plots with 0.1 Hz sampling, the results are practically identical. The way we interpret this is that our cadence is sufficiently high that the results are not affected by it. Cadence effects can be analysed in more detail by binning real exposure meter flux curves to different bin sizes. We expect the result of such an analysis to depend on the signal-to-noise level and other specifics of the particular instrument, and it is therefore beyond the scope of this manuscript.

## 5 DISCUSSION

In the preceding sections, we have demonstrated the dependence and magnitude of a second-order error in the barycentric correction when applied to longer exposures. In the following, we will discuss the possible implications for exoplanet mass measurements, and we will explore some options

Figure 4. Second-order error, v 'h t i' GLYPH&lt;0&gt; h v i , simulated as function of local hour angle for various declinations. The exposure time is 1800 seconds with a uniform flux curve. All curves end at 37 GLYPH&lt;14&gt; altitude, ensuring that all exposures begin and end above 30 GLYPH&lt;14&gt; altitude. The two ? markers indicate the maximum error for the worst declination with this altitude limit ( GLYPH&lt;14&gt; = 11 7 : GLYPH&lt;14&gt; ). The exposure meter flux is uniform.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-a6487d30-8188-48e2-b931-a739d9015d40.jpg)

Figure 5. Second-order error, v 'h t i' GLYPH&lt;0&gt; h v i , simulated as function of local hour angle for various shapes of the exposure meter flux curve. The shapes are a linear ramp, a uniform exposure, a centred V-shape, and a V-shape offset from the centre. The dotted curves in the background are the analytically approximated values, made by multiplying the uniform result in Equation (10) by 0.67, 1.00, 1.29, and 1.50.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-6fbbccb4-f7d9-419b-80ed-c4f05483831b.jpg)

for mitigating errors when observing and/or analysing data, in cases where the exposure meter is broken or unavailable.

## 5.1 Implications for exoplanet mass measurements

As a new generation of PRV spectrographs are aiming at the 10 cms GLYPH&lt;0&gt; 1 and 1 cms GLYPH&lt;0&gt; 1 regimes, hopes are that astronomers will soon be able to detect Earth-like planets around Sunlike stars. If we were to observe the gravitational tug that our own planet exerts on the Sun, the RV semi-amplitude

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-68e96a01-0532-4410-89bb-09d58e660723.jpg)

Mauna Kea (

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-321b6004-d9ac-4a56-8bf2-5c745c7e3cee.jpg)

20

°N,

155

°W), 2017-09-09 10:20:00 (UTC)

La Silla (

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-d146c3d1-7f30-457a-afa9-fb63becd31c6.jpg)

29

°S,

71

°W), 2017-09-09 04:40:00 (UTC)

Figure 6. Second-order error, v 'h t i' GLYPH&lt;0&gt; h v i , simulated as function of sky position at three different observatories. We simulate a 30 minute exposure with uniform flux, observed at local midnight. The colour indicates the size of the error as a function of the target position at the geometric midpoint time; the dashed red lines indicate the altitude limit that ensures the exposure is completed above 30 GLYPH&lt;14&gt; , while solid red lines are contours, drawn for every 0.05 ms 1 GLYPH&lt;0&gt; . The observatories are Happy Jack (EXPRES), Mauna Kea (e.g. HiRes, KPF, MAROON-X, SPIROU), and La Silla (HARPS, CORALIE); for a more comprehensive list of instruments, see e.g. Plavchan et al. (2015), Fischer et al. (2016), and Wright &amp; Robertson (2017).

would be as little as 9 cms GLYPH&lt;0&gt; 1 . We have shown in this paper that the second-order error on the barycentric correction can easily reach and exceed this level. Next-generation PRV instruments must therefore record the exposure meter flux curve for each observation, allowing the observer or an automated pipeline to calculate the photon-weighted barycentric correction. We recommend that the exposure meter data is stored and archived along with the observed spectrum, allowing the photon-weighted barycentric correction to be fully reproduced and checked at a later time. Although the calculation is trivial, mistakes can still happen (e.g. a broken FITS header) and the observer should always have the option to go back and recalculate the barycentric correction from scratch for each data point. As argued by Wright &amp; Eastman (2014) and elaborated by Blackman et al. (2017), it is also recommended that exposure meters split the light into multiple wavelength channels in order to accommodate wavelength dependent fibre coupling and seeing, as it is done at instruments like ESPRESSO (Pepe et al. 2010; Landoni et al. 2014) and EXPRES (Jurgenson et al. 2016; Blackman et al. 2017).

Over the last decade, instruments like HARPS (Mayor et al. 2003) and HARPS-N (Cosentino et al. 2012) have demonstrated their capability to measure RVs with a precision of 1 ms 1 GLYPH&lt;0&gt; (Fischer et al. 2016). One could argue that the barycentric second-order error is mostly negligible at this level, but as we have seen, situations exist where it becomes a substantial fraction of or even more than 1 ms 1 GLYPH&lt;0&gt; . The barycentric second-order errors for a set of RV measurements do not necessarily form a normal distribution - they depend on weather, atmospheric seeing, observing strategies, and decisions (or mistakes) made by the observers. Our concern is that sometimes these errors may alter the measured velocities in ways that change the apparent amplitude of the RV signal. If that happens, it can affect the modelled mass and/or eccentricity.

Most people seek to observe their targets close to the meridian, where the second-order effect diminishes, but sometimes scheduling constraints requires a target to be observed at higher air mass. One particular example of this is when a target is followed throughout its observing season. When it first becomes visible, in the beginning of the observing season, it will be rising in the east, and it will be observed at a large negative hour angle in order to expose at the highest possible altitude before twilight. As time passes, the target becomes visible all night and observable near the meridian. Towards the end of the season, the target is setting in the beginning of the night, and the observer needs to observe at a large positive hour angle. For a target near the celestial equator, with 30 minutes exposure time, the insufficiently corrected velocities could have an apparent GLYPH&lt;6&gt; 25 cms GLYPH&lt;0&gt; 1 'curl' with opposite signs in the beginning and end of the observing season.

A similar situation could happen when a target instead of once per night is visited twice per night. An observer may try to separate the observations in time as much as possible by observing in the beginning and end of the night at the highest feasible air mass. This could lead to systematic barycentric correction errors where the early and late observations each are offset by e.g. 25 cms GLYPH&lt;0&gt; 1 with opposite signs. Sometimes this type of data is fitted with a nightly offset to account for trends due to stellar activity

(Pepe et al. 2013; Howard et al. 2013). Any RV variation within a night would thus be interpreted as Keplerian motion, and the the barycentric error could propagate to the fitted semi-amplitude.

As long as we calculate and apply the barycentric correction with proper photon weights, the RV errors due to barycentric correction will be limited to less than 1 cms GLYPH&lt;0&gt; 1 , as in Wright &amp; Eastman (2014), and the situations described above do not cause any particular problems with respect to barycentric correction. Otherwise, there is a risk that data from current-generation PRV instruments near 1 ms 1 GLYPH&lt;0&gt; precision could sometimes be affected by second-order errors, and we advise caution.

## 5.2 Mitigation strategies

In cases where no flux curve is supplied with an exposure, there are still ways to mitigate the consequent errors in the barycentric correction.

If the observer knows in advance that the exposure meter is broken or for some reason not available, it may be desirable to modify the observing strategy or even switch to other targets. We know from Equations (8) and (9) that targets near the celestial equator are more strongly affected by both errors, so targets at high declination would be better suited for a night like that. We also know that the size of the first-order error can be reduced by shortening the exposure time or observing at larger local hour angle; the latter however comes at the price of higher air mass. In order to keep the flux as uniform as possible, observing through clouds without an exposure meter is not advisable.

Observers may also find themselves at an instrument where the exposure meter only delivers the photon-weighted midpoint time, h i t , and does not preserve the flux curve. In this case, the first-order effect is accounted for by calculating the barycentric correction at h i t , and our only concern is the second-order effect. Since it is proportional to sin ' ' , it can be easily minimised by observing near to the meridian. Decreasing the exposure time or choosing bright targets with shorter exposure times will also help. Fig. 6 gives a visual idea about what areas of the local sky are within a given error margin at three different observatory sites.

When dealing with PRV data that supplies the photonweighted midpoint time and yet no flux curve, the simplest assumption is that the flux was uniformly distributed in time . 2 This allows us to add a second-order correction to the calculated barycentric correction, v 'h i' t . To obtain h v i , we can either subtract the result of Equation (10) or get the difference from a simple numerical simulation as in Section 4. In fact, if we do not apply a second-order correction to v 'h i' t , we implicitly assume that all photons arrived in the same instant at the photon-weighted midpoint time. We argue that the assumption of a uniform flux distribution is much more realistic. Fig. 5 shows how the second-order error

2 One may (rightfully) argue that a uniform flux distribution over the entire exposure time entails that h t i = t 0 . For a more stringent approach, assume a step function set to zero in the beginning or the end of the exposure, such that the photon-weighted midpoint time h t i is reproduced and the exposure time GLYPH&lt;1&gt; t is effectively shortened by 2 j h t i GLYPH&lt;0&gt; t 0 j .

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00991v1-figure-24316731-0fc5-463f-85e0-c9953e0a5ac3.jpg)

Estimated second-order error, uniform flux (cm s

)

Figure 7. For nearly one year of EXPRES observations, we compare the actual second-order error of each measurement with the uniform flux approximation. The multi-channel exposure meter data has been binned to one channel. EXPRES mainly observes bright targets, so the exposure times are short ( GLYPH&lt;20&gt; 1200 seconds), resulting in errors of only a few cms GLYPH&lt;0&gt; 1 . Nevertheless, the data shows that the uniform flux curve is a good approximation if we for some reason do not have access to the exposure meter flux curve when we want to calculate the second-order error.

of three different flux scenarios are all closer to the uniform scenario than zero.

In Fig. 7 we calculate the second-order error using real exposure meter data from 1315 observations with EXPRES (see Blackman et al. 2019), combining all colour channels to one. The data was obtained during commissioning between March 2018 and February 2019. No exposures were longer than 20 minutes, and as a result the second-order errors are relatively small. We plot the real second-order errors against the errors estimated with Equation (10) for a uniform flux curve. The plot shows that within 1 cms GLYPH&lt;0&gt; 1 , the uniform flux approximation is indeed a realistic assumption, and much more realistic than setting the error to zero. There is no indication that longer exposure times would change this picture.

In the rare case of dealing with data that completely lacks information about the temporal flux distribution, there is no way of knowing if the individual observations have large first-order barycentric errors. One can in principle use Equation (8) to estimate the RV error bars, but it requires a qualified guess at the offset 'h i GLYPH&lt;0&gt; t t 0 ' . In our experience, a typical offset is of order 5% of the exposure time, with some larger, negative offsets from when the star was lost on the fibre before the end of the exposure. The actual distribution depends on the instrument and the local sky conditions.

## 6 CONCLUSIONS

This paper describes how to correctly apply a barycentric correction to an exposure using the exposure meter flux curve or an equivalent measure. Rather than following the common practice of computing a single value at the photonweighted midpoint time, we argue why the correct approach is to apply the photon weights to the barycentric correction itself. The systematic difference between the two methods can be described as a second-order effect that leads to systematic RV errors of up to 0.25 ms 1 GLYPH&lt;0&gt; for a 30 minute exposure and 1.0 ms 1 GLYPH&lt;0&gt; for a 60 minute exposure. Certain observing scenarios, when the photon flux drops and rises again during an exposure, can amplify the error by a factor of GLYPH&lt;24&gt; 1.5.

Considering that next-generation PRV spectrographs are aiming for the 10 cms GLYPH&lt;0&gt; 1 regime and beyond, it will be absolutely necessary to account for this effect. Instruments of less precision (1 ms 1 GLYPH&lt;0&gt; ) are likewise affected for some observations, depending on where on the sky and the exposure time. For instruments that already measure the photon-weighted midpoint, it should be somewhat trivial to compute the photon-weighted barycentric correction in order to eliminate unnecessary systematic errors. We also show that it is possible to partially mitigate second-order errors even in archive data that was stored without the flux curve.

## ACKNOWLEDGEMENTS

We would like to thank the referee, Artie Hatzes, for carefully reading the manuscript and supporting its conclusions. We also thank the members of the EXPRES and HARPS-N GTO teams for useful input and fruitful discussions. Our simulations involved an extensive amount of coordinate transformations, made possible with the excellent astropy.coordinates module 3 (Astropy Collaboration et al. 2013; Price-Whelan et al. 2018). All figures were prepared with Matplotlib 4 (Hunter 2007). The Center for Exoplanets and Habitable Worlds is supported by the Pennsylvania State University, the Eberly College of Science, and the Pennsylvania Space Grant Consortium.

## REFERENCES

Astropy Collaboration et al., 2013, A&amp;A, 558, A33

Blackman R. T., Szymkowiak A. E., Fischer D. A., Jurgenson C. A., 2017, ApJ, 837, 18

Blackman R. T., Ong J. M. J., Fischer D. A., 2019, AJ, 158, 40 Cosentino R., et al., 2012, in Ground-based and Airborne Instrumentation for Astronomy IV. p. 84461V, doi:10.1117/12.925738

Fischer D. A., et al., 2016, Publications of the Astronomical Society of the Pacific, 128, 066001

Howard A. W., et al., 2013, Nature, 503, 381

Hunter J. D., 2007, Computing In Science &amp; Engineering, 9, 90 Jurgenson C., Fischer D., McCracken T., Sawyer D., Szymkowiak A., Davis A., Muller G., Santoro F., 2016, in Ground-based and Airborne Instrumentation for Astronomy VI. p. 99086T ( arXiv:1606.04413 ), doi:10.1117/12.2233002

Kanodia S., Wright J., 2018, Research Notes of the American Astronomical Society, 2, 4

- 3 http://astropy.org

4 http://matplotlib.org

Landoni M., Riva M., Pepe F., Conconi P., Zerbi F. M., Cabral A., Cristiani S., Megevand D., 2014, in Ground-based and Airborne Instrumentation for Astronomy V. p. 91478K, doi:10.1117/12.2056406

Mahadevan S., et al., 2012, in Ground-based and Airborne Instrumentation for Astronomy IV. p. 84461S ( arXiv:1209.1686 ), doi:10.1117/12.926102

Mayor M., et al., 2003, The Messenger, 114, 20

Pepe

F.

A., et

al.,

2010, in

Proc.

SPIE.

doi:10.1117/12.857122

Pepe F., et al., 2013, Nature, 503, 377

Plavchan P., et al., 2015, arXiv e-prints, p. arXiv:1503.01770

Price-Whelan A. M., et al., 2018, AJ, 156, 123

Wright J. T., Eastman J. D., 2014, PASP, 126, 838

Wright J. T., Robertson P., 2017, Research Notes of the AAS, 1, 51

This paper has been typeset from a T E X/LT E X file prepared by A the author.

p.

77350F,
René Tronsgaard Lars A. Buchhave Jason T. Wright Jason D. Eastman Ryan T. Blackman
Photon-weighted barycentric correction and its importance for precise radial velocities
When applying the barycentric correction to a precise radial velocity
measurement, it is common practice to calculate its value only at the
photon-weighted midpoint time of the observation instead of integrating over
the entire exposure. However, since the barycentric correction does not change
linearly with time, this leads to systematic errors in the derived radial
velocities. The typical magnitude of this second-order effect is of order 10 cm
s$^{-1}$, but it depends on several parameters, e.g. the latitude of the
observatory, the position of the target on the sky, and the exposure time. We
show that there are realistic observing scenarios, where the errors can amount
to more than 1 ms$^{-1}$. We therefore recommend that instruments operating in
this regime always record and store the exposure meter flux curve (or a similar
measure) to be used as photon-weights for the barycentric correction. In
existing data, if the flux curve is no longer available, we argue that
second-order errors in the barycentric correction can be mitigated by adding a
correction term assuming constant flux.