Microarcsecond VLBI Pulsar Astrometry with PSRπ II. Parallax Distances for 57 Pulsars

An initial list of target pulsars was produced consisting of sources located north of −20° decl. and with a "gated equivalent" flux density¹² at 1400 MHz sufficient to obtain a detection exceeding 35σ within a single PSRπ astrometric observation. This sample consisted of 225 pulsars with a gated equivalent flux density >3.2 mJy (bright enough for observations at the then-available data rate of 512 Mbps) and a further 55 sources with a gated equivalent flux density between 1.6 and 3.2 mJy (bright enough for future observations at a data rate of 2 Gbps). The first phase of PSRπ observations entailed the identification of suitable compact background sources close to the potential target pulsars on the sky that could be used as secondary phase calibrators ("in-beam" calibrators). The astrometric positions of the pulsars are ultimately measured relative to these sources. A pilot program testing the observing strategy was undertaken between 2010 February and May (40 pulsars, 12 hr, Very Long Baseline Array (VLBA) project code BD148), which included some eventual PSRπ targets. In-beam calibrator identification observations for the remaining potential PSRπ targets were undertaken in the main PSRπ observing program (240 pulsars, 85 hr, VLBA project code BD152) between 2010 November and 2011 December. In all cases, all potentially useful candidates within ∼25' of the target pulsar were investigated using the multifield capability of the DiFX software correlator (Deller et al. 2011). The central observing frequency was 1660 MHz, and phase referencing was performed using a nearby calibrator to a grid of four pointing centers arrayed around the target pulsar, with R.A. and decl. offsets of ±10' in each direction. Figure 1 illustrates an example pointing layout, for the target J1136+1551.

Zoom In Zoom Out Reset image size

The candidate sources were taken from the Faint Images of the Radio Sky (FIRST; Becker et al. 1995) catalog where available, and the NRAO VLA Sky Survey (NVSS; Condon et al. 1998) in areas not covered by FIRST. For future campaigns, the ongoing VLA Sky Survey (VLASS; M. Lacy et al. 2019, in preparation) will provide a deeper and higher resolution catalog covering the full NVSS footprint. Each pointing center was visited for approximately 3.5 minutes, which with an observing bandwidth of 64 MHz (dual polarization) gave a typical on-source root-mean-square (rms) noise of 0.3–0.8 mJy beam⁻¹, depending on the location of the candidate in the pointing pattern. In general, this was sufficient to identify (at >6.5σ) any candidates within ∼25' of the target pulsar brighter than 3 mJy (the completeness limit of NVSS) that could potentially serve as useful calibrators or astrometric check sources. The calibration and source detection approach was essentially the same as that used by the mJIVE-20 project (Deller & Middelberg 2014), which was inspired by the procedure undertaken here.

These calibrator search observations for PSRπ served as a survey of over 200 square degrees at milliarcsecond resolution complete to ∼3 mJy, reliably detecting over 1500 sources. Over 90% of the 280 targeted pulsars were found to be located near at least one source suitable as a secondary phase calibrator for high-sensitivity observations (flux density >3 mJy contained within a component of maximum size several milliarcseconds, angular separation <25').

From our initial sample of 280 pulsars, ∼110 met our requirements for astrometric observations with the then-available 512 Mbps recording system on the VLBA capable of recording dual polarization 64 MHz bandwidth. These requirements were: pulsar gated equivalent flux density >3.2 mJy, and at least one compact secondary calibrator within 25' with flux density >6 mJy. We observed each of these 110 sources once using the astrometric scheme described in Section 2.2, before down-selecting to the final sample of 60 pulsars. The initial observations were used to reduce the risk of selecting a target where the final astrometric precision would be insufficient to provide a useful distance constraint. Some targets were rejected due to an unsuitable secondary calibrator, generally due to complicated source structure that was not identified in the initial snapshot search observations due to limited uv coverage. Other targets were rejected because that pulsar's observed flux density was much fainter than the catalog value. Finally, from the remaining viable targets, preference was given to sources that sampled a range of Galactic longitudes and latitudes (which meant primarily discarding sources in the Galactic plane and located toward to inner Galaxy), and sources that did not already have a high-precision parallax distance. The final selection of 60 targets was therefore based on both logistical (range of right ascensions, strong calibrators close on the sky to the target, target flux density) and scientific (range of Galactic heights, predicted distances, individual objects where the distance was a high priority) considerations. The 60 selected sources are summarized in Table 1.

Table 1.  PSRπ Targets

Pulsar	Pulsar	DM	S1.4g	Obs. Freq.	DNE2001	DYMW16
Jname	Bname	(pc cm−3)	(mJy)a	(MHz)b	(kpc)c	(kpc)d
J0040+5716	B0037+56	92.6	4.7	1660	3.05	2.42
J0055+5117	B0052+51	44.1	7.6	1660	1.90	1.94
J0102+6537	B0059+65	65.9	4.9	1660	2.29	1.98
J0108+6608	B0105+65	30.5	5.2	1660	1.42	1.46
J0147+5922	B0144+59	40.1	10.1	1660	2.22	1.58
J0151–0635	B0148–06	25.7	5.4	1660	1.22	25.00
J0152–1637	B0149–16	11.9	9.2	1660	0.51	0.92
J0157+6212	B0154+61	30.2	10.6	1660	1.71	1.39
J0323+3944	B0320+39	26.0	6.5	1660	1.01	1.20
J0332+5434	B0329+54	26.8	1244.9	1660	0.98	1.18
J0335+4555	B0331+45	47.2	4.6	1660	1.64	1.53
J0357+5236	B0353+52	103.7	6.1	1660	2.78	2.02
J0406+6138	B0402+61	65.3	13.5	1660	2.12	1.78
J0601–0527	B0559–05	80.5	11.8	1660	3.93	2.33
J0614+2229	B0611+22	96.9	12.5	1660	2.08	1.74
J0629+2415	B0626+24	84.2	17.9	1660	2.24	1.67
J0729–1836	B0727–18	61.3	8.0	1660	2.90	2.40
J0823+0159	B0820+02	23.7	8.3	1660	1.01	0.81
J0826+2637	B0823+26	19.5	76.4	1660	0.34	0.31
J1022+1001	⋯	10.3	9.5	1660	0.45	0.83
J1136+1551	B1133+16	4.9	181.9	1660	0.34	0.41
J1257–1027	B1254–10	29.6	7.3	1660	1.55	25.00
J1321+8323	B1322+83	13.3	4.3	1660	0.76	0.98
J1532+2745	B1530+27	14.7	4.8	1660	0.83	1.32
J1543–0620	B1540–06	18.4	15.2	1660	0.72	1.12
J1607–0032	B1604–00	10.7	26.2	1660	0.67	0.68
J1623–0908	B1620–09	68.2	4.7	1660	50.00	25.00
J1645–0317	B1642–03	35.7	167.4	1660	1.12	1.32
J1650–1654	⋯	43.2	8.7	1660	1.47	1.05
J1703–1846	B1700–18	49.6	4.8	1660	1.48	1.69
J1735–0724	B1732–07	73.5	9.5	1660	2.26	0.21
J1741–0840	B1738–08	74.9	6.5	1660	2.17	0.22
J1754+5201	B1753+52	35.4	9.3	1660	2.18	4.17
J1820–0427	B1818–04	84.4	37.9	2267	1.94	2.92
J1833–0338	B1831–03	234.5	18.8	2267	5.14	5.17
J1840+5640	B1839+56	26.7	25.0	1660	1.68	2.19
J1901–0906	⋯	72.7	20.4	1660	2.13	2.89
J1912+2104	B1910+20	88.3	7.4	1660	3.96	3.37
J1917+1353	B1915+13	94.5	10.4	2267	3.99	2.94
J1913+1400	B1911+13	145.1	7.3	2267	5.12	5.25
J1919+0021	B1917+00	90.3	5.7	1660	3.06	4.10
J1937+2544	B1935+25	53.2	7.7	1660	3.25	2.87
J2006–0807	B2003–08	32.4	8.8	1660	1.23	1.71
J2010–1323	⋯	22.2	5.8	1660	1.02	1.16
J2046–0421	B2043–04	35.8	12.7	1660	1.75	3.27
J2046+1540	B2044+15	39.8	10.5	1660	2.42	3.34
J2113+2754	B2110+27	25.1	8.9	1660	2.03	1.87
J2113+4644	B2111+46	141.3	62.9	1660	4.53	4.12
J2145–0750	⋯	9.0	21.2	1660	0.57	0.69
J2149+6329	B2148+63	128.0	11.0	1660	5.51	3.88
J2150+5247	B2148+52	148.9	8.6	1660	4.62	3.61
J2212+2933	B2210+29	74.5	3.8	1660	4.20	25.00
J2225+6535	B2224+65	36.1	10.0	1660	1.86	1.88
J2248–0101	⋯	29.1	4.6	1660	1.65	25.00
J2305+3100	B2303+30	49.5	17.2	1660	3.66	25.00
J2317+1439	⋯	21.9	10.3	1660	0.83	2.16
J2317+2149	B2315+21	20.9	6.5	1660	0.95	1.80
J2325+6316	B2323+63	197.4	6.7	1660	8.26	4.86
J2346–0609	⋯	22.5	8.2	1660	0.94	25.00
J2354+6155	B2351+61	94.7	31.6	1660	3.43	2.40

Notes.

^aGated equivalent flux density; calculated from catalog 1.4 GHz flux density scaled by $\sqrt{\mathrm{duty}\,\mathrm{cycle}}$. ^bObserving frequency used for the majority of astrometric observations (see Section 2). ^cDistance estimated from the DM and the NE2001 Galactic electron density distribution (Cordes & Lazio 2002). ^dDistance estimated from the DM and the YMW16 Galactic electron density distribution (Yao et al. 2017).

Download table as:  ASCIITypeset image

Each of the 232 PSRπ astrometric epochs lasted ∼2.4 hr and targeted two pulsars located relatively close to each other on the sky, with typical angular separations of 10°–20°. In each observation, five fields were observed: the two target fields (each of which encompassed both a pulsar and one or more in-beam calibrators), two out-of-beam phase reference sources (one near each pulsar), and one strong "fringe finder" source used to calibrate the instrumental bandpass. The observing sequence was as follows:

• 1.  
  Five 5.5 minutes scans on the first target field interleaved and bracketed by 1.25 minutes scans on the associated primary calibrator
• 2.  
  Five 5.5 minutes scans on the second target field interleaved and bracketed by 1.25 minutes scans on the associated primary calibrator
• 3.  
  A 2-minute scan on the fringe finder
• 4.  
  Five 5.5 minutes scans on the first target field interleaved and bracketed by 1.25 minutes scans on the associated primary calibrator
• 5.  
  Five 5.5 minutes scans on the second target field interleaved and bracketed by 1.25 minutes scans on the associated primary calibrator

This observing sequence was used to ensure that the uv coverage was maximized while keeping the slewing overheads relatively low.

In some cases, a known and suitable VLBI calibrator was separated by less than 25' from the pulsar. In these cases, no nodding calibration was performed, reducing the observing block for that pulsar to a single 27.5 minutes long scan (repeated twice during the observation).

The pointing center for the target fields was typically chosen to be close to the midpoint between the pulsar and the primary in-beam calibrator, although adjustments were made based on the location of additional in-beam calibrator sources in some cases. Figure 2 shows an example astrometric pointing layout for the target PSR J1136+1551, while https://safe.nrao.edu/vlba/psrpi/astrometric_pointings.html shows the pointing layout for all the target pulsars.

Zoom In Zoom Out Reset image size

The astrometric observations were scheduled during the period 2011 January to 2013 December, and were optimized in time to provide maximum sensitivity to annual geometric parallax. Since the VLBA is more extended in the east–west direction than north–south, the synthesized beam is narrower in R.A. than in decl. Accordingly, we scheduled our astrometric observations around the time of the peak parallax signature in R.A.. After the first "check" observation, every successive parallax extremum was sampled with two observations within a ∼20 day period. In total, each pulsar was observed eight or nine times spread over a ∼2 yr period.

The observational setup consisted of four 16 MHz subbands, covering both circular polarizations and sampled at 2 bit precision for a total recording rate of 512 Mbps. In the first astrometric observation of each pulsar, the frequency range chosen was 1624.49–1688.49 MHz, while for the remaining observations the frequency range was shifted slightly to 1627.49–1691.49 MHz, due to strong interference from the Iridium satellite constellation at ∼1625 MHz. For four pulsars (PSRs J1820−0427, J1833−0338, J1913+1400, and J1917+1353), strong scatter-broadening led us to choose observations at higher frequency, 2234.49–2298.49 MHz. In the first observation of each pulsar, a significant portion of the data from the first subband was flagged due to this interference, while in the later observations, a lesser amount of data was flagged in the fourth subband due to interference at ∼1690 MHz.

Correlation was performed in Socorro using the DiFX software correlator (Deller et al. 2011). For each observation, a minimum of three and a maximum of six correlation passes were made, forming 3–6 separate visibility data sets. The first pass correlated all sources, using the position of the primary in-beam calibrator for scans on the target field. We refer to the resultant data set henceforth as the "calibrator" data set. All other passes correlated only the scans on the target fields, and differed in the location of the phase center and the presence or absence of special pulsar processing. The second pass used the position of the target pulsars and employed pulsar gating to boost the signal-to-noise ratio (S/N) by down-weighting time ranges when the pulsar signal is weak or absent (Deller et al. 2007). The pulsar ephemeris and the pulsar gate parameters were obtained from timing observations with the Lovell telescope at Jodrell Bank Observatory. We refer to the output as the "gated" data set. The third pass also used the position of the target pulsars, but used no pulsar gating, to generate the "ungated" data set. Comparison of the S/N of the pulsar detection in the gated and ungated data sets allowed us to check that the correct pulsar ephemeris and gate had been applied. If present, the fourth, fifth, and sixth correlator passes used the position of the additional in-beam calibrator sources; we refer subsequently to these as the "additional" data sets.

In order to generate artifact-free images of the target pulsars that are located in a stable reference frame, careful calibration is needed to remove time- and frequency-dependent corruption of the measured visibilities by instrumental and propagation effects. All calibration was performed using AIPS (Greisen 2003), facilitated by the ParselTongue python interface (Kettenis et al. 2006). AIPS version 31DEC15 was used for the final data processing, and each pulsar was processed independently. Calibration was script-based, with configurable options set using a markup-language control file for traceability. The cumulative calibration derived in the preceding stages is always applied before solving for the next stage of the incremental calibration. We now describe the calibration script stages in detail.

• 1.  
  Load. The visibility data sets were loaded into AIPS using the task FITLD.
• 2.  
  A priori flagging. Logged time ranges when the antennas are slewing, settling, or have otherwise known pointing or recording problems are already recorded in the flag table accompanying the data set. In addition to these existing flags, we also flagged baselines during times when the natural fringe rate was low using the AIPS task UVFLG, as these are susceptible to corruption by radio frequency interference (RFI) and instrumental effects. We also flagged baselines when one or both antennas was pointing below 20° elevation (AIPS task UVFLG). Finally, we applied any user-defined flags that were generated after inspection of final data products of an earlier pipeline run (AIPS task UVFLG).
• 3.  
  Source shifting. If required, the phase center for one or more sources was shifted using the AIPS task CLCOR. This was typically only needed for the first observation, where the pulsar position was sometimes poorly known, and in some cases the in-beam calibrator position was also only poorly constrained after the initial snap-shot observation.
• 4.  
  A priori ionosphere correction. The delay model applied at the correlator does not include any ionospheric contribution. To correct for ionospheric propagation delays, we used the AIPS task TECOR, which makes use of a low-resolution global ionosphere model. While these global models are unable to remove rapid and/or small-scale variations, they do account for bulk ionospheric effects. We used the the final combined analysis models of the International GNSS Service (analysis center code igsg) available from ftp://cddis.gsfc.nasa.gov/gps/products/ionex/.
• 5.  
  EOP corrections. The Earth Orientation Parameters (EOPs) used in the correlator model are often refined after the time of correlation. To update the visibilities and make use of the most accurate available EOPs in order to minimize residual position offsets after phase referencing, we used the AIPS task CLCOR.
• 6.  
  A priori amplitude calibration. The visibilities in the correlator data set are scaled to take the form of pseudo-correlation coefficients. To convert these to an approximate flux density scale in Janskys, the following steps are taken:

  • (a)  
    Quantization correction. Imperfect level-setting in the quantizers biases the amplitude scale; this effect can be detected and corrected by analysis of the station autocorrelations. We used the AIPS task ACCOR to make these corrections.
  • (b)  
    System temperature correction. A continuously operating switched noise diode at each VLBA antenna records the system temperature at that antenna. In combination with an a priori gain curve (which, at low frequencies such as used here, is quite accurate), this can be used to convert the pseudo-correlation coefficients produced by the correlator into Janskys. We use the AIPS task APCAL to generate these corrections.
  • (c)  
    Primary beam correction. In the target pointing, the pulsar and in-beam calibrator(s) are not centered in the primary beam. Accordingly, the amplitude response for each of these sources is attenuated by the primary beam fall-off. We apply a correction based on a simplified model of a uniformly illuminated antenna scaled by the measured parameters for beamwidth and beam squint of VLBA antennas at our observing band, using a custom ParselTongue script described in Deller & Middelberg (2014).
• 7.  
  Instrumental phase calibration. Any instrumental phase variations due to changing propagation through the signal chain are tracked by an injected pulse train, and the measured phases are stored in a table that accompanies the visibilities. We applied these corrections using the AIPS task PCCOR.
• 8.  
  Time-independent delay calibration. Using the AIPS task FRING, we measured the single-band delays for each subband independently on the fringe finder source. A model of the fringe finder source derived from imaging the concatenated data for the source from all 8–9 PSRπ epochs was supplied to FRING. The AIPS task SNSMO was used to apply a median window filter to the resultant delays and automatically exclude any solutions that differed by more than 10 nanoseconds from the median for that subband of that antenna. Rates were zeroed before the delays were applied using CLCAL.
• 9.  
  Time-independent instrumental bandpass calibration. The AIPS task BPASS was used to derive the instrumental bandpass, using the fringe finder scan. As with the preceding step, the model of the fringe finder source was supplied to BPASS. The resultant amplitude corrections were normalized to leave the flux density scale unaffected.
• 10.  
  Time-dependent delay calibration. We used the AIPS task FRING to now derive single-band delays using the phase reference calibrator source. Again, a source model (based on imaging of concatenated PSRπ data sets) was supplied in all cases. In almost all cases, each subband was solved separately, but for several weak phase reference sources, we combined all subbands together to improve the S/N. Solving for all subbands separately is preferred, because the residual ionospheric delays can lead to a (slightly) varying delay between subbands. When this is the case, approximating the delay as constant over all subbands reduces the S/N improvement resulting from the larger bandwidth, and also leaves per-subband phase residuals that must be corrected later. The solutions were median-window filtered using SNSMO to excise values where the delay or phase offsets exceeded 10 ns or 10 mHz, respectively.
• 11.  
  Refinement of time-dependent amplitude calibration. The AIPS task CALIB was used to compute amplitude self-calibration corrections on a per-subband basis for the phase reference calibrator source, usually with a timescale of 20 minutes (shorter in some cases where the calibrator source was particularly strong). As in previous steps, the epoch-averaged calibrator model was employed, which in effect forces the absolute flux density scale to match this model. If the calibrator exhibited substantial flux density variations over the observing period, this could affect the flux density scale of individual observations; however, the absolute flux density scale is not important for the astrometric observables, which depend only on position. SNSMO was used to median-window filter and remove amplitude corrections differing by more than 20% from the median over the whole observation.
• 12.  
  Refinement of time-dependent phase calibration in target direction. The following steps now made use of the strongest source(s) in the target pointing, with the results then being applied to all sources in the target pointing. In most cases, a single in-beam calibrator source was used to derive solutions. For some targets where no strong in-beam calibrator source was available, multiple in-beam sources were used together in a "multi-source selfcal" (Middelberg et al. 2013; Radcliffe et al. 2016). Finally, in several cases (PSR J0332+5434, PSR J1136+1551, and PSR J2113+4644), the pulsar itself was by far the strongest source in the field, and the gated pulsar data set was used to derive these calibrations rather than one of the other in-beam sources.

  • (a)  
    Frequency independent. The data set(s) corresponding to the designated source(s) were split using the AIPS task SPLIT, to apply calibration and flagging and average the data in frequency for speed. For non-pulsar sources, the split data set was then divided by the corresponding source model using the AIPS task UVSUB. If multiple in-beam sources were used in a "multi-source selfcal," these normalized data sets were then combined with the AIPS task DBCON. CALIB was then used to solve for phase corrections on this normalized (and possibly concatenated) data set, summing all subbands and polarizations. The solution interval ranged from 10 s for the brightest sources to 5 minutes for the weakest source, with a median value of 1.25 minutes.
  • (b)  
    Frequency dependent. The preceding step was repeated, but this time treating frequency subbands separately while still summing polarizations. The reduced bandwidth per solution was compensated with a longer solution interval, typically 6 minutes to 30 minutes. This step compensates for the small residual dispersive delays due to the ionosphere between the phase reference calibrator direction and the target direction.
• 13.  
  Refinement of time-dependent amplitude calibration in target direction. For a handful of targets with sufficiently bright in-beam calibrators, we used CALIB to derive further self-calibration amplitude corrections using the in-beam calibrator data. SNSMO was applied to filter out solutions differing by more than 20% from the median.
• 14.  
  Correction of pulsar scintillation. For some nearby pulsars, the scintles produced by diffractive scintillation in the interstellar medium are sufficiently large compared to our visibility resolution (diffractive timescale ≫1 s, diffractive bandwidth ≫1 MHz) that significant amplitude variations are apparent in the pulsar data. When significant diffractive scintillation was present, we derived time-dependent amplitude correction factors using a custom ParselTongue routine described in Deller et al. (2009). The solution interval was empirically determined by inspection of the uncorrected pulsar data.
• 15.  
  Writing calibrated data. The fully calibrated data sets for the pulsar (both gated and ungated), the in-beam calibrator source(s), and the phase reference calibrator source were split and averaged to a single visibility point per frequency subband (excising the two edge channels). The non-pulsar sources were divided by the epoch-averaged source structure model, and all data sets (pulsar gated, pulsar ungated, calibrator sources, calibrator sources divided by model) were written to disk using the AIPS task FITTP.
• 16.  
  Producing log files. While the script was running, statistics on the failure rates of each calibration step were retained. At the conclusion of the script, these are written into a summary webpage, along with any information from the VLBA operator's log and plots of the delay, amplitude, and phase calibration tables generated by the script. Plots showing visibility amplitude as a function of baseline length are also generated for each source and included, to aid in the identification of unflagged RFI. If any evidence of RFI or unsatisfactory calibration was evident, then additional flagging was undertaken and the calibration script was rerun.

After flagging, calibration, and averaging, the visibility data are now in a suitable form for the extraction of our astrometric observables. Given the calibration that has been applied, the pulsar is effectively being determined using the method of differential astrometry with respect to the position of the in-beam calibrators. Among the 73 in-beam calibrators, the positions of 14 were found in the Radio Fundamental Catalog (RFC¹³ ; L. Petrov & Y. Kovalev 2019, in preparation) derived using observations designed to improve the absolute positions of VLBI calibrator sources. The absolute positions of the other in-beam calibrators were determined with respect to the primary out-of-beam phase reference calibrators. For the out-of-beam calibrators, the RFC2019a solution was used, and the absolute positions listed in the RFC have coordinate uncertainties ranging from 0.1 to 0.9 mas (median: 0.18 mas). J2000 coordinates are used for all sources.

Rather than these uncertainties transferred from the out-of-beam calibrators, the absolute position uncertainties of the in-beam calibrators, and hence of the pulsars, is generally dominated by one of the following two factors:

• 1.  
  A substantial angular separation from the out-of-beam calibrators to the target field, (up to 35), leading to biases in position estimates when the calibration is extrapolated. At an observing frequency of 1600 MHz, the uncertainty depends strongly on the ionospheric conditions during the observation(s) as well as the size of the angular separation and the median observing elevation; a few mas is typical (Deller et al. 2016). We estimate this value for each pulsar field by examining the scatter in the in-beam positions obtained when imaging using only phase referencing from the out-of-beam calibrator (i.e., no self-calibration on the in-beam calibrator).
• 2.  
  The frequency dependence of the core position of calibrator sources. The apparent core of calibrator sources, which defines their reference position, is determined from the region of peak brightness in an image. The true jet origin, the region at the jet apex, is invisible to an observer, as it is opaque (optical depth τ ≫ 1) due to synchrotron self-absorption. The apparent core is located where the jet becomes visible further away from the origin, when optical depth reaches τ ≈ 1 at the apparent jet base. The higher the frequency, the closer the observed core is to the jet apex. This effect is called the core-shift; more detailed descriptions can be found in, for example, Marcaide & Shapiro (1984), Lobanov (1998), Kovalev et al. (2008), O'Sullivan & Gabuzda (2009), Voitsik et al. (2018), Pushkarev et al. (2019). According to Sokolovsky et al. (2011), the core-shift at 1600 MHz ranges from 0.6 to 2.4 mas with median 1.1 mas for a specially pre-selected sample of active galactic nuclei (AGNs) jets with significant shifts. Pushkarev et al. (2012) has measured core-shifts for a large complete flux density limited sample between 8 and 15 GHz and has found median values similar to Sokolovsky et al. (2011). Another complication arises with long-term core-shift variability, which can reach significant values of up to 1 mas between 2 and 8 GHz (Plavin et al. 2019). In this work we did not determine the core-shift, and thus the neglected core-shift introduces a ∼mas level bias in the chain of tying the pulsar position to the out-of-beam calibrator via an intermediate in-beam calibrator.

All three sources of uncertainty are added in quadrature when determining the absolute position uncertainty of the pulsar. However, none of these issues affect any shifts relative to the in-beam calibrator, which are relevant for determining proper motion and parallax. The contributions to the relative astrometric error budget are considered in more detail in Section 3.2.

For each visibility data set (gated pulsar, ungated pulsar, calibrator source[s], and calibrator source[s] divided by their average model), we used the following procedure to extract a position, again implemented as a ParselTongue script. We loaded the data set into the difmap package (Shepherd 1997) and inverted the Stokes I visibility data to form a dirty image, using natural weighting. We then shifted the data set to be approximately centered on the peak in the dirty image, and added a single point-source component to the model at the location of the peak. We ran a model fit for 20 iterations, and then wrote the resultant clean image (pixel size 0.75 mas) to disk in FITS format. This clean image was then loaded into AIPS, and the position and position uncertainties were extracted with the task JMFIT to fit an elliptical Gaussian, using a 20 × 20 pixel window centered on the peak. An example image, showing the fitted Gaussian, best-fit position, and uncertainty, is shown in Figure 3.

Zoom In Zoom Out Reset image size

In principle, position information can be extracted directly from the model fit, without the need to form an image and fit an elliptical Gaussian. However, while the best-fit position is easily accessible, extracting a position uncertainty from a model fit is highly dependent on the overall scaling of the visibility weights. In contrast, the position uncertainty resulting from an image plane fit, where the root mean square fluctuations of a residual image can easily be measured, is well defined and robust under conditions that are typically satisfied for radio interferometric images (Condon 1997). It is for this reason that we use JMFIT to extract positions and position uncertainties, although we did cross-check our astrometric results using the model-fit positions and estimated uncertainties, finding results that typically agreed to well within 1σ (where the uncertainty was taken from the image plane fit).

After the source modeling, calibration, and imaging described in the previous section was completed, we were left with a position time series for the target pulsar (taken from the higher S/N gated data sets) and one or more in-beam reference sources. The reference frame in which these positions were measured has been defined by the assumed positions for the source(s) used in step 12 of Section 2.3. This was not necessarily optimal for measuring the time-varying position of the target pulsar—the brightest in-beam calibrator might be located far from the pulsar, or might exhibit structure evolution that result in additional systematic errors. In these cases, we could provide a more stable reference frame by utilizing a weaker and/or more stable background source located at a smaller angular separation to the target pulsar.

By subtracting the position residuals from a given in-beam source (or a weighted average of several in-beam sources), we could transform the reference frame into one defined by an in-beam source (or sources) of our choosing. By selecting the source nearest to the target, it is possible to minimize the systematic position shifts at the target pulsar position introduced by the residual ionosphere. However, this also necessitates adding the formal fit errors for the chosen position reference source(s) in quadrature to the pulsar position fit error, which may not be an acceptable trade-off for weak background sources where these errors are large. Accordingly, for each pulsar we selected the position reference source considering both of these factors, and Table 2 shows the calibrator source(s) and the frame-defining source(s) for each target pulsar, along with the flux density as measured from the combined reference image for each source (at 1670 MHz for most sources, but at 2270 MHz for the calibrators of PSRs J1820−0427, J1833−0338, J1913+1400, and J1917+1353). Finally, any estimate of additional systematic sources of position uncertainty (discussed later) was added to the formal position errors; these were initially set to zero.

This time series of measured positions and estimated uncertainties could then be processed using the pmpar¹⁴ package to perform least-squares minimization and fit for reference position, proper motion, and parallax. Four of the PSRπ targets are a pulsar in a binary system, and where the orbital reflex motion is substantial, additional steps were required in the fitting process. PSR J1022+1001 and PSR J2145−0750 are MSPs that have already been described in Deller et al. (2016), while PSR J0823+0159 is a slow pulsar in a long-period binary. For these pulsars, the orbital period, longitude of periastron, eccentricity, and projected semimajor axis were all well-constrained by pulsar timing, and we therefore were only required to fit for inclination i and longitude of ascending node Ω. For PSR J2317+1439, the orbital reflex motion is negligible compared to our positional uncertainties. For the pulsars where fitting the orbital reflex motion was required, we included these two additional parameters in our least-squares minimization.

In almost all of the PSRπ pulsars, the reduced χ² of the initial least-squares fit exceeded unity, often considerably. This result is not surprising, given that the initial input position uncertainties are purely based on the S/N of the pulsar (and position reference calibrator) images, and does not account for potential systematic position shifts. In most cases, the dominant systematic contribution comes from the residual unmodeled ionosphere, but other possibilities such as source structure evolution in the source(s) defining the reference frame also exist. In general, the distribution (both form and variance) of these systematic errors is extremely difficult to predict a priori, as discussed in Section 3.2, which complicates efforts to accurately estimate the uncertainties on the fitted astrometric parameters. The problem is exacerbated for data sets where the formal position uncertainties vary widely between epochs, as can be the case for pulsars that exhibit significant amplitude variability due to diffractive and/or refractive scintillation. If no adjustment is made to the formal position errors, then the epochs with high-significance detections when the pulsar was "scintillated up" will exhibit a disproportionate impact on the astrometric fit.

In Section 3.2, we investigate different methods for estimating a systematic error that can be added in quadrature to the formal position fit errors in order to mitigate this issue. Unsurprisingly, we find that no method is perfect in all situations, but that the use of an estimator is better than neglecting systematic errors entirely. Our final astrometric solutions therefore make use of the empiral systematic error estimator discussed in Section 3.2.

Once a position time series with final estimated uncertainties is available, best-fit values and uncertainties for the astrometric parameters must be produced. Two options are available:

• 1.  
  A simple least squares fit
• 2.  
  A bootstrap fit

The least-squares fit has the advantage of simplicity, but is sensitively dependent on beginning with a good estimate of the input position errors. If these are underestimated (which will generally result in a reduced χ² that still significantly exceeds unity), then the errors on the astrometric observables will likewise be underestimated. Conversely (but more rarely), overestimating the systematic errors will lead to inflated uncertainties on the astrometric observables. A bootstrap fit (e.g., Efron & Tibshirani 1991) utilizes a large number of trials, where in each trial N position measurements for the input data set for each trial are selected randomly with replacement from the available N astrometric position measurements for that pulsar. In our case, N is usually 8 or 9. For each trial data set, a least-squares fit is made as usual, and the best-fit parameters are saved. After many trials, a cumulative probability distribution for each of the fitted parameters is built, from which the most probable value and a desired confidence interval can be extracted.

A bootstrap fit has the advantage that the uncertainty on the fitted parameters is not determined solely by the uncertainty in the input position measurements, which as we have seen is hard to estimate accurately. However, the bootstrap approach can exacerbate a problem already present for PSRπ and most VLBI astrometry programs: the small sample size. With just eight or nine position measurements, a significant fraction of trials can end up with poor time coverage of one of the desired astrometric quantities, sampling a shorter time range or predominantly one side of the parallax signature. This is especially problematic in cases where the pulsar scintillates and is detected only weakly (or not at all) in some epochs, further reducing the number of useful degrees of freedom. An example is PSR J2317+1439, where non-detections due to unfavorable scintillation were concentrated in the December/January observations and resulted in a poor sampling of the parallax ellipse.

We favor a bootstrap approach for determining the final astrometric uncertainties, as it generally produces the most conservative error estimates (as can be seen in Section 3.2). The results presented here are obtained from a bootstrap with 100,000 trials per pulsar. We highlight the circumstances under which the bootstrap uncertainties may potentially be too conservative in the discussion.

The astrometric results for our 60 target pulsars are shown in Table 3. Asymmetric error bars representing the 68% confidence interval are listed along with the best-fit parameter values. The median parallax uncertainty obtained was 46 μas, with 60% of our targets meeting or exceeded the design goal of 50 μas parallax accuracy. Almost all (53 of the 60) target pulsars have a significant (>95% confidence) parallax measurement, while two-thirds of the sample provide a distance error of 20% or less.

Table 3.  Fitted Astrometric Parameters for All PSRπ Targets

Pulsar	Offset from Reference		Proper Motion		Parallax
	R.A.	Decl.	R.A.	Decl.
	(mas)	(mas)	(mas yr−1)	(mas yr−1)	(mas)
J0040+5716	$-{119735.120}_{-0.045}^{+0.037}$	${783183.975}_{-0.109}^{+0.064}$	${12.399}_{-0.057}^{+0.033}$	$-{5.450}_{-0.082}^{+0.140}$	${0.102}_{-0.025}^{+0.051}$
J0055+5117	$-{332840.251}_{-0.023}^{+0.071}$	$-{302297.263}_{-0.074}^{+0.085}$	${10.490}_{-0.085}^{+0.049}$	$-{17.352}_{-0.204}^{+0.074}$	${0.349}_{-0.055}^{+0.055}$
J0102+6537	${47321.472}_{-0.006}^{+0.074}$	${79472.077}_{-0.045}^{+0.041}$	${9.252}_{-0.081}^{+0.049}$	${1.828}_{-0.206}^{+0.093}$	${0.399}_{-0.045}^{+0.044}$
J0108+6608	$-{140564.314}_{-0.033}^{+0.039}$	${27077.203}_{-0.025}^{+0.050}$	$-{32.754}_{-0.025}^{+0.036}$	${35.162}_{-0.051}^{+0.024}$	${0.468}_{-0.031}^{+0.035}$
J0147+5922	$-{744046.103}_{-0.038}^{+0.100}$	$-{189510.613}_{-0.054}^{+0.093}$	$-{6.380}_{-0.101}^{+0.083}$	${3.826}_{-0.097}^{+0.054}$	${0.495}_{-0.093}^{+0.042}$
J0151–0635	$-{582099.427}_{-0.014}^{+0.164}$	$-{358146.639}_{-0.052}^{+0.066}$	${10.697}_{-0.145}^{+0.094}$	$-{5.373}_{-0.078}^{+0.061}$	${0.217}_{-0.076}^{+0.098}$
J0152–1637	$-{1067819.606}_{-0.014}^{+0.158}$	$-{399877.713}_{-0.330}^{+0.363}$	${0.804}_{-0.201}^{+0.234}$	$-{31.372}_{-0.313}^{+0.424}$	${0.443}_{-0.181}^{+0.214}$
J0157+6212	${814250.347}_{-0.007}^{+0.035}$	${325456.249}_{-0.010}^{+0.037}$	${1.521}_{-0.017}^{+0.105}$	${44.811}_{-0.048}^{+0.034}$	${0.554}_{-0.024}^{+0.039}$
J0323+3944	${401715.413}_{-0.011}^{+0.046}$	$-{189856.218}_{-0.017}^{+0.029}$	${26.484}_{-0.034}^{+0.059}$	$-{30.780}_{-0.015}^{+0.029}$	${1.051}_{-0.040}^{+0.039}$
J0332+5434	$-{158517.350}_{-0.008}^{+0.035}$	$-{328076.504}_{-0.029}^{+0.050}$	${16.969}_{-0.029}^{+0.027}$	$-{10.379}_{-0.036}^{+0.058}$	${0.595}_{-0.025}^{+0.020}$
J0335+4555	${939470.468}_{-0.010}^{+0.073}$	$-{746024.523}_{-0.058}^{+0.066}$	$-{3.638}_{-0.073}^{+0.023}$	$-{0.097}_{-0.105}^{+0.134}$	${0.409}_{-0.027}^{+0.022}$
J0357+5236	$-{61073.169}_{-0.028}^{+0.128}$	$-{744690.138}_{-0.089}^{+0.062}$	${13.908}_{-0.115}^{+0.062}$	$-{10.633}_{-0.058}^{+0.098}$	${0.305}_{-0.077}^{+0.029}$
J0406+6138	$-{41240.544}_{-0.171}^{+0.050}$	${1377882.120}_{-0.093}^{+0.122}$	${12.400}_{-0.085}^{+0.151}$	${22.716}_{-0.060}^{+0.100}$	${0.218}_{-0.057}^{+0.051}$
J0601–0527	$-{765952.497}_{-0.014}^{+0.074}$	${606395.795}_{-0.144}^{+0.097}$	$-{7.348}_{-0.077}^{+0.053}$	$-{15.227}_{-0.105}^{+0.084}$	${0.478}_{-0.045}^{+0.039}$
J0614+2229	${81743.228}_{-0.021}^{+0.050}$	${471984.943}_{-0.014}^{+0.062}$	$-{0.233}_{-0.053}^{+0.036}$	$-{1.224}_{-0.065}^{+0.011}$	${0.282}_{-0.031}^{+0.022}$
J0629+2415	$-{54102.959}_{-0.041}^{+0.179}$	${1069457.665}_{-0.047}^{+0.054}$	${3.629}_{-0.193}^{+0.050}$	$-{4.607}_{-0.153}^{+0.013}$	${0.333}_{-0.054}^{+0.036}$
J0729–1836	${862898.001}_{-0.043}^{+0.143}$	$-{875533.906}_{-0.390}^{+0.400}$	$-{13.072}_{-0.091}^{+0.125}$	${13.252}_{-0.418}^{+0.456}$	${0.489}_{-0.078}^{+0.098}$
J0823+0159	$-{525441.078}_{-0.134}^{+0.150}$	$-{225189.223}_{-0.095}^{+0.080}$	$-{3.797}_{-0.415}^{+0.073}$	${0.171}_{-0.281}^{+0.232}$	${0.376}_{-0.070}^{+0.129}$
J0826+2637	$-{566289.704}_{-0.007}^{+0.007}$	${5345.818}_{-0.012}^{+0.072}$	${62.994}_{-0.007}^{+0.021}$	$-{96.733}_{-0.085}^{+0.045}$	${2.010}_{-0.009}^{+0.013}$
J1022+1001	$-{190717.466}_{-0.016}^{+0.043}$	${25883.813}_{-0.025}^{+0.039}$	$-{14.921}_{-0.033}^{+0.050}$	${5.611}_{-0.035}^{+0.033}$	${1.387}_{-0.028}^{+0.041}$
J1136+1551	$-{87829.834}_{-0.016}^{+0.010}$	$-{73945.169}_{-0.031}^{+0.047}$	$-{73.785}_{-0.010}^{+0.031}$	${366.569}_{-0.055}^{+0.072}$	${2.687}_{-0.016}^{+0.018}$
J1257–1027	$-{135426.468}_{-0.059}^{+0.059}$	$-{182136.895}_{-0.066}^{+0.087}$	$-{7.164}_{-0.105}^{+0.140}$	${12.079}_{-0.110}^{+0.119}$	${0.141}_{-0.092}^{+0.064}$
J1321+8323	${36.330}_{-0.070}^{+0.115}$	${446008.090}_{-0.130}^{+0.173}$	$-{52.674}_{-0.076}^{+0.099}$	${32.373}_{-0.048}^{+0.204}$	${0.968}_{-0.140}^{+0.036}$
J1532+2745	$-{1068597.478}_{-0.029}^{+0.106}$	${646805.042}_{-0.118}^{+0.113}$	${1.542}_{-0.127}^{+0.082}$	${18.932}_{-0.118}^{+0.104}$	${0.624}_{-0.096}^{+0.031}$
J1543–0620	$-{690701.998}_{-0.008}^{+0.041}$	${1328301.506}_{-0.113}^{+0.120}$	$-{16.774}_{-0.063}^{+0.026}$	$-{0.312}_{-0.114}^{+0.147}$	${0.322}_{-0.045}^{+0.028}$
J1607–0032	${1476474.720}_{-0.014}^{+0.060}$	$-{94965.170}_{-0.266}^{+0.163}$	$-{26.437}_{-0.099}^{+0.027}$	$-{27.505}_{-0.200}^{+0.222}$	${0.934}_{-0.047}^{+0.026}$
J1623–0908	$-{1100560.359}_{-0.043}^{+0.175}$	$-{352982.390}_{-0.097}^{+0.182}$	$-{10.769}_{-0.120}^{+0.131}$	${23.509}_{-0.069}^{+0.166}$	${0.586}_{-0.099}^{+0.101}$
J1645–0317	${770570.904}_{-0.016}^{+0.016}$	$-{268070.730}_{-0.136}^{+0.107}$	$-{1.011}_{-0.051}^{+0.003}$	${20.523}_{-0.205}^{+0.147}$	${0.252}_{-0.019}^{+0.028}$
J1650–1654	${172228.769}_{-0.006}^{+0.061}$	${168064.225}_{-0.056}^{+0.097}$	$-{15.024}_{-0.092}^{+0.002}$	$-{6.556}_{-0.131}^{+0.148}$	$-{0.089}_{-0.015}^{+0.031}$
J1703–1846	$-{717146.582}_{-0.025}^{+0.047}$	${712491.475}_{-0.194}^{+0.182}$	$-{0.751}_{-0.056}^{+0.102}$	${16.962}_{-0.230}^{+0.146}$	${0.348}_{-0.047}^{+0.049}$
J1735–0724	${61924.490}_{-0.043}^{+0.076}$	${508913.404}_{-0.093}^{+0.087}$	${0.791}_{-0.029}^{+0.087}$	${20.614}_{-0.046}^{+0.074}$	${0.150}_{-0.035}^{+0.041}$
J1741–0840	${1190194.661}_{-0.069}^{+0.112}$	$-{560681.159}_{-0.087}^{+0.066}$	${0.436}_{-0.126}^{+0.082}$	${6.876}_{-0.066}^{+0.109}$	${0.279}_{-0.058}^{+0.050}$
J1754+5201	$-{341736.402}_{-0.010}^{+0.074}$	$-{1983.221}_{-0.083}^{+0.056}$	$-{3.950}_{-0.046}^{+0.047}$	${1.101}_{-0.059}^{+0.072}$	${0.160}_{-0.022}^{+0.029}$
J1820–0427	$-{169381.766}_{-0.027}^{+0.076}$	$-{63903.822}_{-0.109}^{+0.120}$	$-{7.318}_{-0.055}^{+0.074}$	${15.883}_{-0.069}^{+0.088}$	${0.351}_{-0.055}^{+0.049}$
J1833–0338	${269293.573}_{-0.060}^{+0.119}$	$-{932810.642}_{-0.260}^{+0.182}$	$-{17.409}_{-0.025}^{+0.158}$	${15.038}_{-0.337}^{+0.333}$	${0.408}_{-0.067}^{+0.050}$
J1840+5640	${945467.477}_{-0.011}^{+0.053}$	$-{260360.747}_{-0.029}^{+0.070}$	$-{31.212}_{-0.022}^{+0.033}$	$-{29.079}_{-0.082}^{+0.047}$	${0.657}_{-0.008}^{+0.065}$
J1901–0906	$-{880884.264}_{-0.014}^{+0.071}$	$-{543910.907}_{-0.163}^{+0.113}$	$-{7.531}_{-0.045}^{+0.034}$	$-{18.211}_{-0.159}^{+0.143}$	${0.510}_{-0.042}^{+0.067}$
J1912+2104	$-{170443.421}_{-0.006}^{+0.121}$	$-{180675.561}_{-0.085}^{+0.070}$	$-{11.335}_{-0.097}^{+0.023}$	$-{5.768}_{-0.122}^{+0.092}$	${0.024}_{-0.022}^{+0.171}$
J1913+1400	$-{6337.234}_{-0.018}^{+0.050}$	$-{122121.054}_{-0.033}^{+0.047}$	$-{5.265}_{-0.072}^{+0.040}$	$-{8.927}_{-0.065}^{+0.038}$	${0.185}_{-0.023}^{+0.027}$
J1917+1353	${316294.074}_{-0.016}^{+0.016}$	$-{672694.400}_{-0.025}^{+0.039}$	$-{1.253}_{-0.074}^{+0.022}$	${3.811}_{-0.064}^{+0.057}$	${0.142}_{-0.007}^{+0.068}$
J1919+0021	${893781.316}_{-0.045}^{+0.091}$	$-{7819.701}_{-0.083}^{+0.060}$	${10.167}_{-0.143}^{+0.029}$	$-{4.713}_{-0.073}^{+0.102}$	${0.166}_{-0.042}^{+0.042}$
J1937+2544	$-{865965.163}_{-0.015}^{+0.039}$	${700670.139}_{-0.039}^{+0.029}$	$-{10.049}_{-0.030}^{+0.042}$	$-{13.055}_{-0.039}^{+0.034}$	${0.318}_{-0.029}^{+0.031}$
J2006–0807	$-{522282.331}_{-0.010}^{+0.076}$	${1163967.750}_{-0.126}^{+0.184}$	$-{6.176}_{-0.070}^{+0.035}$	$-{10.616}_{-0.123}^{+0.174}$	${0.424}_{-0.101}^{+0.010}$
J2010–1323	$-{225966.545}_{-0.205}^{+0.305}$	${1202928.433}_{-0.146}^{+0.188}$	${2.358}_{-0.210}^{+0.329}$	$-{5.611}_{-0.303}^{+0.257}$	${0.484}_{-0.120}^{+0.166}$
J2046+1540	${777607.703}_{-0.016}^{+0.099}$	$-{413782.623}_{-0.072}^{+0.041}$	$-{10.455}_{-0.090}^{+0.032}$	${0.681}_{-0.090}^{+0.039}$	${0.310}_{-0.076}^{+0.082}$
J2046–0421	${357981.097}_{-0.037}^{+0.006}$	${848678.492}_{-0.268}^{+0.105}$	${10.760}_{-0.035}^{+0.038}$	$-{4.404}_{-0.076}^{+0.373}$	${0.167}_{-0.042}^{+0.026}$
J2113+2754	$-{109967.706}_{-0.005}^{+0.033}$	${238352.997}_{-0.058}^{+0.047}$	$-{27.981}_{-0.014}^{+0.052}$	$-{54.432}_{-0.096}^{+0.040}$	${0.704}_{-0.022}^{+0.023}$
J2113+4644	$-{704691.944}_{-0.025}^{+0.078}$	${569536.867}_{-0.029}^{+0.066}$	${9.525}_{-0.148}^{+0.068}$	${8.846}_{-0.090}^{+0.076}$	${0.454}_{-0.065}^{+0.077}$
J2145–0750	$-{111529.383}_{-0.016}^{+0.076}$	$-{149818.403}_{-0.080}^{+0.054}$	$-{9.491}_{-0.042}^{+0.052}$	$-{9.114}_{-0.076}^{+0.090}$	${1.603}_{-0.009}^{+0.063}$
J2149+6329	$-{809671.242}_{-0.027}^{+0.101}$	$-{343634.854}_{-0.072}^{+0.184}$	${15.786}_{-0.082}^{+0.131}$	${11.255}_{-0.284}^{+0.092}$	${0.356}_{-0.061}^{+0.072}$
J2150+5247	${1041961.339}_{-0.105}^{+0.176}$	$-{373619.769}_{-0.179}^{+0.163}$	${8.377}_{-0.181}^{+0.226}$	$-{4.427}_{-0.353}^{+0.247}$	${0.034}_{-0.081}^{+0.164}$
J2212+2933	${207796.901}_{-0.052}^{+0.052}$	$-{50710.117}_{-0.083}^{+0.050}$	$-{5.513}_{-0.102}^{+0.067}$	$-{11.322}_{-0.128}^{+0.102}$	${0.265}_{-0.120}^{+0.050}$
J2225+6535	${782263.168}_{-0.065}^{+0.201}$	$-{735020.424}_{-0.058}^{+0.050}$	${147.220}_{-0.223}^{+0.243}$	${126.532}_{-0.115}^{+0.076}$	${1.203}_{-0.204}^{+0.166}$
J2248–0101	${279160.561}_{-0.010}^{+0.151}$	${824047.754}_{-0.083}^{+0.184}$	$-{10.548}_{-0.027}^{+0.117}$	$-{17.407}_{-0.267}^{+0.110}$	${0.256}_{-0.067}^{+0.049}$
J2305+3100	$-{733995.388}_{-0.009}^{+0.042}$	${572881.368}_{-0.068}^{+0.078}$	$-{3.737}_{-0.006}^{+0.082}$	$-{15.571}_{-0.163}^{+0.049}$	${0.223}_{-0.028}^{+0.033}$
J2317+1439	${721838.143}_{-1.295}^{+0.100}$	${259459.178}_{-0.829}^{+0.827}$	$-{1.476}_{-0.065}^{+0.465}$	${3.806}_{-0.704}^{+0.272}$	${0.603}_{-0.241}^{+1.533}$
J2317+2149	${834190.412}_{-0.056}^{+0.040}$	$-{773795.695}_{-0.101}^{+0.093}$	${8.522}_{-0.104}^{+0.035}$	${0.136}_{-0.084}^{+0.192}$	${0.510}_{-0.049}^{+0.057}$
J2325+6316	$-{45049.774}_{-0.012}^{+0.080}$	${16351.365}_{-0.083}^{+0.149}$	$-{5.926}_{-0.073}^{+0.082}$	$-{2.051}_{-0.192}^{+0.188}$	$-{0.010}_{-0.043}^{+0.049}$
J2346–0609	${203020.882}_{-0.006}^{+0.037}$	$-{106887.144}_{-0.052}^{+0.060}$	${37.390}_{-0.042}^{+0.025}$	$-{20.230}_{-0.070}^{+0.107}$	${0.275}_{-0.036}^{+0.021}$
J2354+6155	$-{253697.746}_{-0.005}^{+0.058}$	${1090745.164}_{-0.019}^{+0.035}$	${22.755}_{-0.040}^{+0.056}$	${4.888}_{-0.016}^{+0.033}$	${0.412}_{-0.043}^{+0.031}$

Note. Position offsets are relative to the defined position for the chosen reference source, at epoch MJD 56,000.0; the R.A. offset is calculated at the decl. of the target pulsar.

Download table as:  ASCIITypeset images: 1 2

Detailed results, including astrometric plots and calibrator images, can be found for each pulsar at https://safe.nrao.edu/vlba/psrpi/release.html. As an example, the astrometric plots and bootstrap histograms for PSR J0601−0527 are shown in Figure 4. This is a typical-to-challenging target—the S/N on the target was ∼50, and the in-beam calibrator was slightly resolved, with a total flux density of ∼25 mJy, and separated from the target by 16'. The reduced χ² of the astrometric fit to the position time series using the empirical systematic error estimate discussed in Section 3.2 is 1.3, and the attained parallax precision of ∼40 μas is very close to the median PSRπ value.

Zoom In Zoom Out Reset image size

The fitted parallax and proper motion results can be used to derive distances, Galactic z-heights, and transverse velocities for the target pulsars, or lower limits in the case where the parallax was not measured to sufficient accuracy. Likewise, the fitted offsets from the inbeam calibrators can be combined with an estimate of the in-beam calibrator position uncertainty (comprising contributions from core-shift, phase-referencing to the out-of-beam calibrator, and the out-of-beam calibrator absolute uncertainty added in quadrature as described in Section 2.4), to produce an absolute pulsar position at the reference epoch and associated uncertainty. All of these derived quantities are shown in Table 4. We stress that the absolute positions are of a preliminary nature, since the calibrator positions and core-shifts have not been determined to high precision, and note in particular that the positional uncertainties for PSR J0614+2229, PSR J0629+2415, and PSR J1820–0427 could be substantially underestimated due to the fact that their out-of-beam calibrator source exhibits a compact double structure.

Table 4.  Derived Astrometric Parameters for All PSRπ Targets

Pulsar	R.A.	Decl.	Dist.	z-height	Trans. Vel.
	(J2000)	(J2000)	(kpc)	(kpc)	(km s−1)
J0040+5716	00:40:32.3899(1)	+57:16:24.833(1)	${9.77}_{-3.23}^{+3.13}$	${0.95}_{-0.31}^{+0.30}$	${626.3}_{-221.5}^{+203.2}$
J0055+5117	00:55:45.3981(1)	+51:17:24.601(1)	${2.87}_{-0.39}^{+0.54}$	${0.58}_{-0.08}^{+0.11}$	${275.7}_{-42.5}^{+56.0}$
J0102+6537	01:02:32.9914(1)	+65:37:13.416(1)	${2.51}_{-0.25}^{+0.32}$	${0.12}_{-0.01}^{+0.02}$	${111.9}_{-12.9}^{+16.5}$
J0108+6608	01:08:22.5049(3)	+66:08:34.499(1)	${2.14}_{-0.15}^{+0.15}$	${0.12}_{-0.01}^{+0.01}$	${487.0}_{-35.2}^{+36.4}$
J0147+5922	01:47:44.6434(1)	+59:22:03.284(1)	${2.02}_{-0.16}^{+0.46}$	${0.10}_{-0.01}^{+0.02}$	${71.1}_{-6.6}^{+17.4}$
J0151−0635	01:51:22.7179(2)	−06:35:02.987(2)	${4.60}_{-1.43}^{+2.45}$	${4.17}_{-1.30}^{+2.22}$	${261.0}_{-93.4}^{+116.8}$
J0152−1637	01:52:10.8539(1)	−16:37:53.641(2)	${2.26}_{-0.73}^{+1.56}$	${2.15}_{-0.70}^{+1.48}$	${335.8}_{-124.4}^{+219.4}$
J0157+6212	01:57:49.9434(1)	+62:12:26.648(1)	${1.80}_{-0.12}^{+0.08}$	${0.01}_{-0.00}^{+0.00}$	${383.7}_{-25.8}^{+17.8}$
J0323+3944	03:23:26.6619(1)	+39:44:52.403(1)	${0.95}_{-0.03}^{+0.04}$	${0.24}_{-0.01}^{+0.01}$	${183.2}_{-7.0}^{+7.6}$
J0332+5434	03:32:59.4096(1)	+54:34:43.329(1)	${1.68}_{-0.06}^{+0.07}$	${0.04}_{-0.00}^{+0.00}$	${158.6}_{-5.6}^{+7.5}$
J0335+4555	03:35:16.6416(1)	+45:55:53.452(1)	${2.44}_{-0.12}^{+0.18}$	${0.34}_{-0.02}^{+0.02}$	${42.2}_{-2.4}^{+5.7}$
J0357+5236	03:57:44.8403(2)	+52:36:57.493(1)	${3.27}_{-0.29}^{+1.10}$	${0.03}_{-0.00}^{+0.01}$	${271.7}_{-30.4}^{+93.1}$
J0406+6138	04:06:30.0806(2)	+61:38:41.408(1)	${4.58}_{-0.87}^{+1.63}$	${0.56}_{-0.11}^{+0.20}$	${562.8}_{-107.7}^{+212.6}$
J0601−0527	06:01:58.9752(2)	−05:27:50.871(4)	${2.09}_{-0.16}^{+0.22}$	${0.49}_{-0.04}^{+0.05}$	${167.8}_{-13.4}^{+23.8}$
J0614+2229	06:14:17.0058(1)	+22:29:56.848(1)	${3.55}_{-0.26}^{+0.44}$	${0.15}_{-0.01}^{+0.02}$	${21.0}_{-2.3}^{+4.3}$
J0629+2415	06:29:05.7273(1)	+24:15:41.546(1)	${3.00}_{-0.29}^{+0.57}$	${0.33}_{-0.03}^{+0.06}$	${84.0}_{-11.6}^{+19.0}$
J0729−1836	07:29:32.3369(1)	−18:36:42.244(2)	${2.04}_{-0.34}^{+0.39}$	${0.01}_{-0.00}^{+0.00}$	${179.9}_{-36.5}^{+40.7}$
J0823+0159	08:23:09.7651(1)	+01:59:12.469(1)	${2.66}_{-0.68}^{+0.60}$	${0.96}_{-0.25}^{+0.22}$	${47.8}_{-11.8}^{+18.9}$
J0826+2637	08:26:51.5068(1)	+26:37:21.297(1)	${0.50}_{-0.00}^{+0.00}$	${0.26}_{-0.00}^{+0.00}$	${272.3}_{-1.9}^{+1.3}$
J1022+1001	10:22:57.9957(1)	+10:01:52.765(2)	${0.72}_{-0.02}^{+0.01}$	${0.56}_{-0.02}^{+0.01}$	${54.5}_{-2.0}^{+1.3}$
J1136+1551	11:36:03.1198(1)	+15:51:14.183(1)	${0.37}_{-0.00}^{+0.00}$	${0.35}_{-0.00}^{+0.00}$	${659.7}_{-4.5}^{+4.2}$
J1257−1027	12:57:04.7625(2)	−10:27:05.551(2)	${7.09}_{-2.22}^{+13.18}$	${5.62}_{-1.76}^{+10.44}$	${474.1}_{-160.0}^{+886.4}$
J1321+8323	13:21:45.6315(7)	+83:23:39.432(1)	${1.03}_{-0.04}^{+0.17}$	${0.57}_{-0.02}^{+0.10}$	${302.6}_{-12.0}^{+51.6}$
J1532+2745	15:32:10.3646(1)	+27:45:49.623(1)	${1.60}_{-0.07}^{+0.29}$	${1.31}_{-0.06}^{+0.24}$	${144.0}_{-7.6}^{+28.8}$
J1543−0620	15:43:30.1373(1)	−06:20:45.332(2)	${3.11}_{-0.25}^{+0.51}$	${1.85}_{-0.15}^{+0.30}$	${247.4}_{-22.2}^{+41.3}$
J1607−0032	16:07:12.0598(2)	−00:32:41.527(2)	${1.07}_{-0.03}^{+0.06}$	${0.62}_{-0.02}^{+0.03}$	${193.4}_{-7.0}^{+13.0}$
J1623−0908	16:23:17.6599(1)	−09:08:48.733(2)	${1.71}_{-0.25}^{+0.34}$	${0.78}_{-0.12}^{+0.16}$	${209.4}_{-34.3}^{+44.0}$
J1645−0317	16:45:02.0406(1)	−03:17:57.819(2)	${3.97}_{-0.39}^{+0.33}$	${1.74}_{-0.17}^{+0.14}$	${386.4}_{-43.2}^{+38.7}$
J1650−1654	16:50:27.1694(7)	−16:54:42.282(20)	>3.1	>0.9	>229.8
J1703−1846	17:03:51.0915(2)	−18:46:14.845(6)	${2.88}_{-0.36}^{+0.45}$	${0.67}_{-0.08}^{+0.11}$	${231.5}_{-38.6}^{+40.4}$
J1735−0724	17:35:04.9730(1)	−07:24:52.130(1)	${6.68}_{-1.43}^{+2.03}$	${1.53}_{-0.33}^{+0.47}$	${653.1}_{-144.8}^{+205.8}$
J1741−0840	17:41:22.5629(1)	−08:40:31.711(1)	${3.58}_{-0.55}^{+0.94}$	${0.70}_{-0.11}^{+0.18}$	${116.7}_{-21.8}^{+35.7}$
J1754+5201	17:54:22.9068(1)	+52:01:12.244(1)	${6.27}_{-0.98}^{+1.03}$	${3.10}_{-0.48}^{+0.51}$	${122.0}_{-22.2}^{+22.4}$
J1820−0427	18:20:52.5934(1)	−04:27:37.712(2)	${2.85}_{-0.35}^{+0.52}$	${0.23}_{-0.03}^{+0.04}$	${236.0}_{-30.2}^{+47.0}$
J1833−0338	18:33:41.8945(1)	−03:39:04.258(1)	${2.45}_{-0.27}^{+0.48}$	${0.10}_{-0.01}^{+0.02}$	${266.7}_{-34.0}^{+56.8}$
J1840+5640	18:40:44.5372(1)	+56:40:54.852(1)	${1.52}_{-0.14}^{+0.02}$	${0.62}_{-0.06}^{+0.01}$	${307.9}_{-27.8}^{+5.2}$
J1901−0906	19:01:53.0087(3)	−09:06:11.146(10)	${1.96}_{-0.23}^{+0.17}$	${0.22}_{-0.03}^{+0.02}$	${183.1}_{-24.3}^{+20.4}$
J1912+2104	19:12:43.3391(1)	+21:04:33.926(1)	${41.02}_{-35.90}^{+377.75}$	${3.57}_{-3.12}^{+32.86}$	${2482.4}_{-2184.1}^{+6732.3}$
J1913+1400	19:13:24.3527(1)	+14:00:52.559(1)	${5.42}_{-0.70}^{+0.75}$	${0.15}_{-0.02}^{+0.02}$	${266.3}_{-37.6}^{+39.5}$
J1917+1353	19:17:39.7864(1)	+13:53:57.077(1)	${7.04}_{-2.29}^{+0.38}$	${0.08}_{-0.02}^{+0.00}$	${134.5}_{-47.8}^{+9.7}$
J1919+0021	19:19:50.6715(1)	+00:21:39.722(2)	${6.03}_{-1.21}^{+2.02}$	${0.65}_{-0.13}^{+0.22}$	${319.6}_{-72.3}^{+107.2}$
J1937+2544	19:37:01.2544(1)	+25:44:13.436(1)	${3.15}_{-0.28}^{+0.32}$	${0.12}_{-0.01}^{+0.01}$	${245.9}_{-24.4}^{+25.1}$
J2006−0807	20:06:16.3650(1)	−08:07:02.167(3)	${2.36}_{-0.06}^{+0.73}$	${0.82}_{-0.02}^{+0.25}$	${137.1}_{-8.5}^{+45.4}$
J2010−1323	20:10:45.9211(1)	−13:23:56.083(4)	${2.07}_{-0.53}^{+0.68}$	${0.83}_{-0.21}^{+0.27}$	${58.7}_{-16.0}^{+26.9}$
J2046+1540	20:46:39.3373(1)	+15:40:33.558(1)	${3.22}_{-0.68}^{+1.04}$	${0.93}_{-0.20}^{+0.30}$	${160.0}_{-35.4}^{+48.2}$
J2046−0421	20:46:00.1730(1)	−04:21:26.256(2)	${5.98}_{-0.81}^{+2.00}$	${2.75}_{-0.37}^{+0.92}$	${329.3}_{-50.6}^{+112.3}$
J2113+2754	21:13:04.3506(1)	+27:54:01.160(1)	${1.42}_{-0.04}^{+0.04}$	${0.34}_{-0.01}^{+0.01}$	${412.1}_{-13.8}^{+14.2}$
J2113+4644	21:13:24.3295(1)	+46:44:08.844(1)	${2.20}_{-0.32}^{+0.36}$	${0.05}_{-0.01}^{+0.01}$	${135.8}_{-24.8}^{+18.4}$
J2145−0750	21:45:50.4588(1)	−07:50:18.514(4)	${0.62}_{-0.02}^{+0.00}$	${0.42}_{-0.02}^{+0.00}$	${38.9}_{-1.9}^{+0.5}$
J2149+6329	21:49:58.7033(2)	+63:29:44.277(2)	${2.81}_{-0.47}^{+0.58}$	${0.36}_{-0.06}^{+0.07}$	${258.3}_{-49.4}^{+52.6}$
J2150+5247	21:50:37.7499(1)	+52:47:49.556(1)	>2.4	>0.0	>89.1
J2212+2933	22:12:23.3444(1)	+29:33:05.411(1)	${3.77}_{-0.60}^{+3.14}$	${1.40}_{-0.22}^{+1.16}$	${225.7}_{-38.2}^{+193.8}$
J2225+6535	22:25:52.8627(3)	+65:35:36.371(1)	${0.83}_{-0.10}^{+0.17}$	${0.10}_{-0.01}^{+0.02}$	${765.2}_{-94.5}^{+157.6}$
J2248−0101	22:48:26.8859(1)	−01:01:48.085(1)	${3.90}_{-0.63}^{+1.40}$	${3.02}_{-0.49}^{+1.08}$	${377.1}_{-66.7}^{+149.2}$
J2305+3100	23:05:58.3212(1)	+31:00:01.281(1)	${4.47}_{-0.58}^{+0.65}$	${2.01}_{-0.26}^{+0.29}$	${341.1}_{-47.2}^{+56.4}$
J2317+1439	23:17:09.2364(1)	+14:39:31.265(1)	${1.66}_{-1.19}^{+1.10}$	${1.12}_{-0.80}^{+0.74}$	${31.4}_{-22.9}^{+30.0}$
J2317+2149	23:17:57.8419(1)	+21:49:48.019(1)	${1.96}_{-0.20}^{+0.21}$	${1.16}_{-0.12}^{+0.12}$	${79.3}_{-9.6}^{+8.9}$
J2325+6316	23:25:13.3196(2)	+63:16:52.362(1)	>12.1	>0.4	>327.8
J2346−0609	23:46:50.4978(1)	−06:09:59.899(2)	${3.64}_{-0.26}^{+0.55}$	${3.27}_{-0.23}^{+0.50}$	${732.5}_{-53.1}^{+114.5}$
J2354+6155	23:54:04.7830(1)	+61:55:46.845(1)	${2.42}_{-0.17}^{+0.28}$	${0.01}_{-0.00}^{+0.00}$	${268.0}_{-20.3}^{+32.7}$

Note. Positions are shown at the reference epoch of MJD 56,000.0. 1σ uncertainties are shown on measured values (parentheses show uncertainties on the last digit), and lower limits are shown with 95% confidence.

Download table as:  ASCIITypeset images: 1 2

The most probable distance and the 68% confidence interval were calculated directly from the fitted parallax and confidence interval, without applying any priors based on an assumed pulsar spatial distribution or luminosity distribution (e.g., Verbiest et al. 2012; Igoshev et al. 2016). For high-significance parallax detections, the distance is relatively insensitive to the assumed priors, but we note that for low significance parallax detections (for instance, the 20 PSRπ pulsars with a parallax significance below 5σ), the inferred distance can be substantially dependent on the assumed priors. The most probable transverse velocity was estimated using the most probable distance and most probable proper motion, while the 68% confidence interval was calculated by finding the smallest rectangular cuboid in (parallax, proper motion [R.A.], proper motion [decl.]) space that encompassed 68% of the bootstrap trial results, and taking the highest and lowest transverse velocity from these included trials. Figure 5 shows an example of the transverse velocity estimator for PSR J0601−0527.

Zoom In Zoom Out Reset image size

As shown previously, correctly estimating the total uncertainty of the position measurements used for the astrometric fit is challenging. The following list summarizes the primary contributions to the error budget:

• 1.  
  Thermal noise in the target image. This is the most readily quantified, as it can be easily extracted from the image-plane fitting. It is inversely proportional to the instrumental resolution and inversely proportional to the signal-to-noise in the pulsar image. This term generally dominates for faint sources.
• 2.  
  Systematic offsets introduced by differential propagation effects between the target and calibrator. This is usually the dominant term for bright sources, where the signal-to-noise on the target is not the limiting factor. At 1600 MHz, the ionosphere dominates these path length differences, which vary on a sub-epoch timescale (minutes to hours). The solutions on the calibrator must be extrapolated spatially to the target, and are averaged over a time interval during which the ionosphere can change. Generally, the spatial extrapolation introduces the largest error, meaning this term is most dependent on the calibrator-target separation, along with factors influencing the mean path length through the ionosphere such as the solar activity level, time of day, and antenna elevation. The amount of temporal averaging required depends on the calibrator flux density and structure—for a typical 20 mJy calibrator, a solution interval of order 1.25 minutes was typical.Refraction in the interstellar medium also produces a differential offset between a pulsar and calibrator. For the PSRπ sample, the predicted angular wandering (using the predictions of the NE2001 model and assuming a Kolmogorov density spectrum; Cordes & Lazio 2002) at the observing frequency (1660 or 2270 MHz) for our target pulsars due to refraction has a median value for a given observation of ∼0.05 mas. For most pulsars, this is a negligible component of the astrometric error. However, the predicted scattering disk diameter exceeds 1 mas and the predicted refractive wander exceeds 0.1 mas for four pulsars in our sample: PSR J0601−0527, PSR J1833−0338, PSR J2212+2933, and PSR J2325+6313. In cases such as these, and others where the astrometric precision is extremely high, refractive wander of the pulsar may be a significant component of the error budget. In addition to the pulsar, the refractive wander also affects the in-beam calibrator sources, which is another potential source of error in the target–calibrator separation. Generally, the refractive wander is larger for calibrators than for the pulsars, as the radiation from the calibrators passes by all of the Galactic electrons along the line of sight, leading to a larger scattering disk. For PSR J1833−0338, for instance, the NE2001 model predicts refractive wander with an rms deviation of 0.4 mas for the position reference calibrator source, nearly three times larger than that of the pulsar. However, the refractive wander timescale for the calibrator sources is typically much longer (years), meaning that reference position will be affected more than proper motion, which will itself be affected more than parallax.
• 3.  
  Systematic variations in the image reference frame. An imperfect model of the calibrator source will lead to an offset in the obtained pulsar position. If this is constant in time, it does not impact the measurement of parallax or proper motion, but time variability is an important source of error. Time variability could be the result of evolution intrinsic to the source itself (which is present, at least at a low level, in all compact sources), or from changes in the observing setup (different frequency or uv coverage between observing epochs). This can be the dominant term if the calibrator source is bright and close to the target (minimizing the ionospheric terms) but is a blazar-like source that displays large and rapid variations in the jet structure. Over the 1–2 yr timescale typical for pulsar astrometry programs such as PSRπ, it is often possible to fit a significant component of this reference source offset with a linear function with time, meaning it can corrupt the proper motion measured for the pulsar. However, for most reference sources the likely effect is small compared to our measurement error (the median apparent proper motion seen by Moór et al. 2011 was 19 μas yr⁻¹ versus 106 μas yr⁻¹ for our relative astrometric uncertainty), and parallax (which has a sinusoidal signature with time) is much less affected. As well as effects intrinsic to the source, time-variable position shifts due to the changing Galactic gravitational potential field can be expected (Larchenkova et al. 2017), albeit only at the level of up to ∼10 μas over our timescales, and hence smaller than the reference source structure effects.
• 4.  
  Stochastic noise in the image reference frame. The phase solutions on the calibrator source will have some noise dependent on the signal-to-noise in each solution interval, which is determined by the source flux density, instrumental sensitivity, and calibration interval. Fainter calibrators will generally lead to an increase in this term; although this could be compensated by increasing the solution interval, that would be reduce the ability to compensate for time-variable ionospheric effects. Changing the solution interval can thus shift error between this term and term 2 of this list; we seek to choose a value on a per-source basis that minimizes their sum.

Accordingly, there are five main factors we would expect to influence the total uncertainty of a position measurement in a given epoch:

• 1.  
  pulsar flux density;
• 2.  
  calibrator flux density;
• 3.  
  pulsar-calibrator angular separation;
• 4.  
  ionospheric conditions and average observing elevation; and
• 5.  
  calibrator stability.

We can straightforwardly measure all but the last of these factors, and as noted the calibrator stability is not usually expected to be a significant contributor to the total error budget. We undertook a number of approaches to try and determine the contributions of each of the remaining factors to our error budget. To illustrate the results, we again use the typical-to-challenging source PSR J0601−0527, where the thermal noise errors alone substantially underestimate the total error budget, as can be seen from the reduced ${\chi }^{2}\sim 15$ for a simple least-squares fit to the unmodified position data.

First, we examined the apparent shifts in pulsar position within an observation, by subdividing the data into two halves and imaging each one separately. This approach is particularly useful for identifying observations where short-term ionospheric conditions were unstable—if the two positions differ by considerably more than their formal error bars, it is likely that the mean position over the whole observation also has an underestimated position uncertainty. However, while reliable, this approach is likely not complete, as it would fail to pick up observations with large but relatively stable residual ionospheric "wedges" that lead to a fairly constant position offset over the whole observation duration. Also, when the target source is faint, making a significant measurement of the offset between the two halves of the observation may not be possible.

We evaluated an approach in which we recorded the minimum systematic offset between the two observation halves (accounting for the uncertainty in the position measurements) and set the systematic error contribution to the whole epoch to be half of this value. The R.A. and decl. axes are treated separately. The effect, as expected, was to lower the chi-squared of the resultant astrometric fit, although the position errors remained underestimated (as determined by a χ² value well in excess of 1.0) in many cases. Figure 6 shows the result for PSR J0601−0527. As expected, this approach yields a systematic error estimate that is too low; the reduced χ² remains at 10. Using this refined position set as the input for bootstrap fits results in a change in the best-fit value at the 1σ level, and gives a small reduction in the estimated parameter uncertainty. For most targets, the impact on both best-fit value and uncertainty was smaller than in this example.

Zoom In Zoom Out Reset image size

Second, we processed the data sets multiple times, making use of different ionospheric models, and examined the resultant position scatter for each epoch. If different ionospheric models of comparable quality give widely divergent positions, then the residual error from our chosen ionospheric model is likely high, since we have no way of determining which of the ionospheric models is correct. As with the previous approach comparing the two halves of an observation in time, this method is likely reliable but not necessarily complete, as it will not pick up cases where all models suffer from the same deficiencies.

We investigated all of the products covering our complete observing timespan from ftp://cddis.gsfc.nasa.gov/gps/products/ionex/, and found that the jplg, codg, igsg, and esag models consistently gave the best results, with the lowest residuals on the astrometric fit. The igsg model was chosen for the final data reduction. We therefore trialed an approach in which, for every epoch, we computed the rms scatter in the positions provided by processing using these four models, and used this as the estimate for systematic error for that epoch. As with the epoch-splitting approach given previously, the R.A. and decl. axes were treated separately. We found that, for many sources, a given ionospheric model yielded a statistically significant mean offset across all epochs in addition to random scatter, and so for all sources we subtracted (per model) the mean positional offset from all epochs before computing the rms.

As expected and as with the previous approach splitting the observation in halves, we typically capture some but not all of the systematic error with this technique. Figure 7 shows the result for PSR J0601−0527; the reduced χ² is still 8. The best-fit value for parallax changes by approximately 0.2σ, and the estimated parallax uncertainty from the bootstrap fit is reduced by ∼10%.

Zoom In Zoom Out Reset image size

Our third approach to estimate systematic error made use of the residual position errors across our entire data set of 60 pulsars. As shown previously, we expect the dominant error sources to depend on the calibrator-target separation, the observing elevation, and the calibrator brightness. As a simplification, we consider that the systematic error should be proportional to two quantities: the mean "deprojected" calibrator-target separation (calculated as the angular separation multiplied by the cosecant of the observing elevation, averaged over all antennas and all scans in the observation) and the S/N achieved on the inbeam calibrator source(s). For each epoch, we added a systematic error estimate given by

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{\mathrm{sys}}=A\times \displaystyle \frac{s\times {\displaystyle \sum }_{a=1}^{N}{\displaystyle \sum }_{o=1}^{M}\mathrm{cosec}\left({\mathrm{el}}_{a,o}\right)}{M\times N}+B/S,\end{eqnarray} \tag{ 1 }$$

where Δ_sys is the systematic error estimate in fractions of a synthesized beam; s is the calibrator-target separation in arcminutes; el_a,o is the observing elevation for antenna a in scan o on the target pulsar; N and M are the number of antennas and target scans, respectively; and S is the S/N on the calibrator source (added in quadrature if multiple sources were used). We conducted a brute-force grid search for the optimal values of the coefficients A and B, seeking the values that gave the tightest grouping of reduced χ² values around 1.0 for our ensemble of 60 target pulsars.

The optimal values were found to be A = 0.001, B = 0.6. When this estimate of systematic error is added for all pulsars, the 25%, median, and 75% values of reduced ${\chi }^{2}$ across all pulsars become 0.64, 1.08, and 1.69, compared to 1.78, 4.81, and 11.96, when no estimate of systematic error is added. The spread in reduced χ² values is comparable to that expected, given the typical number of degrees of freedom (∼11) in the astrometric fits. The results for our example pulsar PSR J0601−0527 are shown in Figure 8; the typical systematic error contribution is 110 μas in R.A. and 280 μas in decl. at each epoch, and the revised reduced χ² is 1.35.

Zoom In Zoom Out Reset image size

We summarize the results of our different estimates of systematic error for PSR J0601−0527 in Table 5. Three things are immediately apparent:

• 1.  
  The inclusion of systematic error (however estimated) pushes the estimated parallax in one direction. When the position errors are severely underestimated (as they are initially), individual discrepant epochs can have overly large effects on the fit. This diminishes once a more realistic error is applied.
• 2.  
  The fitted parameters and their uncertainties remain relatively unchanged regardless of the systematic error estimate used when estimated using a bootstrap.
• 3.  
  The fitted parameters and their uncertainties exhibit good agreement between a bootstrap and a simple least-squares fit when the systematic errors are reasonably well estimated (as appears to be the case for the empirically estimated values).

For our quoted results, we choose to use the bootstrap fit to the data set including empirically estimated systematic errors, which is generally the most conservative (and we believe) correct error estimate we can make with our available information. We do, however, note that in some cases, generally when the parallax has been poorly sampled due to non-detections, this bootstrap error estimate may be overly conservative (because many trials have effectively no sensitivity to parallax). In these cases, better constraints could be obtained by interpreting, with caution, the least-squares fit to the data set incorporating empirically estimated systematic errors. We highlight this for individual pulsars in the discussion that follows in Section 4.1.

Finally, since the unmodeled ionosphere dominates the error budget in many cases, we should expect that the level of solar activity should significantly impact the results obtained. While the solar cycle peaking in 2013 was not particularly active by historical standards, our observations were nevertheless made near the solar maximum, and accordingly we would expect that the same observations repeated a half-decade later would yield better results. Likewise, the precise results seen in, for example, Deller et al. (2012, 2013) might have been more difficult to obtain at the time of the observations presented here. Importantly, our empirical estimates of systematic error should be used with caution when applied to observations in different observing conditions. An even larger observing program might consider including a measure of ionospheric activity as a parameter in the empirical fit to account for this.

From our sample of 60 pulsars, three sources display discrepancies that indicate a potentially biased parallax estimation. We consider the results for these sources in detail and estimate the probability that any of the remaining 57 sources have comparable but undetected errors.

4.1.1. PSR J1650–1654

PSR J1650–1654 has a significant negative parallax of $-{0.089}_{-0.019}^{+0.030}\,\mathrm{mas}$ (where the uncertainty denotes the 68% confidence interval from the bootstrap fit). Since a negative parallax is unphysical, this indicates that the obtained value is incorrect by at least 3σ, but possibly more as the NE2001 distance obtained using the pulsar's DM is just 1.5 kpc.

Two in-beam calibration sources are present, allowing us to cross-check results. The relative position separation of the in-beam phase calibrator source J165133–170928 from the second in-beam source J165015–165730 displays quite a large scatter (0.2 mas rms in R.A., 0.25 mas rms in decl., as seen in Figure 9), and an astrometric fit to this relative separation time series gives a best-fit parallax of 0.05 mas with an uncertainty of ∼0.1 mas. This is not unexpected, given the differential ionosphere across the ∼20' separation. However, as shown in Figure 9, the position reference J165015–165730 is separated by just 4' from the pulsar, meaning that differential ionospheric effects should be smaller than between the two calibrators. An alternative is that the position reference source J165015–165730 exhibits substantial structure variations that lead to position offsets up to a few tenths of a milliarcsecond between epochs, and a substantial amount of the power this introduces into the pulsar's position time series can be fit by the parallax term. The structure of the two calibrator sources can be seen in Figure 9—both exhibit prominent milliarcsecond-scale jets, meaning structure evolution is likely at some level.

Zoom In Zoom Out Reset image size

At present, we have no definitive explanation for the discrepant parallax for PSR J1650–1654, but the results indicate that the parallax uncertainty is underestimated by a factor of several for this pulsar.

4.1.2. PSR J1820–0427

PSR J1820–0427 has two in-beam calibrator sources (J182043–042412 and J182103–042633), which are close to each other and to the pulsar on the sky (angular separations 3'–5'). PSR J1820–0427 was observed at 2.3 GHz due to the strong scattering along this low Galactic latitude line of sight.

Figure 10 shows the clear structure in the residuals between the two in-beam calibrators. The large offset in the first observation is likely due to the different observing frequency (1650 MHz versus 2360 MHz), which would result in a different observed source structure. Neglecting this, however, it is clear that a significant parallax is measured between the two calibrators. The higher observing frequency, along with the small angular separation, should have given relatively small ionospheric errors (which in any case should not yield a parallax-like signature). Both in-beam calibrators do, however, exhibit significant milliarcsecond-scale jet structure, and as with the PSR J1650–1654 calibrators, evolution in this structure over time is a potential source of error. PSR J1820–0427 is also at low Galactic latitude and the NE2001 predictions for scattering disk size and refractive wander of the in-beam calibrators are among the largest in the PSRπ sample, meaning refractive effects may also be a contributing factor.

Zoom In Zoom Out Reset image size

The possibility that one of the two in-beam calibrators is actually Galactic cannot be definitively ruled out, but both are relatively bright (tens of mJy) and stable, arguing strongly against possibilities such as a microquasar or magnetically active protostar.

Since the discrepancy between the two calibrator sources is larger than parallax uncertainty calculated for the pulsar, there is a substantial possibility that the distance uncertainty of PSR J1820–0427 is underestimated. If J182043–042412 had been used as a position reference instead of J182103–042633, the best–fit distance would change from ∼3 to ∼9 kpc. Arguing in favor of using J182103–042633 is the fact that it is the closer source to the pulsar on the sky; this is the criterion applied throughout the PSRπ sample. And since the pulsar–calibrator separation is smaller than the calibrator–calibrator separation, if the offsets are due to a differential term such as the ionosphere (as opposed to an offset created by one or the other calibrator source that is independent of angular separation), then the effect on the pulsar should be smaller than those seen here between the two calibrators.

4.1.3. PSR J2325+6316

4.1.4. Implications for the Remainder of the PSR$\pi$ Sample

4.1.5. PSR J2317+1439

4.1.6. Pulsars with Previous VLBI Astrometry

Three pulsars from the PSRπ sample have previously been the subject of VLBI astrometry: PSR J0332+5434 and PSR J1136+1551 (Brisken et al. 2002), and PSR J0826+2637 (Gwinn et al. 1986). Table 6 shows the PSRπ results for these pulsars compared against the previous results, for proper motion in R.A. and decl. (μ_α, μ_δ) measured in mas yr⁻¹ and parallax in mas. While agreement is good in most cases, with five of the nine measured parameters agreeing to better than 1σ, three values have a discrepancy exceeding 2σ, which is not expected statistically.

Table 6.  Comparison against Previous VLBI Astrometry

Astrometric Quantity	PSR J0332+5434	PSR J0826+2637	PSR J1136+1551
μα (PSRπ)	${16.969}_{-0.029}^{+0.027}$	${62.994}_{-0.007}^{+0.021}$	$-{73.785}_{-0.010}^{+0.031}$
μα (previous)	17.00 ± 0.27	62.6 ± 2.4	−73.95 ± 0.38
μδ (PSRπ)	$-{10.379}_{-0.036}^{+0.058}$	$-{96.733}_{-0.085}^{+0.045}$	${366.569}_{-0.055}^{+0.072}$
μδ (previous)	−9.48 ± 0.37	−95.3 ± 2.4	368.05 ± 0.28
Parallax (PSRπ)	${0.595}_{-0.025}^{+0.020}$	${2.010}_{-0.009}^{+0.013}$	${2.687}_{-0.016}^{+0.018}$
Parallax (previous)	0.94 ± 0.11	2.8 ± 0.6	2.8 ± 0.16

Download table as:  ASCIITypeset image

For PSR J0332+5434, the parallax and proper motion in decl. measured by Brisken et al. (2002) differ from the more precise PSRπ values by 2.5–3σ. In Brisken et al. (2002), only four position measurements were made for this pulsar, meaning the resultant astrometric fit had only three degrees of freedom and would be susceptible to larger errors induced by poor fits in one or two epochs. The final and most discrepant case is proper motion in decl. for PSR J1136+1551, where the Brisken et al. (2002) value (based on five observations) differs from the more precise PSRπ value by around 5σ. While the shorter timespan (12 months) and small number of observations in the Brisken et al. (2002) program would make it more susceptible to potential biases to proper motion such as calibrator source structure evolution, such a large discrepancy remains difficult to explain.

The PSRπ parallax sample allows estimation of the line-of-sight (LoS) average electron density and will provide important input to the next generation Galactic electron density model.¹⁵

The mean electron density (cm⁻³) for a given LoS is n_e = DM(pc cm⁻³) × parallax(mas)/1000. The mean across the sample is $\langle {n}_{{\rm{e}}}\rangle =0.020\,{\mathrm{cm}}^{-3}$, and the rms value is ${\sigma }_{{n}_{{\rm{e}}}}=0.016\,{\mathrm{cm}}^{-3}$. The median distance ∼2.5 kpc implies that the electron densities are representative of the solar region in the Galaxy. However, the sample includes a few objects that are well above any realistic scale height (≲2 kpc) for the electrons, so the sample-mean density is biased lower than the mid-plane value. Restricting the sample to six objects within 1 kpc of the Sun, we obtain a larger value, $\langle {n}_{{\rm{e}}}\rangle =0.026\,{\mathrm{cm}}^{-3}$.

The distribution of LoS electron densities is positively skewed. Figure 11 shows the probability density function (PDF) for the electron density. It is calculated as the sum of rectangle functions, each centered on the nominal value of n_e and having a width equal to the 68% confidence interval given in Table 3. For the three pulsars with lower limits on the distance, the rectangle extends from n_e = 0 to its upper bound. The PDF is normalized to unit area. Its width reflects the wide variation of electron densities between LoSs caused by Galactic structure on both small and large scales.

Zoom In Zoom Out Reset image size

Figure 12 shows the perpendicular component ${\mathrm{DM}}_{z}\,=\mathrm{DM}\,\sin \,| b|$ plotted against the z distance of each pulsar from the Galactic plane. Dashed lines show the expected values for constant densities of 0.01 and 0.1 cm⁻³. Most of the points are between these two lines except at large values of z, where there appears to be a maximum in DM_z, even though about eight of the PSRπ pulsars extend from 2 to 5 kpc above the plane. The PSRπ sample by itself therefore verifies that the Galaxy must have an electron density component similar to the thick disk component in the NE2001 (Cordes & Lazio 2002) and YMW16 (Yao et al. 2017) models.

Zoom In Zoom Out Reset image size

Distances of PSRπ pulsars projected onto the Galactic plane (Figure 13) show that the sample extends to large distances parallel to the plane and sample several spiral arms (as defined in the NE2001 model).

Zoom In Zoom Out Reset image size

Galactic electron density models are based on a wide variety of measurements but ultimately require independently obtained distances for as many pulsars as possible. The NE2001 model used distance constraints on 112 pulsars, of which only 14 were parallax measurements; the remainder were mostly from H i-absorption constrained distances and from pulsars in globular clusters.¹⁶ For the YMW16 model, 73 parallax measurements were used, many from pulsar timing, but 29% had a constrained distance range (maximum to minimum ratio) of 1.5:1 and 13% had more than a 2:1 range.

The performance of the NE2001 and YMW16 models can be compared against the PSRπ distances as a "blind" test because the PSRπ sample was not used in the construction of either model. The comparison is particularly useful for more distant pulsars; the median distance of the PSRπ sample is 2.5 kpc, while the median distance to pulsars with previously published VLBI parallaxes is 1.1 kpc.

Figure 14 shows the ratio d_π/d_model for the NE2001 and YMW16 models in the top and bottom panels, respectively. Circle sizes indicate values of this ratio while colors indicate approximate distances. Comparison of the two figures indicates that both models show large errors for some objects, with the NE2001 model doing better on some objects and the YMW16 model on others. The YMW16 model performs somewhat better on a few high-latitude pulsars than the NE2001 model. The distances of two pulsars with the most negative latitudes (the red circles between 60° and 90° longitude) are overestimated in the YMW16 model and underestimated by the NE2001 model. The median distance ratio exceeds unity for both models (1.5 and 1.1 for NE2001 and YMW16, respectively), and the rms values of this ratio are 1.1 and 0.8, respectively. These results demonstrate that the new PSRπ sample will be extremely valuable for the next-generation Galactic electron density distribution model.

Zoom In Zoom Out Reset image size

Yao et al. (2017) claim that in 95% of cases, the YMW16 predicted distance d_YMW will fall within the range $0.1\,\times {d}_{\mathrm{actual}},1.9\times {d}_{\mathrm{actual}}$. This claim can be examined using the PSRπ data set. We restrict ourselves to pulsars whose VLBI parallax significance is at least 5σ, of which there are 42. Of these 42 pulsars, 6 have a predicted YMW16 distance that falls outside the range $\left(0.1\times {d}_{{\pi }_{\min }},1.9\times {d}_{{\pi }_{\max }}\right)$, where ${d}_{{\pi }_{\min }}$ and ${d}_{{\pi }_{\max }}$ are the values given by inverting the 95% confidence interval for parallax. This would be inconsistent with 95% of the YMW16 distances falling within the range $\left(0.1\times {d}_{\mathrm{actual}},1.9\times {d}_{\mathrm{actual}}\right)$, if the PSRπ sample was representative of the entire pulsar population. However, the PSRπ sample is explicitly not an unbiased sample of pulsars, and in particular, pulsars at high Galactic latitudes are intentionally over-represented in order to help constrain the Galactic scale height. Of the six discrepant YMW16 predictions, five are at moderate to high Galactic latitudes ($| b| \gt 20^\circ$), and the YMW16 model places them beyond the edge of the Galaxy, while the median Galactic latitude of the 42 pulsars with a significant parallax distance is 14°. The sixth source, at b = 113, is underpredicted by an order of magnitude. Accordingly, based on the PSRπ sample, we advise that all DM-based distance estimates be used with caution, especially for high Galactic latitude pulsars.

Finally, as noted in Section 3.2, refractive wander is potentially a significant contributor to the differential astrometry error budget for some lines of sight, particularly at low Galactic latitude. While beyond the scope of this work, future analysis using the PSRπ sample could refine the scattering disk and refractive wander predictions of the NE2001 model, leading to improved estimates of systematic error contributions for future studies.

Galactic pulsars have a larger scale height than their progenitor massive stars, leading to the early inference that they have high velocities (Gunn & Ostriker 1970). Individual high velocity objects such as PSR B1508+55 (v_⊥ ∼ 1000 km s⁻¹; Chatterjee et al. 2005) and PSR J2225+6535, the Guitar Nebula pulsar (v_⊥ ∼ 800 km s⁻¹; this work; Cordes et al. 1993), establish stringent constraints on natal kicks and the minimum asymmetry requirements in simulations of supernova core collapse (e.g., Fryer 2004), and the overall population velocities inform models for neutron star birth, supernova explosions, and the evolution of close binary systems. The pulsar velocity distribution has thus been a topic of continued interest (e.g., Lyne & Lorimer 1994; Hansen & Phinney 1997; Cordes & Chernoff 1998; Arzoumanian et al. 2002; Hobbs et al. 2005; Verbunt et al. 2017).

The parallax and proper motion measurements presented here provide model-independent estimates of pulsar distances and transverse velocities (Table 4) and thus mitigate a key uncertainty in deriving the pulsar velocity distribution. We note, however, that the astrophysically relevant quantity is the three-dimensional birth velocity for the entire pulsar population. For individual pulsars, their uncertain age limits the precision of any extrapolation in the Galactic gravitational potential, and their radial velocity is unknown, rendering a full three-dimensional birth velocity unmeasurable. While such uncertainties can be addressed statistically, the inference of population parameters is further affected by biases in the sample of objects with precise astrometry (see, e.g., Arzoumanian et al. 2002). The targets in this work were selected based on flux density and calibrator availability (Section 2.1), and inferring the properties of the population would require the addressing of selection effects in the original detection surveys. Such selection effects are not trivial: for example, many pulsar surveys focus on the Galactic plane where the stellar progenitors of pulsars are concentrated,¹⁷ but on average, the higher velocity neutron stars will spend less time near the Galactic plane compared to the lower velocity objects. Thus the high velocity tail of the pulsar population will be suppressed in a typical survey yield.

We defer the treatment of the full range of selection effects to future work, but here we present histograms of the measured transverse velocities for each pulsar in our sample, along with a normalized probability distribution summing across all pulsars (Figure 15). Note that the four millisecond (recycled) pulsars in our sample (PSRs J1022+1001, J2010−1323, J2145−0750, and J2317+1439) are excluded from the cumulative distribution, since recycled pulsars are an older population with a lower characteristic velocity distribution (e.g., Cordes & Chernoff 1997). As a comparison for our sample of transverse velocities, we plot a recent velocity distribution model for young pulsars (Verbunt et al. 2017), scaled to two dimensions. The distribution of transverse velocities of the sample of young pulsars presented here is broadly compatible with previous published models, and a detailed treatment of selection effects is required before we can usefully discriminate between the models.

Zoom In Zoom Out Reset image size

As an aside, we note that the present work provides the first model-independent distance and velocity estimate for PSR J2225+6535, the Guitar Nebula pulsar (v_⊥ ∼ 800 km s⁻¹; Table 4). That is the highest well-measured velocity in the current sample, although somewhat lower than previous estimates (e.g., Chatterjee & Cordes 2004); further analysis and comparison to long-term optical monitoring of the time-evolution of the Hα bow shock nebula is underway.

As highlighted in Section 1, having independent measurements of pulsar distances and astrometric parameters is extremely valuable for several pulsar science cases. While pulsar timing can provide the pulsar DM, extracting a distance from this measurement is dependent on having an accurate model of the Galactic electron density (see Section 4.2). Moreover, multi-frequency (or wide-band) observations are required to obtain a good DM measurement.

Astrometric terms in pulsar timing models are covariant with red noise in timing data that arises from fluctuations in the pulsar spin-down and/or propagation delays through the ISM that are not completely captured in the pulsar ephemeris. As shown in Deller et al. (2016) and Madison et al. (2013), red noise can lead to substantial errors in the values and underestimates of the uncertainty for timing-derived parameters. These errors are especially large for unrecycled pulsars with surface field strengths ∼10¹¹–10¹³ G, but are present at some level for all pulsars.

The highest precision measurements from pulsar timing are obtained with MSPs due to their frequent, short pulses and stable rotation. In our sample we have four MSPs: PSR J1022+1001, PSR J2010−1323, PSR J2145−0750, and PSR J2317+1439. The VLBI results for PSR J1022+1001 and PSR J2145−0750 were already presented and discussed in Deller et al. (2016) and compared with the most recent pulsar timing measurements. For completeness, we summarize the results here for all MSPs in the PSRπ sample. Each of these pulsars is observed by at least one pulsar timing array (PTA) in the search for low-frequency gravitational waves (e.g., Verbiest et al. 2016).

PSR J1022+1001 is located extremely close to the ecliptic plane (ecliptic latitude β = −006), which in pulsar timing leads to suboptimal measurements of the position and proper motion, due to the components being covariant in equatorial coordinates. Moreover, low-ecliptic latitude pulsars have their LoS passing close to the Sun every year, which leads to annual increases of DM that may not be modeled optimally in pulsar timing (e.g., Tiburzi & Verbiest 2018).

Although PSR J2010−1323 is also located relatively close to the ecliptic plane (β = 649), pulsar timing has been able to measure the proper motion of this pulsar with relatively high accuracy, as shown in Table 7. This pulsar is observed by the European Pulsar Timing Array (EPTA) and the North American Nanohertz Observatory for Gravitational Waves (NANOGrav), but not the Parkes Pulsar Timing Array (PPTA). Figure 16 shows that the uncertainties from their timing programs are comparable to VLBI, and no significant discrepancies are seen. It is noteworthy that the 11 yr NANOGrav data set values (Arzoumanian et al. 2018) changed significantly from the 9 yr values (Matthews et al. 2016), which were inconsistent at the ≳2σ level with the VLBI and EPTA values (Desvignes et al. 2016).

Zoom In Zoom Out Reset image size

Table 7.  Comparison between Proper Motion Measurements

Pulsar	VLBI	EPTA	NANOGrav	PPTA
J1022−1001 μα	$-{14.92}_{-0.03}^{+0.05}$	−18.2(64)	⋯	−17.09(3)
J1022−1001 μδ	${5.61}_{-0.04}^{+0.03}$	−3(16)	⋯	⋯
J2010−1323 ${\mu }_{\alpha }$	${2.36}_{-0.21}^{+0.33}$	2.53(9)	2.59(5)	⋯
J2010−1323 μδ	$-{5.61}_{-0.30}^{+0.26}$	−5.7(4)	−6.0(2)	⋯
J2145−0750 μα	$-{9.49}_{-0.04}^{+0.05}$	−9.58(4)	−10.0(2)	−9.59(8)
J2145−0750 μδ	$-{9.11}_{-0.08}^{+0.09}$	−8.86(10)	−8.0(5)	−8.9(3)
J2317+1439 μα	$-{1.43}_{-0.08}^{+0.08}$	−1.19(7)	−1.36(2)	⋯
J2317+1439 μδ	${3.74}_{-0.18}^{+0.18}$	3.33(13)	3.49(4)	⋯

Note. All values are given in mas yr⁻¹. For PSR J2317+1439, the VLBI value presents the results of the least-squares fit. Timing references as mentioned in the caption of Figure 16. Uncertainties on the timing parameters refer to the last digit(s) quoted.

Download table as:  ASCIITypeset image

PSR J2145−0750 is observed by all three PTAs, and is a good example where VLBI measurements give up to an order-of-magnitude improved accuracy compared to some timing measurements. The ecliptic latitude of this pulsar is 5.31 degrees, and this pulsar is known to show DM variations that affect the timing observations (e.g., Arzoumanian et al. 2018). As was the case for PSR J2010–1323, the proper motion obtained for PSR J2145−0750 from the NANOGrav 11 yr data set is considerably less discrepant with other measurements (our VLBI results, and also the PPTA and EPTA timing results) than the proper motion from the NANOGrav 9 yr data set was.

As described in Section 4.1, the VLBI astrometry of PSR J2317+1439 resulted in relatively poor constraints due to failed observations, which all fell on the same side of the parallax extrema. This also resulted in very conservative and skewed uncertainties for the proper motion parameters when using the bootstrapping method (green error ellipse in Figure 16), and all timing-derived measurement are consistent with the VLBI values. When using the least-squares fitting method, the uncertainties are reduced (black curve in Figure 16), and the EPTA timing value is offset by about 1σ from the VLBI value. Additional observations for PSR J2317+1439 could greatly reduce the VLBI proper motion uncertainty and provide a much more stringest comparison against timing.

Overall, Figure 16 shows that although for MSPs in some cases the timing measurements of proper motion parameters are comparable to the accuracy of VLBI measurements, the actual values can differ significantly between PTAs. As discussed in Deller et al. (2016), there could be multiple explanations, such as contamination by annual DM variations, systematic instrumental noise, the use of different versions of solar system ephemeris (SSE), or including different levels of noise modeling in the timing solutions. Arzoumanian et al. (2018) find the effect of using a different SSE on the proper motion to be on the order of 10 μas yr⁻¹ or less, which is insignificant compared to the current uncertainty levels.

Besides further comparison between timing models and, for example, the effect of including different types of noise modeling and DM correction, an extended set of independent and improved VLBI measurements of MSP proper motions will be extremely useful to find the underlying causes of any discrepancies between measured values.

Table 8 shows that in order to derive an independent distance measurement, VLBI observations can play an essential role to improve pulsar timing. In general, measuring a parallax signature in pulsar timing data is dependent on having a long baseline of observations, and as seen in Table 8 the uncertainties are between a factor of 2 and an order of magnitude larger compared to what is achieved with VLBI within 2 yr.

Table 8.  Comparison between Parallax Measurements

Pulsar	VLBI	EPTA	NANOGrav	PPTA
J1022−1001	${1.39}_{-0.03}^{+0.04}$	0.72(20)	⋯	1.1(3)
J2010−1323	${0.48}_{-0.12}^{+0.17}$	⋯	0.3(1)	⋯
J2145−0750	${1.60}_{-0.01}^{+0.06}$	1.53(11)	1.6(4)	1.84(17)
J2317+1439	${0.65}_{-0.07}^{+0.07}$	0.7(3)	0.50(8)	⋯

Note. All values are given in mas. For PSR J2317+1439, the VLBI value presents the results of the least-squares fit. References for the timing values are identical to those given in the caption of Figure 16. Uncertainties on the timing parameters refer to the last digit(s) quoted.

Download table as:  ASCIITypeset image

Finally, a comparison was made between the proper motions of the non-MSPs in our sample and the proper motions derived from timing observations. However, none of these pulsars had sufficiently significant detections from timing (Hobbs et al. 2004; Zou et al. 2005; Li et al. 2016) to make a useful comparison with our VLBI-derived values. When compared to previously derived interferometric measurements (see Table 9), we find that al most all are consistent within the given (generally low precision) error bounds; the exceptions are PSR B0329+54 and PSR B1133+16, for which the previous VLBI proper motion measurements of Brisken et al. (2002) are discrepant at the ≳2σ level in the decl. coordinate.

Three factors contribute to the accuracy of pulsar absolute positions that we can obtain from phase-referenced VLBI observations:

• 1.  
  The accuracy of the off-beam calibrator absolute position
• 2.  
  A contribution from the frequency-dependent core-shift of the off-beam calibrator
• 3.  
  The accuracy of the determination of position offsets with respect to the off-beam calibrators

We now consider each of these effects in turn for the PSRπ sample.

Among 60 sources used as off-beam calibrators, the absolute position accuracy (as recorded in the RFC; http://astrogeo.org/rfc/) ranged from 0.10 to 0.37 mas, with a median of 0.17 mas.

We have not measured the core-shift of any of the off-beam calibrators in our observations. Therefore, we can present only a rough estimate of its unaccounted contribution. Sokolovsky et al. (2011) presented results of multifrequency observations of core-shift. The core-shift of 17 AGNs at 1.6 GHz varied from 0.4 to 2.2 mas with the median 1.1 mas. We can take this estimate and assume it is typical for the sources used by the PSRπ sample.

According to Table 3, uncertainties in the position offsets from the pulsars to their position reference source range from 0.04 to 1.1 mas, with a median of 0.09 mas. The offset from the position reference source to the out-of-beam calibrator is typically be an order of magnitude greater, given the typical angular separations (14' median separation from position reference to target, versus 19 for off-beam calibrator to target).

Assuming these sources of errors to be independent, the overall uncertainties in the PSRπ pulsar absolute positions range from 0.4 to ∼10 mas, with a median of 1.4 mas (Table 4). The unaccounted core-shift in the off-beam calibrator and the position offset from the off-beam calibrator to the position reference source contribute roughly equally in most cases, and generally dominate over the absolute position uncertainty of the off-beam calibrator.

The absolute pulsar position does not play a significant role in the context of our study, but this does not mean it is not valuable at all. Comparison of the VLBI pulsar positions with positions determined with pulsar timing provides important information. First, a determination of the net rotation of pulsar positions determined with timing against the positions determined with VLBI can be used to improve three parameters that describe the orientation of the Earth orbit in the inertial space (i.e., the position of the ecliptic pole and the point of vernal equinox). Second, analysis of the residual position differences between timing and VLBI after the removal of net rotation gives us a measure of possible systematic errors in VLBI and/or timing. Timing and VLBI position estimates are to a great extent independent, and therefore their intercomparison provides us a unique opportunity to make an assessment of their accuracy.

Dedicated VLBA observations are able to bring the absolute accuracy of calibrator positions down to at least 0.1 mas. As Fey et al. (2015) demonstrated via decimation tests, a level of 0.05 mas for random position errors can even be reached if sources are observed long enough. Reaching that level of accuracy requires additional observations in the mode of absolute astrometry, similar to the regular geodesy "RDV" program conducted on the VLBA (Petrov et al. 2009).

The most important step for improvement of accuracy of pulsar VLBI absolute astrometry is determination of the core-shift. This requires multi-frequency dedicated observations and analysis, which we have not utilized for PSRπ. Sokolovsky et al. (2011) describes the technique of such observations. As it was shown in this and following works, frequency dependence of observed core-shifts obeys the power law ν^−1/r, with r close to 1. It was shown theoretically by Lobanov (1998) that in a case if (i) the plasma is in the state of equipartition with the magnetic field, (ii) the dominating absorption mechanism is synchrotron self-absorption, and (iii) the jet has a conical shape, then r = 1.

Determination of the core-shift requires significant observational resources, and it is not practical to do it for every pulsar. The PSRπ program is presently being continued with a second sample consisting exclusively of MSPs: MSPSRπ, with observations of 18 MSPs taking place during the period 2015–2018 and first results reported in Vigeland et al. (2018). To facilitate a comparison between high-precision pulsar timing and VLBI absolute positions, we have decided to focus on a list of 17 MSPs drawn from PSRπ and MSPSRπ that have timing positions accurate to 1 mas or better. In 2018 March we commenced a VLBA program targeting their off-beam calibrators, aiming to improve the absolute positions of these sources to the 0.1 mas level, as well as measuring the core-shift at 1.6 GHz. As of 2018 August, 28% of planned observations have been conducted.

For pulsars in binary systems, an optically visible companion offers the possibility of comparing VLBI positions against Gaia positions for truly point-like sources, avoiding the problems of systematic VLBI–Gaia offsets that are located preferentially along the jet direction of AGN (Kovalev et al. 2017). Unfortunately, pulsars with bright optical companions are rare; of the four binary pulsars in PSRπ, none have a companion above the Gaia magnitude limit. Given that more than 9000 common sources can be identified in VLBI and Gaia catalogs (Petrov et al. 2019), the sheer weight of numbers means that an ensemble comparison of radio/optical AGN will provide a more accurate alignment of VLBI and Gaia positions than will be possible with pulsar companions.

As discussed in Section 3.2, the astrometric accuracy achieved for a given target is expected to be influenced by the target brightness, the calibrator brightness, the target–calibrator separation, and the average observing conditions (principally the magnitude of the ionospheric gradients and the observing elevation). In Figures 17–19, we explore this hypothesis by examining the parallax uncertainty achieved for PSRπ pulsars as a function of each potential influence separately, before attempting to find an parameterized function that predicts the achieved accuracy given the known observing parameters. In all plots, we have excised PSR J2317+1439, where the parallax uncertainty is artificially inflated by the non-detections of the pulsar discussed in Section 4.1.5, and the red line shows a best-fit linear regression. The parallax uncertainty shown is the width of the 68% confidence interval in milliarcseconds (i.e., the addition of the uncertainties in the positive and negative directions).

Zoom In Zoom Out Reset image size

Figure 17 shows the parallax uncertainty plotted against the average S/N achieved on the target. Because the brightest pulsars achieve very high S/N, a log scale is used for the x-axis. As expected, the faintest targets tend to have higher uncertainties, but systematic errors dominate in most cases.

Figure 18 shows the parallax uncertainty plotted against separation to the position reference source. A weak trend toward larger uncertainties at larger separations is seen, but with a large scatter (r² = 0.08, where r is the correlation coefficient).

Zoom In Zoom Out Reset image size

Figure 19 shows the parallax uncertainty plotted against the average S/N achieved on the inbeam calibrator, or the quadrature addition of the S/N if multiple calibrators were used. Because the brightest calibrators achieve very high S/N, a log scale is used for the x-axis. The pulsars with the faintest calibrators tend to have higher uncertainties, but calibrator brightness is no more dominant than calibrator–target separation.

Zoom In Zoom Out Reset image size

Finally, we use a multiparameter estimation including calibrator S/N, elevation–weighted calibrator–target separation, and target S/N to attempt to predict the parallax accuracy achieved in a PSRπ observing setup. The predicted parallax error ${\pi }_{{\rm{p}}}$ in mas is given by

$$\begin{eqnarray}&&{\pi }_{{\rm{p}}}=\sqrt{{\left(\displaystyle \frac{A}{{\rm{S}}/{{\rm{N}}}_{t}\times \sqrt{{N}_{\mathrm{obs}}}}\right)}^{2}+\ {\left(\displaystyle \frac{B}{{\rm{S}}/{{\rm{N}}}_{c}}\right)}^{2}+{\left(C\times {\rm{\Delta }}\theta \right)}^{2}},\end{eqnarray} \tag{ 2 }$$

where A, B, and C are constants which we fit from our data set and find values A = 9.0, B = 4.5, C = 0.0028; ${\rm{S}}/{{\rm{N}}}_{t}$ is the average S/N on the target; N_obs is the number of observations; S/N_c is the average S/N on the in-beam calibrator; and Δθ is the angular separation of the pulsar and in-beam calibrator in arcminutes divided by the sine of the average observing elevation.

Figure 20 shows the actual parallax uncertainty π_obs versus the predicted value given by Equation (2), with the red line highlighting the expected 1:1 relationship. More than 80% of the pulsars fall within the range 0.5 × π_p < π_obs < 2×π_p, indicated by the gray lines on the plot.

Zoom In Zoom Out Reset image size

The values of the constants fitted in Equation (2) can be used with caution to estimate the probable astrometric accuracy of a future 1600 MHz VLBA astrometric campaign of comparable duration (8 epochs). For instance, in the limit of a very bright target and calibrator, then for an observation with typical observing elevation of 45°, the expected parallax uncertainty is 4 μas per arcminute of separation to the in-beam calibrator. Alternatively, given a target accuracy, these results can be used to estimate the characteristics of the in-beam calibrator that will be required. For a bright target, if the desired parallax accuracy is 20 μas, the in-beam calibrator should be separated by no more than 5 arcmin from the target, and should be at least ∼20 mJy (in order to achieve the necessary S/N of 225). Obviously, extending the number of observations in the campaign could be undertaken to lower these expected limits.