# A detailed study of the barium central star of the planetary nebula Abell 70

David Jones${}^{1,2}\href https://orcid.org/0000-0003-3947-5946{}$ , Henri M.
J. Boffin${}^{3}\href https://orcid.org/0000-0002-9486-4840 {}$, Alex J.
Brown${}^{4}\href https://orcid.org/0000-0002-3316-7240{}$, Jiri
Zak${}^{3}\href https://orcid.org/0000-0001-9416-9007{}$, George
Hume${}^{4,5}\href https://orcid.org/0000-0002-8518-4974{}$, James
Munday${}^{6}\href https://orcid.org/0000-0002-1872-5398{}$, and Brent
Miszalski${}^{7}\href https://orcid.org/0000-0003-2561-6306{}$
1Instituto de Astrofísica de Canarias, E-38205 La Laguna, Tenerife, Spain
2Departamento de Astrofísica, Universidad de La Laguna, E-38206 La Laguna,
Tenerife, Spain
3European Southern Observatory, Karl-Schwarzschild Strasse 2, 85748 Garching,
Germany
4Department of Physics and Astronomy, University of Sheffield, Sheffield, S3
7RH, UK
5Isaac Newton Group of Telescopes, Apartado de Correos 368, E-38700 Santa Cruz
de La Palma, Spain
6Department of Physics, Gibbet Hill Road, University of Warwick, Coventry CV4
7AL, United Kingdom
7Australian Astronomical Optics - Macquarie, Faculty of Science and
Engineering, Macquarie University, North Ryde, NSW 2113, Australia
E-mail<EMAIL_ADDRESS>
(Accepted xxxx xxxxxxxx xx. Received xxxx xxxxxxxx xx; in original form xxxx
xxxxxxxx xx)

###### Abstract

We present a detailed study of the barium star at the heart of the planetary
nebula Abell 70. Time-series photometry obtained over a period of more than
ten years demonstrates that the barium-contaminated companion is a rapid
rotator with temporal variability due to spots. The amplitude and phasing of
the photometric variability changes abruptly, however there is no evidence for
a change in the rotation period (P = 2.06 d) over the course of the
observations. The co-addition of 17 high-resolution spectra obtained with VLT-
UVES allow us to measure the physical and chemical properties of the
companion, confirming it to be a chromospherically-active, late G-type sub-
giant with more than +1 dex of barium enhancement. We find no evidence of
radial velocity variability in the spectra, obtained over the course of
approximately 130 d with a single additional point some 8 years later, with
the radial velocities of all epochs approximately $-$10 km s-1 from the
previously measured systemic velocity of the nebula. This is perhaps
indicative that the binary has a relatively long period (P $\gtrsim$ 2 yr) and
high eccentricity ($e\gtrsim$ 0.3), and that all the observations were taken
around radial velocity minimum. However, unless the binary orbital plane is
not aligned with the waist of the nebula or the systemic velocity of the
binary is not equal to the literature value for the nebula, this would imply
an unfeasibly large mass for the nebular progenitor.

###### keywords:

planetary nebulae: individual: PN A66 70, PN G038.1$-$25.4 – planetary
nebulae: general – stars: AGB and post-AGB – stars: chemically peculiar –
accretion, accretion discs

††pagerange: A detailed study of the barium central star of the planetary
nebula Abell 70–A.5††pubyear: 2022††pagerange: A detailed study of the barium
central star of the planetary nebula Abell 70–A.5

##  1 Introduction

It is now clear that binary stars hold the key to understanding the formation
of many planetary nebulae (PNe; Jones & Boffin, 2017a). However, with much of
the recent focus being placed on common envelope evolution (Boffin & Jones,
2019), the role of wider binaries remains almost completely unconstrained
(Tyndall et al., 2013). To date, only three wide binary central stars of PNe
have confirmed periods (Jones et al., 2017) with other systems either being so
wide that they are resolved (Ciardullo et al., 1999) or displaying composite
spectra with giant companions (Tyndall et al., 2013). Four central star
systems are known to be barium stars (Bond et al., 2003; Miszalski et al.,
2012; Miszalski et al., 2013; Löbling et al., 2019) with only one having a
known orbital period (LoTr 5; Jones et al., 2017; Aller et al., 2018). Barium
stars are mostly – but not only: see, e.g., Escorza et al. (2019) – giant
stars of spectral type G-K which display enhanced abundances of carbon and
s-process elements like barium and strontium, and are now known to all be
long-period binaries (100d $\lesssim$ Porb $\lesssim$ 104d; McClure et al.,
1980; Jorissen et al., 2019a). The chemical contamination of barium stars is
believed to be due to accretion of chemically enriched material from an
evolved binary companion, likely through wind (Boffin & Jorissen, 1988) or
wind Roche lobe overflow (WRLOF; Theuns & Jorissen, 1993; Nagae et al., 2004).
WRLOF will occur when the acceleration radius of the asymptotic giant branch
(now a white dwarf in these barium star systems) star’s stellar wind is
comparable or greater than its Roche lobe radius (Mohamed & Podsiadlowski,
2012), meaning that the wind itself is strongly influenced by the binary
potential and can be accreted onto the companion via the first Lagrangian
point. This accreted material chemically contaminates the surface of the
companion as well as increases its spin rate due to the conservation of
angular momentum (Theuns et al., 1996). Intriguingly, the amount of material
required to account for the observed contamination is often in excess of that
which would be expected to result in critical rotation rates in the companion,
providing an indication that significant (and as yet not understood) angular
momentum losses must be experienced in these systems (Matrozis et al., 2017).
PNe with barium central stars offer an important window into this process as
the presence of the short-lived nebula ($\tau\leq$ 30 000 yr; Jacob et al.,
2013) means that the mass transfer has occurred too recently for significant
changes in the companion’s spin rate, and moreover the nebula itself traces
the mass-loss history of the system as it is formed from the material that has
escaped the binary potential (Jones, 2018).

In this paper, we present a detailed study of the central star of the
planetary nebula Abell 70 (PN G038.1$-$25.4, or A66 70, hereafter A 70), which
was revealed by Miszalski et al. (2012) to be a binary comprised of a hot
white dwarf, revealed by GALEX photometry, and a roughly G8IV-V barium star
secondary, which has been heavily polluted by s-process material from the
nebular progenitor ([Ba/Fe] overabundance $\sim$0.5 dex). A 70 itself is a
southern PN most well known for its striking appearance, earning it the
nickname the “diamond ring" due to its ring-like form and superposition with a
background elliptical galaxy (which forms the diamond of the diamond ring,
Miszalski et al., 2012). The kinematical study of Tyndall et al. (2013) showed
A 70 to most likely comprise a faint bipolar shell, the waist of which is
encircled by a bright toroid, the deprojection of which results in an
estimated nebular inclination of 30°$\pm$ 10°. They also determined a systemic
heliocentric velocity for the nebula of $-$73$\pm$4 km s-1 and a kinematical
age of approximately 15,000 years for the distance of 6.04$\pm$2.13 kpc
derived by Frew et al. (2016).

The central star of Abell 70 appears in the Gaia Data Release 3 (DR3; source
id: 6907822573352460032, Gaia Collaboration, 2022), where it was not
identified as binary nor was it bright enough for the radial velocity data to
be released. However, Gaia DR3 does contain stellar parameters for the central
star of A 70 (Creevey et al., 2022; Fouesneau et al., 2022), derived using the
General Stellar Parameterizer from Photometry (GSP-Phot) module of the
Astrophysical parameters inference system (Apsis). GSP-Phot simultaneously
fits the BP/RP spectra, parallax and apparent magnitude, using isochrone
models to constrain the stellar properties. The corresponding parameters, as
well as the geometric parallax from EDR3 (Bailer-Jones et al., 2021), for the
central star of A 70 are summarised in Tab. 1, however these should be used
with caution as will be discussed in Sec. 4.

##  2 Photometry

###  2.1 New i-band photometry

We supplement the literature photometry of Bond & Ciardullo (2018) with
photometric observations of the central star of A 70 made in the $i$-band at
various facilities between August 2010 and July 2022.

The first of the supplementary observations were acquired on the nights of
2010 August 11–13 and 15–17 with the 1.9-m Radcliffe Telescope of the South
African Astronomical Observatory (SAAO) and the SAAO CCD instrument. The STE4
CCD was employed along with the Bessell I filter, providing a FOV of
approximately 2.5′$\times$2.5′ and a binned pixel scale of 0.28″ pixel-1.

Further observations were then acquired on the nights of 2013 June 3–7 with
the European Southern Observatory’s 3.6-m New Technology Telescope (ESO-NTT)
and the ESO Faint Object Spectrograph and Camera v.2 (EFOSC2) instrument
(Buzzoni et al., 1984; Snodgrass et al., 2008). EFOSC2 was employed in its
imaging mode along with a Gunn $i$-band filter (ESO ID #705) to provide a FOV
of 4.1′$\times$4.1′ with a binned pixel scale of 0.24″ pixel-1.

Service mode observations were acquired with the Las Cumbres Observatory
Global Telescope (LCOGT; Brown et al., 2013) on 2015 April 1, May 14, and June
12, 14, 18, 19, 20, 21, 22, 25 and 27 with the 1m network, and on 2015 June 24
and 27 on the 2m network. All observations were taken with a Cousins I-band
filter, with 1m network observations employing SBIG cameras for a
15.8′$\times$15.8′ FOV and an unbinned pixel scale of 0.23″ pixel-1, while the
2m network observations employed Spectral cameras for a 10.5′$\times$10.5′ FOV
and a binned pixel scale of 0.30″ pixel-1.

Finally, observations were taken with the 2-m Liverpool Telescope (LT) and the
IO:O instrument on the nights of 2022 May 25, 27, 28, 29 and 31, June 1, 11,
12, 17 and 18, and July 4, 7 and 16. All observations were taken unbinned
(10′$\times$10′ FOV and 0.15″ pixel-1) and with a Sloan $i^{\prime}$ filter.

All data were bias and flat field corrected using standard reduction routines
or the instrument pipelines (where available). Differential aperture
photometry of the central star of A 70 was then performed using the astropy-
affiliated photutils python
package111https://photutils.readthedocs.io/en/stable/.

Bond & Ciardullo (2018) identify a periodicity of 2.06 d in their data, with a
similar period clearly evident in the newly acquired data. However, no single
period and time of first minima could be found that fits all the data (perhaps
not unexpected given that the variability is clearly not constant, see Sec.
2.3 for further discussion). In order to better constrain the period, as well
as the changes between observing epochs, we focus on the periodicitiy of all
the data taken between 2011 and 2015 (shown in the lower left panel of Fig.
1). Given the nature of the variability, we measure the period via the reduced
$\chi^{2}$ of a sinusoidal fit as a function of frequency (ideally suited to
sinusoidal variability, e.g., Horne et al., 1986; Munday et al., 2020). The
resulting periodogram is shown in Fig. 2. Two strong minima in the reduced
$\chi^{2}$ are found, the strongest at 2.061$\pm$0.005 d (consistent with the
period of Bond & Ciardullo, 2018) and a second slightly shallower at
1.938$\pm$0.005 d. The data from 2010, combining that of Bond & Ciardullo
(2018) and our own from SAAO, phases well on the first period but not the
second, lower period, likely indicating that this is an alias. As such, we
conclude that the rotation period is indeed $\approx$2.06 d and that there is
no strong evidence for a change in period before and after the phase of
negligible variability in early 2011. There is, however, a clear shift in
phase between the 2010 data and the later data used to derive the period (see
Fig. 2). There is no statistically significant evidence for a shift between
the 2011–2015 data and the later data from 2022, although the lack of data in
the intervening seven years as well as the uncertainty on the period mean that
we cannot exclude a phase shift between these two data sets. The 2022 data
does, however, appear to have a larger photometric amplitude, more consistent
with the 2010 data. Nonetheless, if the variability did remain constant
throughout this time, we can refine the period to P=2.060561$\pm$0.000002 d.

Figure 1: The combined $i$-band photometry (including the data from Bond &
Ciardullo, 2018) demonstrating the constant periodicity (P=2.061 d), but
changing phasing and amplitude, of the variability. All data are phased
according to the ephemeris derived from the 2010 data using the period derived
from the 2011b, 2013 and 2015 data.

###  2.2 Survey photometry

The field containing A 70 has been observed over relatively long periods of
time by a number of all-sky surveys. Unfortunately, the relatively faint
central star and bright nebula combined with the large pixel scales employed
by the majority of these surveys means their photometry is not particularly
useful in directly probing the variability in A 70 – the 2 d period would not
be uniquely identified in any of the datasets. Nonetheless, they can be used
to probe the apparent changes in phase between the epochs presented in Fig. 1.
The details of this analysis are outlined in Appendix A, with the main
conclusion being that we cannot rule out an additional shift in phase between
the 2015 and 2022 epochs of photometry.

###  2.3 Origin of the changing variability

The observed photometric period – presumably the rotation period of the less
evolved component of the binary nucleus of A 70 – is consistent with the
periodicities observed in other fast-rotating companions to central stars of
PNe. Miszalski et al. (2013) found a period of 5.5 d for Hen 2-39, while
similar periods have been found for giant components of WeBo 1 (4.7 d; Bond et
al., 2003), LoTr 1 (6.4 d; Tyndall et al., 2013) and LoTr 5 (5.95 d; Aller et
al., 2018). The study of Miszalski et al. (2012) found that the companion in A
70 is either a main sequence or sub-giant G-type star, both of which would be
below break-up for a rotational period of 2 d (although a sub-giant would be
fairly close).

However, no other rapidly rotating companion has been shown to present with
changing variability. Within observing periods, where the data is sufficiently
precise (e.g., the 2022 LT data or the 2010 data from Bond & Ciardullo, 2018),
there seems to be some evidence for additional scatter beyond the over-arching
sinusoidal periodicity. This may simply be underestimated uncertainties which
do not account for the contamination from the surrounding nebula, however the
sinusoidal variability itself also seems to come and go as well as change in
phasing. In early 2011, the variability drops below a detectable level, and
then returns in late 2011 with a different phasing but similar amplitude. The
new data obtained in 2013 and 2015 phase well with the late 2011 data,
indicating that the variability was relatively stable over this four year
period. The long time span between the 2015 and 2022 observations means that
we cannot conclusively determine whether there has been a phase shift during
that time (small uncertainties in the rotation period equate to significant
changes in phase). Nonetheless, similar shifts in phase of photometry
variability are observed in other magnetically active stars (not at the
centres of PNe), for example the giant component of IT Com where the change
was associated with a periastron passage (Oláh et al., 2013). The fact that
the central star of A 70 is the first of the rotationally variable PN central
stars to be found to show changes in its photometric variability is more
likely a consequence of the lack of long-term monitoring for other central
stars rather than a special property of the central star of A 70.

Figure 2: The reduced $\chi^{2}$ of a sinusoidal fit as a function of period for the 2011b, 2013 and 2015 data (shown folded on this period in Fig. 1. The best fitting period (P=2.061$\pm$0.005d) is marked by the vertical red line. Table 1: The Gaia DR3 parameters of the central star of A 70. Parameter | A 70
---|---
$T_{\mathrm{eff}}$ (K) | 5144${}^{+19}_{-17}$
log $g$ | 4.784${}^{+0.012}_{-0.015}$
[Fe/H] | $-2.73^{+0.26}_{-0.32}$
GSP-Phot Aeneas distance (kpc) | 1.47${}^{+0.02}_{-0.01}$
GSP-Phot Aeneas radius (R⊙) | 0.506${}^{+0.05}_{-0.04}$
Parallax, $\varpi$ (mas) | 0.25$\pm$0.12
Geometric parallax distance (kpc) | 3.6${}^{+1.5}_{-1.0}$
Bp magnitude | 18.13$\pm$0.03
G magnitude | 17.66$\pm$0.01
Rp magnitude | 17.08$\pm$0.02

##  3 Spectroscopy

The central star of A 70 was observed 16 times between 2011 May 28 and
September 4, and once more on 2019 October 17 (for exact dates see Tab. 2)
using the Ultraviolet and Visual Échelle Spectrograph (UVES; Dekker et al.,
2000) on ESO’s Very Large Telescope (VLT) Unit Telescope 2 (UT2), also known
as Kueyen. All exposures were of 3000 s through a 1″ slit. The dichroic #2 was
employed, with the red arm centred at 7600Å (with a BK7_5 blocking filter) and
the blue arm centred at 4370Å (with a HER_5 below-slit filter). This set-up
provides continuous spectral coverage in the ranges 3800-5000Å (blue-arm) and
5700-9400Å (red-arm, with a small gap in coverage between CCDs around 7600Å),
with a resolution of approximately 40 000. The data were reduced using the
standard ESO pipeline.

Miszalski et al. (2012) presented intermediate-dispersion spectroscopy
obtained with the Gemini Multi-Object Spectrograph (GMOS; Hook et al., 2004)
on the Gemini South Telescope, and with the FOcal Reducer and low-dispersion
Spectrograph (FORS2; Appenzeller et al., 1998) on ESO’s VLT Unit Telescope 1
(UT1, Antu). In order to analyse them with the same methodology as the
aforementioned UVES spectra, these data (as well as a further two FORS2
spectra both obtained on 2012 July 8 with the same instrumental setup) were
downloaded from their respective archives and re-reduced.

###  3.1 Radial velocities

The reduced spectra were continuum subtracted, corrected to heliocentric
velocity and cross-correlated against a G8IV template (in accordance with the
spectral type determined by Miszalski et al., 2012) created using spectrum222
http://www.appstate.edu/~grayro/spectrum/spectrum.html (Gray & Corbally,
1994). Each cross-correlation function (CCF) displayed a clear peak of width
consistent with a high $v$ sin $i$ ($\sim$30–50 km s-1) as expected from the
rotational period derived in section 2. For the highest signal-to-noise
spectra the shape of the CCF could be seen to deviate from a single gaussian
peak, perhaps indicative of spots crossing the surface of the star (as
expected from the observed photometric variability). This hypothesis is
difficult to test given the relatively sparse phase sampling, as well as the
low signal-to-noise of the individual spectra and, particularly, individual
lines - thus preventing a detailed doppler tomography of the stellar surface.
In any case, the shape of the CCF does not seem to correlate with the
photometric phase of the observation (see Fig. 3).

Figure 3: The CCFs of three spectra obtained close to photoemtric minimum (the
MJDs of the spectra are in the legend, but each is $\approx$0.1 in phase away
from the photometric minimum), highlighting that the shape of the CCF does not
appear to correlate with the photometric phase.

In order to minimise the effects of the possible spot(s) on the shape of the
CCF, the final radial velocities were measured using a bisector of the entire
peak, although some “jitter” may still be present. The resulting heliocentric
radial velocities (RVs) are presented in Tab. 2. Note that just as in
Miszalski et al. (2012), the FORS2 spectra showed significant signs of
flexure, and so the derived RVs were corrected using gaussian fits to the [O
iii] 4959Å and 5007Å emission lines and assuming a nebular systemic velocity
of $-$73 km s-1 (Tyndall et al., 2013).

The resulting RVs show no obvious periodic variability, with the UVES
measurements clustering around a heliocentric velocity of between $-$85 and
$-$82 km s-1 (See Fig. 4. Given that the typical periods of barium stars are
in the range of hundreds to tens of thousands of days (Jorissen et al., 1998,
2019b), one would perhaps expect to see some signs of variability in the UVES
RVs. The 2011 data spans the shorter end of the barium star period
distribution (P$\sim$100 days), such that periods of a year or shorter can be
ruled out (unless the orbital inclination is particularly low). However, with
only a single data point some 8 years after the initial 2011 data, a number of
orbital solutions are still possible. If one reasonably assumes that the
nebular systemic velocity, $-$73 km s-1 (Tyndall et al., 2013), is also the
systemic velocity of the binary, then the semi-amplitude of the binary would
be roughly 10 km s-1 with all data being taken close to quadrature (i.e., at
the RV minimum). The lack of observed variability would then be a consequence
of a relatively long period and high eccentricity. Restricting to just the
UVES data, periods of approximately 2, 4 and 8 years would all give plausible
fits to the data. However, including the generally more uncertain GMOS and
FORS2 radial velocities clearly favours a longer period (see Fig. 5). This
longer period would also be supported by the previously discussed photometric
variability where the variability was found to be stable for at least four
years (see Sec. 2). Note, however, that the similar changes in photometric
variability of other magnetically-active stars may have been associated with
periastron passages, perhaps implying that the time of pericenter passage for
the central stars of A 70 should be roughly coincident with the 2011 period of
negligible photometric variability. This is not the case for the fit presented
in Fig. 5 and no reasonable fit could be obtained constraining the time of
pericenter passage. Ultimately, given the low precision of the GMOS and FORS2
RVs, and the non-variability of the UVES RVs, any fit is highly speculative.

Figure 4: Heliocentric radial velocity measurements of the central star of A 70. The right-hand axis shows the deviation of the stellar heliocentric radial velocity from the heliocentric systemic radial velocity of the nebula as measured by Tyndall et al. (2013). Table 2: A table of measured heliocentric radial velocities versus heliocentric Julian date for the polluted giant secondary of A 70. Heliocentric | Radial velocity | Instrument | Approx. S/N
---|---|---|---
MJD (days) | ( km s-1) |  | of spectrum
54944.32710 |  $-$73.93 |  $\pm$9.27 | GMOS | 9.3
55026.41908 |  $-$77.83 |  $\pm$12.25 | FORS2 | 4.0
55364.42193 |  $-$77.16 |  $\pm$7.73 | FORS2 | 10.5
55366.41773 |  $-$74.03 |  $\pm$7.57 | FORS2 | 12.8
55679.35887 |  $-$83.45 |  $\pm$0.71 | UVES | 5.4
55680.34523 |  $-$84.36 |  $\pm$0.84 | UVES | 7.1
55700.36627 |  $-$82.13 |  $\pm$3.65 | UVES | 2.8
55700.40375 |  $-$83.76 |  $\pm$3.24 | UVES | 2.7
55726.23343 |  $-$83.09 |  $\pm$1.19 | UVES | 4.3
55726.29200 |  $-$83.37 |  $\pm$0.73 | UVES | 4.0
55765.28905 |  $-$83.28 |  $\pm$1.00 | UVES | 2.6
55765.32868 |  $-$85.63 |  $\pm$0.31 | UVES | 2.3
55775.11918 |  $-$82.46 |  $\pm$1.15 | UVES | 3.4
55775.15412 |  $-$84.28 |  $\pm$0.08 | UVES | 4.8
55783.09244 |  $-$82.10 |  $\pm$2.84 | UVES | 3.2
55783.13296 |  $-$81.44 |  $\pm$2.28 | UVES | 4.6
55791.22318 |  $-$82.64 |  $\pm$1.03 | UVES | 6.0
55791.26339 |  $-$83.06 |  $\pm$1.06 | UVES | 5.7
55808.06681 |  $-$82.65 |  $\pm$1.50 | UVES | 6.4
55808.10472 |  $-$84.15 |  $\pm$0.34 | UVES | 5.9
56116.40453 |  $-$78.08 |  $\pm$7.73 | FORS2 | 17.4
56116.41831 |  $-$77.25 |  $\pm$6.11 | FORS2 | 15.9
58773.12292 |  $-$82.77 |  $\pm$2.27 | UVES | 3.1
Figure 5: A speculative orbital solution with a period of approximately 8
years and an eccentricity of 0.3. The horizontal dashed line is the nebular
velocity derived by Tyndall et al. (2013), taken to be the systemic velocity
of the binary for each of the fits.

###  3.2 Spectral emission features

Ca H and K are in emission in the UVES spectra (see Fig. 6), consistent with
chromospheric activity. Similarly, just as in other rapidly rotating
companions (e.g. LoTr 1 and LoTr 5; Tyndall et al., 2013; Jasniewicz et al.,
1994), H$\alpha$ shows broad emission with strong absorption superimposed (see
Fig. 7). The full-width at half-maximum of the emission ($\sim 130$ km s-1,
although this is challenging to measure accurately given the strength of the
absorption feature) means that it is unlikely to originate from the
chromosphere of the companion but may, instead, be associated with a disc-like
structure similar to those found in symbiotic stars (Tyndall et al., 2013, and
references therein). Ultimately, the signal-to-noise of the spectra is
insufficient to draw stronger conclusions other than to say that the central
star of A 70 is spectrally very similar to others of the class.

Figure 6: Calcium H and K in emission in the co-added UVES spectrum of the
central star of A 70.  Figure 7: The H$\alpha$ emission/absorption feature in
the co-added UVES spectrum of the central star of A 70.

###  3.3 Spectral modelling

To determine the properties of the bright component of A70, we co-added all
our spectra (after correcting for the derived radial velocities) and
normalised the resulting spectrum. The iSpec software package (Blanco-Cuaresma
et al., 2014; Blanco-Cuaresma, 2019)333https://www.blancocuaresma.com/s/iSpec
was then used to model the spectrum and to probe stellar parameters and
abundances. The stellar parameters ($T_{\mathrm{eff}}$, log $g$, $v$ sin $i$,
[M/H]) were derived following the workflow described in Blanco-Cuaresma et al.
(2014) using a $\chi^{2}$-minimization comparison with synthetic spectra
calculated using the spectrum software (Gray & Corbally, 1994) and MARCS model
atmospheres (Gustafsson et al., 2008). The initial effective temperature and
surface gravity chosen to begin the fitting process were those derived by
Miszalski et al. (2012, $T_{\mathrm{eff}}=5300$ K and $\log g$=3.9,
respectively) with solar metallicity, and a $v$ sin $i$ consistent with the
spectral type of Miszalski et al. (2012) and the rotation period measured in
Sec. 2 (40 km s-1, assuming an inclination, $i$, comparable to that of the
nebula as determined by Tyndall et al. 2013). The fitting procedure was then
allowed to undergo as many iterations as required to converge (in this case
seven iterations).

It is important to note that only iron lines (Fe i and Fe ii) were used to
derive the stellar parameters, thus ensuring that any nebular contamination
(or poor nebular subtraction by the UVES pipeline) did not affect the results,
and that any over-abundances in s-process elements did not affect the
derivation of the metallicity (i.e., [M/H] $\equiv$ [Fe/H]). The surface
temperature ($T_{\mathrm{eff}}\sim$ 5000 K), surface gravity (log $g\sim 3.5$)
and metallicity ([Fe/H]$\sim-0.4$) are in line with those expected based on
the similarity with the G8 IV-V star HD 24616 (as shown by Miszalski et al.,
2012)444As HD 24616 is characterised by an effective temperature $T_{\rm
eff}\approx 4900$ K, $\log g\approx 3.3$, and metallicity [Fe/H]$\approx-0.2$
(Steinmetz et al., 2020).. The Gaia Apsis surface gravity ($\log g\approx
4.8$) is significantly larger than the determined surface gravity, but results
in an appreciably poorer fit to the spectrum especially at the extremely low
Apsis metallicity, [Fe/H]=$-2.73$. Here, it is also worth noting that our
decision to fix the surface gravity does not dramatically impact the other
parameters derived (beyond the mass) nor our subsequently derived chemical
abundances. More generally, it is important to recognise that given the
relatively large rotational velocity of the star, resulting in broad lines and
a lot of blending, and the relatively low signal-to-noise ratio of the co-
added spectrum (S/N$\sim$20), it is not possible to establish any parameters
with high precision. In particular, the gravity is very poorly constrained.

Using the best model, we also tried to estimate the abundance of several
chemical elements using the equivalent width method, including s-process
elements. Unfortunately, the abundance of only 2 elements could be derived
with reasonable precision (see Fig. 8): Ca and Ba, while it appears that there
is no need to enhance the carbon abundance to get a good fit (nor the
chromium, nickel or titanium, although these elemental abundances are less
well constrained). The barium abundance is based on the barium lines at 6141
and 6496 Å. Using the line at 4554 Å results in a lower abundance of
$\approx$0.8. However, it is well known that the strong line at 4554 Å is
severely affected by NLTE effects, which are not taken into account in our
modeling, and this value should thus be considered a lower limit. With +1 dex
of Ba overabundance, the central star of A 70 is classed as a strong barium
star – these can have periods of up to 104 days but are generally found to
cluster in the range of 1000–4000 days (Jorissen et al., 2019b), compatible
with the range of periods previously inferred from the UVES RVs.

The final parameters and abundances as well as their uncertainties are listed
in Tab. 3.

Figure 8: Portion of the co-added UVES spectrum of the central star of A 70 showing the region around the 6141Å Ba line (in grey). The blue line shows the best fit model with [Fe/H]=$-0.5$, [Ba/Fe]=$+1.2$ and [Ca/Fe]=$+0.5$, while the red line shows the same model without enhanced barium or calcium.  Table 3: Most probable range of the parameters of the bright secondary of A 70. Parameters | A 70
---|---
$T_{\mathrm{eff}}$ (K) | 5000$\pm$250
log $g$ | 3.5 (fixed)
$v$ sin $i$ ( km s-1) | 37 $\pm$ 5
[Fe/H] | $-0.5\pm 0.2$
Microturbulence ( km s-1) | 1.0 $\pm$ 0.5
[Ba/Fe] | 1.2 $\pm$ 0.2
[Ca/Fe] | 0.5 $\pm$ 0.1
[C/H] |  $\approx$0
[Ti/Fe] |  $\approx$0
[Cr/Fe] |  $\approx$0
[Ni/Fe] |  $\approx$0

##  4 Distance, mass and radius of the companion

Assuming the companion rotates in the same plane as the binary orbit and the
binary orbit is coincident with the waist of the nebula (i.e. $i$=30°; Tyndall
et al., 2013), the measured rotation period of 2.06 d and $v$ sin $i$ of
$\approx$37 km s-1 imply a radius of $\sim$3 R⊙. For a surface gravity of log
$g$=3.5, this radius would imply a mass for the companion of $\approx$1 M⊙.
These parameters, as well as the previously derived temperature, are
consistent with a sub-giant on the verge of joining the first giant branch
(Pols et al., 1998). Irrespective of the precise inclination of its rotation
axis, the radius of the companion star must be $\gtrsim$1.5 R⊙ (for the
derived rotation period and stellar parameters). This is still a strong
indication that the star has evolved off the main sequence unless it has been
inflated following the mass transfer (as seen in some post-common-envelope
central stars; Jones et al., 2015, 2022) which resulted in its contamination
with the products of AGB nucleosynthesis from the nebular progenitor (Ca and
Ba).

Using the best-fitting synthetic spectrum from section 3 and assuming an
extinction of E(B$-$V)=0.05 (Frew et al., 2016), we calculated synthetic
magnitudes in the three Gaia bands (see Tab. 1) allowing the distance to be a
free parameter. The best-fitting distance was 4.6 kpc, in reasonable agreement
with both the Gaia geometric parallax distance (3.6${}^{+1.5}_{-1.0}$ kpc) and
the nebular surface-brightness-radius relationship distance (6.0$\pm$2.1 kpc)
of Frew et al. (2016). The Gaia Apsis distance is much closer at $\approx$1.5
kpc, however several other parameters derived by the Gaia Apsis also appear to
be strongly incongruent. While the Apsis extinction agrees well with that
measured from the nebula by Frew et al. (2016), the Apsis radius is not
compatible with the Apsis distance and the observed Gaia magnitudes (with the
radius required to match the observed photometry being approximately double
the Apsis radius). Similarly, the Apsis metallicity of [Fe/H]$\approx-$2.7 is
entirely incompatible with the spectra presented in section 3. The strong
disagreement is hard to explain, however GSP-Phot is known to provide
unreliable (and underestimated) distances for sources further than 2 kpc
(Fouesneau et al., 2022), as well as to systematically underestimate
metallicities (Andrae et al., 2022). Ultimately, the discrepancies could be
related to the GSP-Phot module’s fixing of the available parameter space using
isochrones, which may not be entirely compatible with a star that has
experienced appreciable accretion as is the case here (Andrae et al., 2022).

The very strong barium enhancement we derive implies that the secondary has
accreted a large amount of polluted material from the former AGB. Indeed,
according to Stancliffe (2021), one would require a dilution factor of the
order of 0.25 to account for such an s-process excess. If the secondary is
about 1 M⊙, this means it must have accreted about 0.25 M⊙ from its AGB
companion. This in turn means that, assuming the primary was originally 2-3 M⊙
(typical for AGB progenitors of the barium star’s companion; Karinkuzhi et
al., 2018), the mass transfer efficiency must have been of the order of 10-20%
– compatible with wind Roche lobe overflow (Abate et al., 2013).

##  5 Primary mass

As the primary is not visible in our spectra, we cannot place any further
constraints on its temperature or radius beyond those already presented by
Miszalski et al. (2012). However, the radial velocities presented in section 3
can be used to probe its mass. The mass function of an eccentric binary is
defined as:

$f(\mathrm{M}_{1})=\frac{\mathrm{M}_{1}^{3}\mathrm{sin}^{3}i}{(\mathrm{M}_{1}+\mathrm{M}_{2})^{2}}=\frac{\mathrm{P}~{}\mathrm{K}_{2}^{3}}{2\pi
G}(1-e^{2})^{\frac{3}{2}}$ (1)

where M1 and M2 are the masses of the central star and companion,
respectively. P is the orbital period, K2 is the radial velocity semi-
amplitude of the companion, and $e$ is the orbital eccentricity.

Assuming that the systemic velocity of the binary is, indeed, the same as the
systemic velocity of the surrounding nebula, we can place limits on the
majority of these parameters. The orbital period must be in the range 2 yr
$\lesssim$ P $\lesssim$ 8 yr (shorter periods can be ruled out by the lack of
variability in the 2011 UVES RVs), the semi-amplitude of the companion must be
K${}_{2}\gtrsim$10 km s-1 (the difference between the nebular systemic
velocity and the UVES RVs), and the orbital eccentricity must be e$\gtrsim$0.3
(lower eccentricities would also not reproduce the observed lack of
variability in the 2011 UVES RVs). This leads to a mass function of 0.05–0.1
M⊙ which, accounting for the derived mass of the companion (M${}_{2}\approx$ 1
M⊙) and assuming the binary orbital plane corresponds to the waist of the
nebula, implies an problematically large mass for the primary (close to, or
even exceeding, the Chandrasekhar mass). Stellar evolution and nucleosynthesis
models indicate that the pollutors of Ba stars should almost all have had
initial masses M$\sim$ 2–3 M⊙ (Karinkuzhi et al., 2018). This places a rough
upper limit of 0.75 M⊙ for the remnant based on the initial-final mass
relation (but with large uncertainties; Cummings et al., 2018).

While misalignment between rotation axis and orbital planes is not impossible,
it would be particularly strange that the rotation axis of the G-type star is
coincident with the nebular symmetry axis while the binary orbital plane is
not. Given that the mass loss should be focused in the orbital plane one would
expect that for a misaligned system it would be the binary plane that is
aligned with the nebular symmetry axis rather than the stellar rotation axis.
However, such a misalignment between binary and nebula is not unheard of -
there is some suggestion that the same is true for another barium star
planetary nebula LoTr 5 (Jones et al., 2017; Aller et al., 2018, although,
given the uncertainties, there are still doubts as to whether this is proven),
perhaps indicating that the expected alignment between binary plane and
nebular symmetry axis (Hillwig et al., 2016, and references therein) does not
hold for such long period systems.

It is possible that the G-type star observed at the projected centre of A 70
is simply a chance alignment and is not related to the nebula. Indeed, other
curious chance alignments have been discovered in other PNe, for example, SuWt
2 (Jones & Boffin, 2017b) and M 2-3 (Boffin et al., 2018) both of which have
bright main-sequence binaries close to their projected centres. However, a
chance alignment is highly unlikely in the case of A 70 given that the star is
contaminated in s-process elements, indicating that it must be in a binary
with a more evolved companion from which it accreted s-process rich material.
Similarly, while the lack of observed radial velocity variability in the
G-type star may support the chance alignment hypothesis, a similar conclusion
was drawn for NGC 1514 before more extensive observations revealed that the
orbit was simply too long and eccentric to have been detected in a single
observing season (Jones et al., 2017).

Alternatively, the systemic velocity of the binary may not be the same as the
systemic velocity of the PN. Assuming a moderate eccentricity of e$\approx$0.3
and that the orbital inclination is the same as the nebular inclination, the
radial velocity semi-amplitude should be 4–5 km s-1 in order for the primary
mass to be approximately 0.6 M⊙ (a relatively canonical value for such white
dwarfs). Based on the UVES RVs, this would mean that the binary systemic
velocity should be approximately $-$79 km s-1 – some 6 km s-1 (although only
still less than two uncertainties) removed from the nebular systemic velocity
of $-$73$\pm$4 km s-1 (Tyndall et al., 2013).

A difference in systemic velocity between binary and nebula could, in
principle, be due to a “kick” which also helped impart eccentricity into the
Ba star orbit (Izzard et al., 2010). In this case, the central star might be
expected to be slightly offset from the centre of the nebula but, given the
clumpy nature of the nebular ring, this is hard to assess using existing
imagery (see, e.g., Fig. 8 of Tyndall et al., 2013). For the same reason, it
is also possible that the uncertainty on the systemic velocity of the nebula
has been underestimated. There is precedent for such a conclusion with Jones
et al. (2010) showing that, unless a detailed spatio-kinematical modelling is
performed, the brightness variations across a PN can skew the measurement of
its systemic velocity by several km s-1. This hypothesis is particularly
plausible for A 70 given the clear brightness variations across the nebula
where, for example, the slit number 1 of Tyndall et al. (2013) shows bright
emission from the ring centred at a heliocentric velocity of approximately
$-$80 km s-1 but significantly fainter emission extending redwards to $-$30 km
s-1, while their slit position 3 shows the converse with bright emission
centred at approximately $-$60 km s-1 and fainter emission extending bluewards
up to $-$120 km s-1.

Ultimately, without further observations which definitively constrain the
binary period and radial velocity amplitude, it is not possible to constrain
the mass of the primary, nor the possibility of misalignment between the
nebular symmetry axis and binary orbital plane.

##  6 Conclusions

The rotation period of the G-type companion to the central star of A 70 is
confirmed to be roughly 2 d based on photometric variability, while the
spectroscopically derived parameters ($T_{\mathrm{eff}}$, log $g$ and $v$ sin
$i$) are consistent with a G8IV-type star rotating with the same period at an
inclination equal to that found for the nebula by Tyndall et al. (2013). The
G-type star is chromospherically active, presenting Ca H and K in emission,
consistent with having been spun-up by the accretion event which led to its
enhancement in barium. A broad H$\alpha$ emission feature is visible in the
spectrum of the G-type star, the origin of which is unclear but which has been
observed in similar stars (Tyndall et al., 2013). The spots due to the
chromospheric activity are found to change phase (i.e., location on the
stellar surface) in 2011, associated with a period of negligible activity,
followed by a period of at least four years of stable variability.

The radial velocity measurements of the G-type star, however, do not appear to
be entirely consistent with the previously derived properties of the nebula –
either indicating a higher inclination or a slightly different systemic
velocity. Both discrepancies could potentially be resolved by continued radial
velocity monitoring in order to conclusively constrain the orbital period and
semi-amplitude. For now, we conclude that there is insufficient evidence to
support the hypothesis that the nebula and binary are misaligned (as has been
suggested in another long-period central star, LoTr 5; Aller et al., 2018),
and we encourage further detailed study of both the central star and nebula of
A 70. In particular, if the highly speculative $\sim$8 year period fit shown
in Fig. 5 is representative, the next radial velocity maximum will occur in
late-2029 – as such this will be a critical time to obtain further radial
velocity measurements.

## Acknowledgments

The authors would like to thank Howard E. Bond and Robin Ciardullo for making
the data in their Research Note publicly available, and the anonymous referee
for their insightful report.

DJ acknowledges support from the Erasmus+ programme of the European Union
under grant number 2020-1-CZ01-KA203-078200. DJ also acknowledges support
under grant P/308614 financed by funds transferred from the Spanish Ministry
of Science, Innovation and Universities, charged to the General State Budgets
and with funds transferred from the General Budgets of the Autonomous
Community of the Canary Islands by the Ministry of Economy, Industry, Trade
and Knowledge. JM and AJB acknowledge the support of STFC in the form of
studentships. This research was supported by the Erasmus+ programme of the
European Union under grant number 2017-1-CZ01-KA203-035562.

Based on observations made with ESO Telescopes at the La Silla Paranal
Observatory under programme IDs 083.D-0654, 085.D-0629, 087.D-0174,
087.D-0205, 091.D-0475, 0103.D-0186 and 0104.D-0326. This work makes use of
observations from the LCOGT network, and observations made at the South
African Astronomical Observatory (SAAO). The CRTS survey is supported by the
U.S. National Science Foundation under grants AST-0909182 and AST-1313422.
Based on observations obtained with the Samuel Oschin 48-inch Telescope at the
Palomar Observatory as part of the Zwicky Transient Facility project. ZTF is
supported by the National Science Foundation under Grant No. AST-1440341 and a
collaboration including Caltech, IPAC, the Weizmann Institute for Science, the
Oskar Klein Center at Stockholm University, the University of Maryland, the
University of Washington, Deutsches Elektronen-Synchrotron and Humboldt
University, Los Alamos National Laboratories, the TANGO Consortium of Taiwan,
the University of Wisconsin at Milwaukee, and Lawrence Berkeley National
Laboratories. Operations are conducted by COO, IPAC, and UW

This research has made use of NASA’s Astrophysics Data System Bibliographic
Services; the SIMBAD database, operated at CDS, Strasbourg, France; the VizieR
catalogue access tool, CDS, Strasbourg, France; Astropy, a community-developed
core Python package for Astronomy (Astropy Collaboration et al., 2013; Price-
Whelan et al., 2018); PyAstronomy; SciPy (Jones et al., 2001); NumPy (Van Der
Walt et al., 2011).

## Data availability

Individual photometric and radial velocity measurements are made available
either in the paper or as online material via CDS. The raw ESO, LCOGT and LT
data are available from the respective online archives, while the raw SAAO are
available upon reasonable request to the corresponding author. Survey data are
available as part of their respective online data releases.

## References

  * Abate et al. (2013) Abate C., Pols O. R., Izzard R. G., Mohamed S. S., de Mink S. E., 2013, A&A, 552, A26
  * Aller et al. (2018) Aller A., Lillo-Box J., Vučković M., Van Winckel H., Jones D., Montesinos B., Zorotovic M., Miranda L. F., 2018, MNRAS, 476, 1140
  * Andrae et al. (2022) Andrae R., et al., 2022, arXiv e-prints, p. arXiv:2206.06138
  * Appenzeller et al. (1998) Appenzeller I., et al., 1998, The Messenger, 94, 1
  * Astropy Collaboration et al. (2013) Astropy Collaboration et al., 2013, A&A, 558, A33
  * Bailer-Jones et al. (2021) Bailer-Jones C. A. L., Rybizki J., Fouesneau M., Demleitner M., Andrae R., 2021, AJ, 161, 147
  * Blanco-Cuaresma (2019) Blanco-Cuaresma S., 2019, MNRAS, 486, 2075
  * Blanco-Cuaresma et al. (2014) Blanco-Cuaresma S., Soubiran C., Heiter U., Jofré P., 2014, A&A, 569, A111
  * Boffin & Jones (2019) Boffin H. M. J., Jones D., 2019, The Importance of Binaries in the Formation and Evolution of Planetary Nebulae. Springer Nature, doi:10.1007/978-3-030-25059-1
  * Boffin & Jorissen (1988) Boffin H. M. J., Jorissen A., 1988, A&A, 205, 155
  * Boffin et al. (2018) Boffin H. M. J., et al., 2018, A&A, 619, A84
  * Bond & Ciardullo (2018) Bond H. E., Ciardullo R., 2018, Research Notes of the American Astronomical Society, 2, 143
  * Bond et al. (2003) Bond H. E., Pollacco D. L., Webbink R. F., 2003, AJ, 125, 260
  * Brown et al. (2013) Brown T. M., et al., 2013, PASP, 125, 1031
  * Buzzoni et al. (1984) Buzzoni B., et al., 1984, The Messenger, 38, 9
  * Ciardullo et al. (1999) Ciardullo R., Bond H. E., Sipior M. S., Fullton L. K., Zhang C.-Y., Schaefer K. G., 1999, AJ, 118, 488
  * Creevey et al. (2022) Creevey O. L., et al., 2022, arXiv e-prints, p. arXiv:2206.05864
  * Cummings et al. (2018) Cummings J. D., Kalirai J. S., Tremblay P. E., Ramirez-Ruiz E., Choi J., 2018, ApJ, 866, 21
  * Dekker et al. (2000) Dekker H., D’Odorico S., Kaufer A., Delabre B., Kotzlowski H., 2000, in M. Iye & A. F. Moorwood ed., Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series Vol. 4008, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series. pp 534–545
  * Drake et al. (2009) Drake A. J., et al., 2009, ApJ, 696, 870
  * Escorza et al. (2019) Escorza A., et al., 2019, A&A, 626, A128
  * Fouesneau et al. (2022) Fouesneau M., et al., 2022, arXiv e-prints, p. arXiv:2206.05992
  * Frew et al. (2016) Frew D. J., Parker Q. A., Bojičić I. S., 2016, MNRAS, 455, 1459
  * Gaia Collaboration (2022) Gaia Collaboration 2022, A&A, p. in press
  * Gray & Corbally (1994) Gray R. O., Corbally C. J., 1994, AJ, 107, 742
  * Gustafsson et al. (2008) Gustafsson B., Edvardsson B., Eriksson K., Jørgensen U. G., Nordlund Å., Plez B., 2008, A&A, 486, 951
  * Hillwig et al. (2016) Hillwig T. C., Jones D., Marco O. D., Bond H. E., Margheim S., Frew D., 2016, The Astrophysical Journal, 832, 125
  * Hook et al. (2004) Hook I. M., Jørgensen I., Allington-Smith J. R., Davies R. L., Metcalfe N., Murowinski R. G., Crampton D., 2004, PASP, 116, 425
  * Horne et al. (1986) Horne K., Wade R. A., Szkody P., 1986, MNRAS, 219, 791
  * Izzard et al. (2010) Izzard R. G., Dermine T., Church R. P., 2010, A&A, 523, A10
  * Jacob et al. (2013) Jacob R., Schönberner D., Steffen M., 2013, A&A, 558, A78
  * Jasniewicz et al. (1994) Jasniewicz G., Acker A., Mauron N., Duquennoy A., Cuypers J., 1994, A&A, 286, 211
  * Jayasinghe et al. (2019) Jayasinghe T., et al., 2019, MNRAS, 485, 961
  * Jones (2018) Jones D., 2018, in Boffin H. M. J., Beccari G., eds, The Importance of Binaries in Stellar Evolution. (arXiv:astro-ph/0108426)
  * Jones & Boffin (2017a) Jones D., Boffin H. M. J., 2017a, Nature Astronomy, 1, 0117
  * Jones & Boffin (2017b) Jones D., Boffin H. M. J., 2017b, MNRAS, 466, 2034
  * Jones et al. (2001) Jones E., Oliphant T., Peterson P., et al., 2001, SciPy: Open source scientific tools for Python, http://www.scipy.org/
  * Jones et al. (2010) Jones D., Lloyd M., Mitchell D. L., Pollacco D. L., O’Brien T. J., Vaytet N. M. H., 2010, MNRAS, 401, 405
  * Jones et al. (2015) Jones D., Boffin H. M. J., Rodríguez-Gil P., Wesson R., Corradi R. L. M., Miszalski B., Mohamed S., 2015, A&A, 580, A19
  * Jones et al. (2017) Jones D., Van Winckel H., Aller A., Exter K., De Marco O., 2017, A&A, 600, L9
  * Jones et al. (2022) Jones D., et al., 2022, MNRAS, 510, 3102
  * Jorissen et al. (1998) Jorissen A., Van Eck S., Mayor M., Udry S., 1998, A&A, 332, 877
  * Jorissen et al. (2019a) Jorissen A., Boffin H. M. J., Karinkuzhi D., Van Eck S., Escorza A., Shetye S., Van Winckel H., 2019a, A&A, 626, A127
  * Jorissen et al. (2019b) Jorissen A., Boffin H. M. J., Karinkuzhi D., Van Eck S., Escorza A., Shetye S., Van Winckel H., 2019b, A&A, 626, A127
  * Karinkuzhi et al. (2018) Karinkuzhi D., et al., 2018, A&A, 618, A32
  * Löbling et al. (2019) Löbling L., Boffin H. M. J., Jones D., 2019, A&A, 624, A1
  * Masci et al. (2019) Masci F. J., et al., 2019, PASP, 131, 018003
  * Matrozis et al. (2017) Matrozis E., Abate C., Stancliffe R. J., 2017, A&A, 606, A137
  * McClure et al. (1980) McClure R. D., Fletcher J. M., Nemec J. M., 1980, ApJL, 238, L35
  * Miszalski et al. (2012) Miszalski B., Boffin H. M. J., Frew D. J., Acker A., Köppen J., Moffat A. F. J., Parker Q. A., 2012, MNRAS, 419, 39
  * Miszalski et al. (2013) Miszalski B., et al., 2013, MNRAS, 436, 3068
  * Mohamed & Podsiadlowski (2012) Mohamed S., Podsiadlowski P., 2012, Baltic Astronomy, 21, 88
  * Munday et al. (2020) Munday J., et al., 2020, MNRAS, 498, 6005
  * Nagae et al. (2004) Nagae T., Oka K., Matsuda T., Fujiwara H., Hachisu I., Boffin H. M. J., 2004, A&A, 419, 335
  * Oláh et al. (2013) Oláh K., Moór A., Strassmeier K. G., Borkovits T., Granzer T., 2013, Astronomische Nachrichten, 334, 625
  * Pols et al. (1998) Pols O. R., Schröder K.-P., Hurley J. R., Tout C. A., Eggleton P. P., 1998, MNRAS, 298, 525
  * Price-Whelan et al. (2018) Price-Whelan A. M., et al., 2018, AJ, 156, 123
  * Ricker et al. (2015) Ricker G. R., et al., 2015, Journal of Astronomical Telescopes, Instruments, and Systems, 1, 014003
  * Shappee et al. (2014) Shappee B. J., et al., 2014, ApJ, 788, 48
  * Snodgrass et al. (2008) Snodgrass C., Saviane I., Monaco L., Sinclaire P., 2008, The Messenger, 132, 18
  * Stancliffe (2021) Stancliffe R. J., 2021, MNRAS, 505, 5554
  * Stassun et al. (2019) Stassun K. G., et al., 2019, AJ, 158, 138
  * Steinmetz et al. (2020) Steinmetz M., et al., 2020, AJ, 160, 83
  * Theuns & Jorissen (1993) Theuns T., Jorissen A., 1993, MNRAS, 265, 946
  * Theuns et al. (1996) Theuns T., Boffin H. M. J., Jorissen A., 1996, MNRAS, 280, 1264
  * Tonry et al. (2018) Tonry J. L., et al., 2018, PASP, 130, 064505
  * Tyndall et al. (2013) Tyndall A. A., et al., 2013, MNRAS, 436, 2082
  * Van Der Walt et al. (2011) Van Der Walt S., Colbert S. C., Varoquaux G., 2011, Computing in Science and Engineering, 13, 2

##  Appendix A Survey photometry

###  A.1 Catalina Real-Time Transient Survey

Figure 9: The CRTS photometry of the central star of A 70.

The central star of A 70 was also observed 212 times by the Catalina Real-Time
Transient Survey (CRTS) between October 2005 and November 2013 (Drake et al.,
2009). The combined data originating from all three telescopes used in the
survey does not present with a single consistent period, with observations
taken with the larger (and higher spatial resolution) 1.5m Mount Lemmon Survey
(MLS) and 0.7m Catalina Sky Survey (CSS) telescopes presenting with a large
($>$1 magnitude) and aperiodic scatter. However, the 32 points observed using
the 0.5m Siding Springs Survey (SSS) telescope do present variability
consistent with that derived from the photometry presented above (period
$\sim$ 2 days, semi-amplitude $\sim$ 0.1 magnitudes). Given that the
photometric extraction performed by the CRTS is not optimised for the central
stars of planetary nebulae, and that the nebular emission of A 70 is
appreciable in the V-band employed by the survey, we conclude that the high
amplitude, aperiodic variability found in the photometry from MLS and CSS are
most likely spurious and related to variations in seeing (and thus variable
nebular contribution within the aperture; Jones et al., 2015). This offers a
natural explanation as to why the data from the lower spatial resolution SSS
telescope is more consistent with the 2-day period derived above, as the
larger pixels (pixel scale $\gtrsim$ seeing) would be less susceptible to
variable nebular contamination in the aperture. A sinusoidal fit to the early
SSS data is shown in the left panel of Fig. 10. In spite of the small number
of measurements, the period and amplitude are seemingly consistent with the
variability shown in Fig. 1 but with a different phasing to the 2010 data.
This is perhaps indicative of another period of negligible variability shortly
before the first observations by Bond & Ciardullo (2018). However, given the
quality of the fit and small number of highly uncertain data points, we stress
that this conclusion is at best speculative. Fitting a sinusoid to all the
CRTS data after mid-2011 (obtained with all three telescopes: MLS, CSS, and
SSS), and assuming a 2.061 d period, reveals a best fit phasing more or less
consistent with the 2011-2015 ephemeris albeit with a very large $\chi^{2}$
(right panel, Fig. 10).

Figure 10: Left: The phase-folded CRTS photometry of the central star of A 70
obtained with the SSS telescope prior to 2009. Right: The phase-folded CRTS
photometry obtained after mid-2011. All data are folded on the same ephemeris
as in Fig. 1, with a sinusoidal fit shown underlaid in red.

###  A.2 Asteroid Terrestrial-impact Last Alert System

The central star of A 70 was observed 1643 times in the $o$-band and 478 in
the $c$-band by the Asteroid Terrestrial-impact Last Alert System (ATLAS;
Tonry et al., 2018) between July 2015 and July 2022. The data are rather
imprecise for our purposes, with a median uncertainty in both bands of
approximately 0.6 magnitudes (see Fig. 11). Perhaps unsurprisingly, there are
no significant peaks in the periodograms and the data do not phase well on the
previously measured photometric period.

Figure 11: The ATLAS photometry of the central star of A 70. Note that data
points with uncertainties larger than 1 mag have been discarded.

###  A.3 All-Sky Automated Survey for Supernovae

The central star of A 70 was observed by the All-Sky Automated Survey for
Supernovae (ASAS-SN; Shappee et al., 2014; Jayasinghe et al., 2019) 338 times
in the $V$-band between November 2012 and November 2018, and a further 517
times in the $g$-band between October 2017 and July 2022. The magnitudes from
ASAS-SN are significantly brighter ($\sim$4 magnitudes) than the values from
both CRTS and ATLAS (see Fig. 12), strongly indicative that the photometry is
heavily contaminated by the surrounding nebula. Furthermore, the $g$-band data
separates in to two separate regimes, similarly indicative of an issue with
the photometry. Sinusoid fitting to the $V$-band data results in a phasing
consistent with the 2011-2015 ephemeris, but with an amplitude that is too low
(only 0.01 mag) and a very large reduced $\chi^{2}$ ($>$10). As one might
expect from their bi-modal nature, the $g$-band data do not phase well with
any ephemeris. Therefore, just as for the ATLAS data, the ASAS-SN data are not
useful in constraining the variability of the central star of A 70.

Figure 12: The ASAS-SN photometry of the central star of A 70.

###  A.4 Zwicky Transient Facility

The central star of A 70 was observed 45 times in the $g$-band and 90 times in
the $r$-band by the Zwicky Transient Facility (ZTF; Masci et al., 2019)
between May 2018 and December 2019. Very much like the CRTS data, the ZTF data
shows a relatively large ($\sim$0.4 magnitudes) and seemingly aperiodic
scatter in both bands (see Fig. 13). Both show peaks in their periodograms at
periods close to 2.06 d (although neither is significant). Simultaneous
sinusoid fitting of both $g$\- and $r$-bands, as shown in Fig. 14, results in
a phasing which is not consistent with the 2011-2015 and/or 2022 data
presented in Fig. 1 (but rather more consistent with the 2010 epoch). While
not conclusive (given the uncertainty on the photometric period), this is
seemingly indicative that the variability was not stable over the intervening
time period.

Figure 13: ZTF photometry of the central star of A 70.  Figure 14: Phase-
folded ZTF photometry of the central star of A 70. The sinusoidal fit to the
data (assuming the same time of minimum but different amplitudes) are shown
underlaid in grey.

###  A.5 Transiting Exoplanet Survey Satellite

The field containing A 70 was observed by the Transiting Exoplanet Survey
Satellite (TESS; Ricker et al., 2015) in Sector 54 (9 July 2022 – 5 August
2022), however at the time of submission that data has yet to be released. In
any case, at an approximate TESS magnitude of 17.1 (as listed in the TESS
input catalog; Stassun et al., 2019), it is highly unlikely that TESS will
return useful photometry of the central star system. Furthermore, the TESS
pixel size of 21″ is comparable to the angular size of the nebula, so the TESS
photometry of the central star is likely be heavily contaminated by nebular
emission.


Implications of exclusive $J/\psi$ photoproduction in a tamed collinear
factorisation approach to NLO

Chris A. Flett1$\star$, Alan D. Martin2, Misha G. Ryskin and Thomas Teubner3

1 Department of Physics, University of Jyväskylä, P.O. Box 35, 40014
University of Jyväskylä, Finland and Helsinki Institute of Physics, P.O. Box
64, 00014 University of Helsinki, Finland

2 Institute for Particle Physics Phenomenology, Durham University, Durham, DH1
3LE, U.K.

3 Department of Mathematical Sciences, University of Liverpool, Liverpool, L69
3BX, U.K.

*Speaker<EMAIL_ADDRESS>
51st International Symposium on Multiparticle Dynamics (ISMD2022)
Pitlochry, Scottish Highlands, 1-5 August 2022 10.21468/SciPostPhysProc.?

## Abstract

We discuss exclusive $J/\psi$ photoproduction, initially in conventional
collinear factorisation at NLO and then subsequently in a refined approach
with a programme of low $x$ resummation and implementation of a crucial low
$Q_{0}$ subtraction included. We compare and contrast predictions in both
frameworks and remark about the possibility to constrain and ultimately
determine the low $x$ and low scale gluon PDF, emphasising the significance of
this for future global PDF analyses.

##  1 Introduction

The multitude of inclusive Deep Inelastic Scattering (DIS) experimental data
from colliders past and present do not well constrain the global parton
distribution function (PDF) analyses at small momentum fractions $x$ and
factorisation scales $\mu^{2}$.

Figure 1: Left: Comparison of the absolute values of the NLO gluon PDFs
obtained in various global analyses [1, 2, 3] at the scale $\mu$ = 1.55 GeV.
Right: Comparison of ratios of the NLO gluon PDFs normalised to the CT18A set,
emphasising the large uncertainties present at low- and high-$x$.

Fig. 1 (left) shows, for example, the gluon PDFs obtained in the MSHT20 [1],
NNPDF3.1 [2] and CT18A [3] NLO analyses as a function of the parton momentum
fraction $x$ at a fixed low scale $\mu$ = 1.55 GeV. At this scale, it is clear
that as $x$ gets smaller, the uncertainties attributed to the gluon PDF become
larger, and the extractions are consistent with an arguably unphysical
decreasing gluon density as well as a vanishing one. This is primarily due to
the lack of experimental data constraints in the global fits in this small $x$
and small $\mu$ domain.

As is evident in Fig. 1 (right), there also exists a sizeable uncertainty for
larger $x>0.1$ in the region where the valence quark PDFs dominate. This is
not visible in the log-linear display of Fig. 1 (left) but becomes so if we
take ratios of the sets normalised to e.g. the CT18A set, as shown in Fig. 1
(right). This large uncertainty has implications for quark PDF flavour
separation and Beyond-The-Standard-Model searches but, in these proceedings,
we will be concerned with only the small $x$ and small $\mu$ domain. In
particular, we study exclusive $J/\psi$ photoproduction, $\gamma p\rightarrow
J/\psi p$, and explain why the data for this process from HERA and LHC can be
used in the global PDF analyses to determine the gluon PDF down to $x\sim
3\times 10^{-6}$ and factorisation scales $\mu^{2}=O(M_{J/\psi}^{2}/4)$, where
$M_{J/\psi}$ is the mass of the $J/\psi$. As we will describe, this is
complicated theoretically as, due to the off-forward, or skewed, parton
kinematics, the process is sensitive to the generalised parton distribution
function (GPD) rather than the usual forward one and, moreover, within the
perturbative QCD (pQCD) approach to next-to-leading order (NLO), there is a
large factorisation scale dependence exhibited by the $\gamma p\rightarrow
J/\psi p$ hard matrix element. Nevertheless, we will show that this scale
dependence can be largely alleviated by supplementing standard collinear
factorisation with resummation of a class of large double logarithms [4] and
implementation of a crucial low $Q_{0}$ subtraction [5], where $Q_{0}$ is the
PDF input scale. We refer to this as a ‘tamed’ collinear factorisation
approach in the following. Furthermore, the GPD can be related to the PDF via
the so-called Shuvaev integral transform [6], which admits an accuracy of
$O(x)$ at NLO in pQCD and so, on a practical level, becomes a reliable means
to ascertain the behaviour of skewed partons at small $x$.

After describing our general framework in Section 2, we compare and contrast
exclusive $J/\psi$ photoproduction to NLO in the conventional and tamed
collinear factorisation approaches in Section 3. Then, in Section 4, we
discuss the implications of this tamed approach for the $\overline{\text{MS}}$
gluon PDF and end, in Section 5, with concluding remarks.

##  2 Framework

The exclusive $J/\psi$ photoproduction process, $\gamma p\rightarrow J/\psi
p$, proceeds via the fluctuation of a quasi-real photon into a $c\bar{c}$
pair, which hadronises into the final state $J/\psi$ meson through a two-
parton exchange mechanism, see Fig. 2. In the Bjorken limit at leading twist,
the net momentum fraction along the lightcone direction $P^{+}$ is $2\xi$,
where $\xi$ is the skewedness parameter, and so the parton momentum fractions
can be parametrised as $X+\xi$ and $X-\xi$. The $c\bar{c}$ to $J/\psi$
transition is described within LO Non-Relativistic QCD (NRQCD), with the
relative quark velocity in the outgoing $c\bar{c}$ pair set to zero and with
$M_{J/\psi}=2m_{c}$, where $m_{c}$ is the mass of the charm quark. The first
relativistic correction has been studied and shown to be small once the LO
non-perturbative matrix element is normalised to the leptonic decay width of
the $J/\psi$ [7].

Figure 2: Examples of a LO contribution (left) and NLO quark contribution
(right) to the exclusive process $\gamma p\to J/\psi p$. Here, the momentum
$P\equiv(p+p^{\prime})/2$ and $l$ is the loop momentum. The momentum fractions
of the left and right input partons are $x=X+\xi$ and $x^{\prime}=X-\xi$
respectively; for the gluons coupled directly to the on-shell heavy-quark
pair, we have $x^{\prime}\ll x$ and so $x\simeq 2\xi$.

The amplitude for exclusive $J/\psi$ photoproduction in the collinear
factorisation framework at NLO can therefore be written as

$A(\xi)\sim\sqrt{\langle
O_{1}\rangle_{V}}\int_{-1}^{1}\frac{\text{d}X}{X}\left[C_{g}(X,\xi,\mu_{R},\mu_{F})F_{g}(X,\xi,\mu_{F})+C_{q}(X,\xi,\mu_{R},\mu_{F})F_{q}(X,\xi,\mu_{F})\right],$
(1)

where $\langle O_{1}\rangle_{V}$ is the LO non-perturbative NRQCD matrix
element, $C_{g,q}$ are the NLO gluon and quark coefficient functions [8] and
$F_{g,q}$ are the NLO gluon and quark singlet GPDs. The renormalisation and
factorisation scales are denoted by $\mu_{R}$ and $\mu_{F}$, respectively. At
LO, there is only the gluonic term while at NLO there is a contribution from
the quark sector, too. The GPDs at sufficiently small values of $x,\xi$ are
reliably obtained via the Shuvaev transform [6], with the forward PDFs taken
from the LHAPDF6 interface [9].

##  3 Exclusive $J/\psi$ photoproduction in collinear factorisation

In the conventional collinear factorisation approach to NLO in the
$\overline{\text{MS}}$ scheme, the exclusive $J/\psi$ photoproduction suffers
from a large factorisation scale dependence and the NLO result is larger and
of opposite sign to the LO one, suggesting a poor perturbative stability in
$\alpha_{s}$. This is illustrated in Fig. 3, which shows Im $A/W^{2}$ as a
function of $W$, the $\gamma p$ centre-of-mass energy. We use state-of-the-art
CT18ANLO input partons and show the LO and NLO scale variation bands for
$\mu_{F}^{2}=\mu_{R}^{2}=Q_{0}^{2},m_{c}^{2},2m_{c}^{2}$. It is clear that the
factorisation scale dependence worsens at NLO and that in general this process
exhibits a strong sensitivity with respect to scale variations. However,
inspection of the amplitude in the high-energy limit, $W^{2}\gg
M_{J/\psi}^{2}$, shows that logarithmically enhanced contributions at small
$x$ can spoil the perturbative convergence.

Figure 3: LO (black) and LO+NLO (blue) contributions to Im $A/W^{2}$ generated
using CT18ANLO [3] global partons with
$\mu_{F}^{2}=\mu_{R}^{2}=Q_{0}^{2},m_{c}^{2}$ and $2m_{c}^{2}$ (from bottom to
top). Im $A$ is the imaginary part of the amplitude.

In particular, as the first ingredient of our tamed collinear factorisation
approach, we consider the resummation of a class of large double logarithms of
the form $\alpha_{s}\ln(1/\xi)\ln(\mu_{F}^{2}/m_{c}^{2})$ [4]. In the high-
energy limit, the NLO amplitude becomes directly proportional to a logarithm
of this kind and so the judicious factorisation scale choice $\mu_{F}=m_{c}$,
motivated by the structure of the amplitude in the high-energy asymptotics,
eliminates this large NLO contribution, as well as simultaneously resuming
this important class of double logarithmically enhanced contributions.111Note
that the alternate $k_{t}$-factorisation scheme intrinsically resums the same
class of large logarithms, but in this framework only a subset of NLO
corrections belonging to the equivalence class of gluon-ladder diagrams are
considered. The setting of the factorisation scale translates into a shifting
of terms between the higher-order contributions. For example, at some
factorisation scale $\mu$, the quark initiated contribution is part of the NLO
hard matrix element while at another factorisation scale $\mu\textquoteright$,
it is absorbed into the LO result. Our ideology is to find a so-called
‘optimal’ scale that removes the largest contribution from the NLO correction
and that, as explained above, is given by the choice $\mu_{F}=m_{c}$. The
effect of the scale change from a given $\mu_{F}$ to $\mu_{f}$ is driven by
DGLAP evolution, in which generalised (skewed) splitting functions $V$ are
present, and the change in the LO amplitude is given by

$A^{(0)}(\mu_{f})=\left(C^{(0)}+\frac{\alpha_{s}}{2\pi}\ln\left(\frac{\mu_{f}^{2}}{\mu_{F}^{2}}\right)C^{(0)}\otimes
V\right)\otimes F(\mu_{F}),$ (2)

while to NLO it is

$A(\mu_{f})=C^{(0)}\otimes F(\mu_{F})+C^{(1)}(\mu_{F})\otimes F(\mu_{f}).$ (3)

Figure 4 (left) shows the effect of this double logarithmic resummation. Here,
the factorisation scale $\mu_{F}=m_{c}$ is fixed numerically and variations
are made with respect to the residual factorisation scale $\mu_{f}$,
introduced above. As expected, the scale dependence improves noticeably but
still not satisfactorily as, for example, the NLO result is still of the
opposite sign to the LO one. The variations can be eased further by
considering an additional class of contributions which take the shape of power
corrections $O(Q_{0}^{2}/\mu_{F}^{2})$ to the standard coefficient functions
[5]. These constitute the other ingredient in the tamed approach and arise
through the implementation of a crucial $`Q_{0}$’ subtraction which cures an
inherent double counting in the standard collinear factorisation prescription
within the $\overline{\text{MS}}$ scheme. Explicitly, in this procedure, the
contribution from scales $l^{2}<Q_{0}^{2}$ circulating in the Feynman diagram
virtual loops are subtracted from the full NLO coefficient function. This
avoids the double counting in the convolution of the coefficient functions
with the input GPDs at some $Q_{0}$. Figure 4 (right) shows the effect of the
double logarithmic resummation together with the $Q_{0}$ subtraction. The
scale variation has dramatically improved and the results are indicative of a
perturbative convergence, suitable for a reliable comparison with experimental
data. We now provide some further discussion of the $Q_{0}$ subtraction.

Figure 4: As in Fig. 3 but now with $\mu_{F}=m_{c}$ fixed before (left) and
after (right) the double counting correction has been implemented, see text.
Scale variations are now made with respect to $\mu_{f}.$

##  4 Sensitivity to the $\overline{\text{MS}}$ scheme gluon PDF

It is important to emphasise that we remain in the $\overline{\text{MS}}$
scheme while working with $Q_{0}$ subtracted coefficient functions to NLO
accuracy. Specifically, the subtraction does not affect the infrared collinear
or ultraviolet divergence renormalisation procedures and the soft singularity
at $l=0$ is removed by subtracting off the LO part of the NLO coefficient
function before the integral over the loop momentum from 0 to $Q_{0}$ is
performed. It is precisely this finite contribution from $l<Q_{0}$ that is
subtracted to avoid the double counting inherent within the
$\overline{\text{MS}}$ scheme. To NLO this subtraction procedure does not give
rise to a new set of splitting functions and PDF evolution that otherwise a
new scheme would entail. We remark that this subtraction is fundamentally
ubiquitous and should be employed to correct for the double counting in all
$\overline{\text{MS}}$ coefficient functions, see [10] for the procedure
applied to inclusive DIS and Drell-Yan production. However, as these
additional corrections are of the form $O(Q_{0}^{2}/\mu_{F}^{2})$, when the
typical factorisation scale of the process $\mu_{F}\gg Q_{0}$, the effects of
this subtraction are not visible and only become relevant for lower scale
processes, like exclusive $J/\psi$ photoproduction discussed here.222For
exclusive $J/\psi$ photoproduction, $\mu_{F}=O(m_{c})$. In our approach, as
discussed in Section 2, $\mu_{F}=m_{c}$ and so the ratio $Q_{0}/m_{c}$ = O(1)
and the subtraction is therefore crucial.

##  5 Conclusion and final remarks

In these proceedings, we have presented a comparison of approaches of
exclusive $J/\psi$ photoproduction to NLO within collinear factorisation and
shown that a programme of resummation together with consideration of
additional crucial power corrections solves the long known problem of the
large scale sensitivity exhibited by the exclusive $J/\psi$ photoproduction
process. In an earlier analysis [11], we showed that analogous predictions to
those presented in Section 3 gave cross-section predictions that agreed
favourably with the HERA data for exclusive $J/\psi$ photoproduction, and with
good scale stability. The large difference between the predictions obtained
using a variety of different global PDFs in the larger $W$ LHC regime,
together with their reconciliation in the lower $W$ HERA regime, motivated the
extraction of a low $x\sim 3\times 10^{-6}$ and low $\mu\sim m_{c}$ gluon PDF
via two approaches which were shown to be consistent [12]. We emphasised that,
as the global fit PDF errors are so much greater than the exclusive $J/\psi$
experimental data uncertainties and the underlying scale uncertainty
associated with the theory prediction in our tamed collinear factorisation
framework, the exclusive $J/\psi$ data are in a position to readily improve
the global PDF analyses. In this way, the exclusive $J/\psi$ data at large $W$
from LHC will constrain the low $x$ and low $\mu$ behaviour of the gluon PDF,
while the data at lower $W$ from HERA can provide constraints in the regime
where there is the aforementioned multitude of data constraints from
inclusive, as well as diffractive, DIS experiments spanning a wide range of
momentum fractions $x$ and scales $\mu$.

## Acknowledgements

CAF is supported by the Helsinki Institute of Physics core funding project
QCD-THEORY. The work of TT is supported by the STFC Consolidated Grant
ST/T000988/1.

## References

  * [1] S. Bailey, T. Cridge, L. A. Harland-Lang, A. D. Martin and R. S. Thorne,  _Parton distributions from LHC, HERA, Tevatron and fixed target data: MSHT20 PDFs_ ,  Eur. Phys. J. C 81(4), 341 (2021),  10.1140/epjc/s10052-021-09057-0,  2012.04684.
  * [2] R. D. Ball _et al._ ,  _Parton distributions from high-precision collider data_ ,  Eur. Phys. J. C 77(10), 663 (2017),  10.1140/epjc/s10052-017-5199-5,  1706.00428.
  * [3] T.-J. Hou _et al._ ,  _New CTEQ global analysis of quantum chromodynamics with high-precision data from the LHC_ ,  Phys. Rev. D 103(1), 014013 (2021),  10.1103/PhysRevD.103.014013,  1912.10053.
  * [4] S. P. Jones, A. D. Martin, M. G. Ryskin and T. Teubner,  _Exclusive $J/\psi$ and $\Upsilon$ photoproduction and the low $x$ gluon_,  J. Phys. G 43(3), 035002 (2016),  10.1088/0954-3899/43/3/035002,  1507.06942.
  * [5] S. P. Jones, A. D. Martin, M. G. Ryskin and T. Teubner,  _The exclusive $J/\psi$ process at the LHC tamed to probe the low $x$ gluon_,  Eur. Phys. J. C 76(11), 633 (2016),  10.1140/epjc/s10052-016-4493-y,  1610.02272.
  * [6] A. D. Martin, C. Nockles, M. G. Ryskin, A. G. Shuvaev and T. Teubner,  _Generalised parton distributions at small $x$_,  Eur. Phys. J. C 63, 57 (2009),  10.1140/epjc/s10052-009-1087-y,  0812.3558.
  * [7] P. Hoodbhoy,  _Wave function corrections and off forward gluon distributions in diffractive $J/\psi$ electroproduction_,  Phys. Rev. D 56, 388 (1997),  10.1103/PhysRevD.56.388,  hep-ph/9611207.
  * [8] D. Y. Ivanov, A. Schafer, L. Szymanowski and G. Krasnikov,  _Exclusive photoproduction of a heavy vector meson in QCD_ ,  Eur. Phys. J. C 34(3), 297 (2004),  10.1140/epjc/s2004-01712-x,  [Erratum: Eur.Phys.J.C 75, 75 (2015)],  hep-ph/0401131.
  * [9] A. Buckley, J. Ferrando, S. Lloyd, K. Nordström, B. Page, M. Rüfenacht, M. Schönherr and G. Watt,  _LHAPDF6: parton density access in the LHC precision era_ ,  Eur. Phys. J. C 75, 132 (2015),  10.1140/epjc/s10052-015-3318-8,  1412.7420.
  * [10] A. D. Martin and M. G. Ryskin,  _Comments on Global Parton Analyses_ ,  Acta Phys. Polon. B 51, 1079 (2020),  10.5506/APhysPolB.51.1079,  1912.09304.
  * [11] C. A. Flett, S. P. Jones, A. D. Martin, M. G. Ryskin and T. Teubner,  _How to include exclusive $J/\psi$ production data in global PDF analyses_,  Phys. Rev. D 101(9), 094011 (2020),  10.1103/PhysRevD.101.094011,  1908.08398.
  * [12] C. A. Flett, A. D. Martin, M. G. Ryskin and T. Teubner,  _Very low $x$ gluon density determined by LHCb exclusive $J/\psi$ data_,  Phys. Rev. D 102, 114021 (2020),  10.1103/PhysRevD.102.114021,  2006.13857.