C M Espinoza

Lucas Guillemot

O Celik

P Weltevrede

B W Stappers

D A Smith

M Kerr

V E Zavlin

Ismaël Cognard

R P Eatough

Ö C ¸elik

P C C Freire

G H Janssen

F Camilo

G Desvignes

J W Hewitt

X Hou

S Johnston

M Keith

M Kramer

A Lyne

R N Manchester

S M Ransom

P S Ray

R Shannon

G Theureau

N Webb

Six millisecond pulsars detected by the

come

could be attributed to caustic radio emission produced at a range of different altitudes in the magnetosphere. We note that most of these characteristics are also observed in the Crab pulsar, the only other radio pulsar known to exhibit phase-aligned radio and gamma-ray emission.

Key words: pulsars: general -gamma-rays: general -X-rays: general.

I N T RO D U C T I O N

Pulsed gamma-ray emission from more than 100 pulsars has been detected by the Large Area Telescope (LAT) aboard the Fermi Gamma-ray Space Telescope1 (Abdo et al. 2010b;Abdo et al., in preparation). About a third of these are millisecond pulsars (MSPs), fast rotators that spin-down very slowly and steadily. MSPs are thought to have been spun-up through accretion of material from an evolved companion star (Alpar et al. 1982;[START_REF] Radhakrishnan | [END_REF]. During accretion, the system appears as a low-mass X-ray binary (LMXB) and no pulsed radio emission is observed. Although other possibilities have been considered (e.g. Ruderman 1991), it is commonly believed that accretion is responsible for the low surface magnetic fields observed amongst MSPs. In this model, the original field has been buried by the infalling material, decreasing its surface strength and decreasing the magnetic torque acting on the neutron star [START_REF] Romani | [END_REF]Bhattacharya & van den Heuvel 1991).

Most MSPs are found in binary systems. Among those in binaries, there are a growing number in very tight binary systems with lowmass companions, the so-called black widows (e.g. [START_REF] Roberts | Radio Pulsars: An Astrophysical Key to Unlock the Secrets of the Universe[END_REF]. It is believed that most solitary MSPs came through the black widow evolutionary path, in which the companion is ablated by the pulsar wind, either during the X-ray phase or once the pulsar turns on in the radio (Eichler & Levinson 1988;Ruderman et al. 1989).

The gamma-ray emission mechanism for MSPs and normal pulsars is believed to be the same (Abdo et al. 2009). MSPs exhibit the lowest gamma-ray luminosities among gamma-ray pulsars, a consequence of their low spin-down energy rates ( Ė), compared to young pulsars (Abdo et al. 2010b). Most observed gamma-ray pulse profiles consist of two dominant, sharp peaks, which suggest the emission is caustic in nature (Abdo et al. 2010b). Among the different models for radio and gamma-ray emission of pulsars, these properties tend to favour models in which gamma-rays are generated in the outer magnetosphere (e.g. Romani & Yadigaroglu 1995;Dyks & Rudak 2003) and radio emission is produced at lower altitudes (e.g. Rankin 1993). This is supported by the fact that most gamma-ray detected radio pulsars exhibit gamma-ray pulses out of phase with the radio pulsations (Abdo et al. 2010b). However, this scheme has been challenged by the increasing number of MSPs found to exhibit phase-aligned radio and gamma-ray pulse profiles, suggesting that, at least in these cases, they are both produced at a similar location in the magnetosphere (Abdo et al. 2010a;Freire et al. 2011a;Guillemot et al. 2012b). It has been noted that some of the MSPs showing phase-aligned radio and gamma-ray emission exhibit very low levels of radio linear polarization (Guillemot et al. 2012b), which is predicted to some extent by caustic emission models (Venter, Johnson & Harding 2012).

Before these discoveries, the Crab pulsar was the only pulsar known to exhibit phase-aligned radio and gamma-ray emission (Kuiper et al. 2003). MSPs present radio polarization properties similar to normal radio pulsars (Xilouris et al. 1998;Yan et al. 2011;Keith et al. 2012) but often much wider pulse profiles. Interpulses and bridge emission between different peaks are more common among MSPs (Yan et al. 2011). MSPs present radio spectra with a slope similar to normal radio pulsars, although arguably slightly steeper (Kramer et al. 1998;Toscano et al. 1998;Maron et al. 2000).

As well as detecting many new radio MSPs in gamma-rays (Abdo et al. 2009), Fermi has also pointed the way to many MSP discoveries via their gamma-ray emission properties (e.g. Cognard et al. 2011;Keith et al. 2011;Ransom et al. 2011;Guillemot et al. 2012a;Kerr et al. 2012). Here, we report on another six MSPs detected in gamma-rays by the Fermi LAT, all of which were previously known radio pulsars (PSRs J0610-2100, J1024-0719, J1600-3053, J1713+0747, J1741+1351 and J2051-0827). This is the first time that pulsed gamma-rays have been detected from these sources, except for PSR J2051-0827, for which a 4σ detection was reported by Wu et al. (2012).

By combining these six MSPs with 24 previously reported gamma-ray MSPs, we are able to start evaluating trends in the emission properties of the gamma-ray MSP sample. In particular, we use this sample to study the pulse profile properties and other connections between the radio and gamma-ray characteristics of these objects.

In Section 2, we describe the methods used to detect and analyse the gamma-ray pulsations from the six new gamma-ray MSPs and in Section 3, we comment on their multiwavelength properties. In Section 4, we present a search for orbital modulation in the gammaray emission of the five MSPs in binary systems. We present a study of the radio and gamma-ray properties of most known gamma-ray MSPs in Section 5 and we discuss this study in Section 6. We provide a summary of our work in Section 7.

M E T H O D O L O G Y: D E T E C T I N G G A M M A -R AY P U L S AT I O N S

The detection of pulsed gamma-ray emission from known MSPs is possible through the use of precise rotational ephemerides obtained through frequent radio observations, which are used to assign a rotational phase to each gamma-ray photon (Smith et al. 2008). Photons are then binned in the pulse phase to create a histogram that represents the light curve of the pulsar.

Radio analysis

We use pulse times of arrival (TOAs) obtained from radio observations with the Arecibo Observatory (AO; Freire et al. 2011b), the Nanc ¸ay Radio Telescope (NRT; Theureau et al. 2005), the Jodrell Bank Observatory (JBO; Hobbs et al. 2004), the Parkes Observatory (PKS; Weltevrede et al. 2010) and the Westerbork Synthesis Radio Telescope (WSRT; Voûte et al. 2002;[START_REF] Karuppusamy | [END_REF]. See Section 2.1.1 for more details on the use of the different data sets.

To study the phase-alignment between gamma and radio pulses, it is necessary to account for the delay suffered by the radio waves in their passage through the interstellar medium. In order to do this precisely, accurate estimates of the dispersion measure (DM)

Table 1. Properties of the radio ephemerides used to fold the gamma-ray data. The columns show the observatories involved in the observations, the total number of TOAs used (N TOAs ), the root mean square variation (RMS) of the timing residuals and the time range of validity for the timing solution. The last two columns show the DM values and the uncertainty of the phases assigned to the gamma-ray photons caused by the uncertainties on the DM, δφ. The DM uncertainties on the last quoted digits are in parentheses. These uncertainties correspond to the formal errors given by TEMPO2 and might under-represent the real uncertainties. The photon phase errors (δφ) caused by the uncertainty on the DM were calculated using these values.

Pulsar

Observatories are essential, especially for rapidly rotating pulsars (Smith et al. 2008). With the aim of having accurate and updated DM values, we measured the DM using multifrequency radio data taken during the time of the Fermi mission for five of the six gamma-ray MSPs (Table 1).

The best available radio TOAs were used to produce radio ephemerides valid from at least three months before the start of the Fermi mission until 2011 October. Ephemerides were produced by analysing the TOAs with the timing software TEMPO2 (Hobbs, Edwards & Manchester 2006). Because of the use of different pulse-profile templates to match the observations -hence different points of reference -and the existence of small clock divergences between observatories, the different data sets sometimes require phase-alignment. Thus, in addition to the spin frequency and its first derivative, binary parameters, position, proper motion and DM, we also fitted for phase delays between the data sets from different observatories. More information on the timing solutions built for each pulsar can be found in Table 1 and Section 2.1.1.

While pulsar timing can provide most of the known parameters describing pulsars, distances can be obtained by different techniques. If no better distance estimate was available (e.g. via parallax), the DM was used to calculate a distance based on the Galactic Free Electron Model (NE2001; Cordes & Lazio 2002). The best distance estimates (d) for the six MSPs are listed in Table 2. We note that the error bars of the DM-based distances could potentially be largely underestimated, affecting the error bars quoted for gamma-ray luminosities and efficiencies (Section 2.2.1).

Because MSPs present low rotational period derivatives ( Ṗ ), the Shklovskii effect (Shklovskii 1970) can be significant for those objects with a large proper motion (such as PSR J1024-0719; see below). We have corrected all the Ṗ values for this effect and the corrected values can be found in Table 2. Other radio properties, such as position, rotational and orbital periods and proper motion, can also be found in Table 2.

Building timing solutions, dispersion measures and distances

For PSR J0610-2100, most observations during the last four years have been carried out at ∼1.4 GHz, and so it is not possible to perform a sensitive analysis of the DM during the time of the Fermi mission. None the less, using a few ∼2-GHz TOAs and the small spread of observing frequencies around 1.4 GHz of the other avail-able TOAs, we are confident that the DM has not changed by more than 0.02 per cent with respect to the value quoted in the radio discovery paper (Burgay et al. 2006). This ensures an accuracy better than 0.006 rotations on the phases assigned to the gamma-ray photons.

For PSR J1024-0719, we measure a DM of 6.488 ± 0.001 pc cm -1 , which indicates a distance of 380 ± 40 pc, based on the NE2001 model. However, the Lutz-Kelker bias-corrected distance is 500 ± 100 pc (Hotan, Bailes & Ord 2006;Verbiest et al. 2012). As discussed below, considerations of the transverse movement and the Shklovskii effect are important before adopting any of these values. There are slightly different published proper motion values for this pulsar (see Bailes et al. 1997;Hotan et al. 2006). Using five years of NRT data at ∼1.4 GHz, we fitted for position and proper motion and obtained values entirely consistent with those of Hotan et al. (2006). The large total proper motion measured (∼59 mas yr -1 ) implies a significant contribution to the observed Ṗ value by the Shklovskii effect. In fact, it seems that most of the observed spin-down comes from this effect. The parallax distance indicates a negative intrinsic Ṗ and the maximum possible distance, in order to obtain a positive Ṗ , is 410 pc. Using the lowest possible distance given by the parallax uncertainties (i.e. 400 pc), we calculate the upper limit Ṗ ≤ 0.05 × 10 -20 . We use this Ṗ upper limit and adopt the Lutz-Kelker bias-corrected distance calculated by Verbiest et al. (2012).

For PSR J1600-3053, the timing solution was obtained using data at two main frequencies (1.4 and 2.1 GHz), which allowed us to calculate an accurate and up-to-date DM. We have also included the first time derivative for the DM in the fit (i.e. it was allowed to vary linearly with time), measuring a value consistent with that of You et al. (2007).

For PSR J1713+0747, a distance of 0.9 ± 0.1 kpc has been obtained from the NE2001 model using our DM estimate (Table 2). However, using 12 yr of radio timing, Splaver et al. (2005) measured this pulsar's yearly parallax and calculated a distance of 1.1 ± 0.1 kpc, which is the value we adopt here.

For PSR J1741+1351, the timing solution was obtained with observations carried out at the AO using the L-wide receiver, recording data between 1.1 and 1.6 GHz, which offers enough bandwidth to monitor DM variations. Based on the NE2001 model and the DM obtained from these observations, the distance to this pulsar would be 0.9 ± 0.1 kpc. The same data set contains the clear effects of a yearly parallax, suggesting a similar distance of ∼1.08 ± 0.05 kpc (this is a preliminary value; Freire et al., in preparation), which is the value we use.

For PSR J2051-0827, the orbital parameters are changing rapidly. Stappers et al. (1998) have reported the measurement of the first time derivative of the orbital period, Doroshenko et al. (2001) have detected up to three time derivatives of the orbital period and Lazaridis et al. (2011) have detected variations of the orbital period and the projected semimajor axis. Almost four years of WSRT data at 1.4 GHz have been used to fit for position, proper motion and binary parameters, using the ELL1 binary model (Lange et al. 2001) implemented in TEMPO2, and including one derivative of the orbital period. A second WSRT data set taken at 0.35 GHz and JBO TOAs at 1.5 GHz were then used together with the 1.4-GHz data to fit for the DM and its first derivative, keeping all orbital and astrometric parameters fixed.

These timing solutions will be made available through the Fermi Science Support Center.2

Gamma-ray analysis

To study the gamma-ray emission from these six MSPs, we selected Fermi LAT data taken between 2008 August 4 and 2011 August 4 using the Fermi Science Tools (STs). 3 We restricted the data set to events with energies between 0.1 and 100 GeV, reconstructed directions within 15 • of the pulsar locations, zenith angles smaller than 100 • , and belonging to the 'source' class of the P7_V6 instrument response functions (IRFs). We rejected the data collected when the LAT rocking angle exceeded 52 • , when the instrument was not operating in the science observation (or configuration) mode or when the data quality flag was not set as good.

Spectral properties

The gamma-ray spectral properties of the pulsars were determined using a binned maximum likelihood method, as implemented in the pyLikelihood PYTHON module of the Fermi STs. This method fits a model representing the point sources in the selected region of interest (ROI) and the diffuse emission to the data, and finds the bestfitting parameters to optimize the likelihood function describing the data. Our models included all sources from the Fermi LAT Two-Year Source Catalogue (2FGL; Nolan et al. 2012) found within 20 • of the pulsar positions. The parameters of sources within 8 • of the pulsars were left free in the fit, while the parameters of sources more than 8 • away were fixed at the values listed in the 2FGL. The diffuse Galactic emission was modelled using the gal_2yearp7v6_v0 map cube. The residual instrument background and the diffuse extragalactic emission were modelled using the iso_p7v6source template. 4 The normalization of the diffuse components was left free in the fits. The first step of the spectral analysis involved modelling the pulsar spectra with exponentially cut-off power-law (ECPL) shapes of the form:

dN dE = N 0 E 1 GeV - exp - E E c . (1)
Here, N 0 is a normalization factor, is the photon index and E c is the cut-off energy of the spectrum. The spectral parameters for PSRs J0610-2100 and J1713+0747 found with this model are listed in Table 3, along with the derived integrated photon fluxes F and energy fluxes G above 0.1 GeV. These pulsars spectra have also been fitted with a simple power-law model to test the validity of the ECPL model, by comparing the goodness of fit for the two models.

For both pulsars, the power-law model was rejected with more than 3σ significance. For the other four pulsars, we were unsuccessful in fitting the spectra with all three parameters of the ECPL spectral shape left free, and we obtained unsatisfactory results, indicating bad convergence issues. To estimate the cut-off energy and the integrated energy flux G for these pulsars, a second spectral fit was performed using an ECPL model with two free parameters only, fixing to a value of 1.3, which is the average seen for the 32 strongest MSPs in the second catalogue of LAT pulsars (Abdo et al., in preparation). The results of these fits are also listed in Table 3. In this table, the first and second quoted uncertainties are statistical and systematic, respectively. For PSRs J0610-2100 and J1713+0747, the latter uncertainties were estimated by running the fitting procedure using bracketing IRFs, where the effective area was perturbed by the estimated uncertainties ±10 per cent at 0.1 GeV, ±5 per cent near 0.5 GeV and ±10 per cent at 10 GeV, using linear extrapolations (in log space) in between (Ackermann et al. 2012a). For the

- - Photon index, 1.2 ± 0.4 +0.1 -0.1 1.3 b 1.3 b Cut-off energy, E c (GeV) 1.6 ± 0.8 +0.3 -0.2 2.2 ± 0.7 +4.1 -0.9 5.0 ± 1.8 +42.4 -2.4
Photon flux, F (≥0.1 GeV, 10 -9 ph cm -2 s -1 ) 7 . 8± 2.5 +0.2 -0.2

4.0 ± 1.1 +3.4 -1.4 3.5 ± 0.9 +3.0 -1.3 Energy flux, G (≥0.1 GeV, 10 -12 erg cm -2 s -1 ) 6.6± 1.1 +0.4 -0.3 3.8 ± 0.7 +0.8 -0.5 5.2 ± 1.0 +0.7 -0.4 Luminosity, L γ /f (≥0.1 GeV, 10 33 erg s -1 ) 1 0± 6 +6 -6 0.11 ± 0.05 +0.05 -0.05 4 ± 3 +3 -3 Efficiency, η/f (≥0.1 GeV) 11 ± 27 +27 -11 >0.8 ± 0.3 +0.3 -0.3 0.5 ± 0.4
- - Photon index, 1.6 ± 0.3 +0.1 -0.2 1.3 b 1.3 b Cut-off energy, E c (GeV) 2.7 ± 1.2 +0.3 -0.3 3.1 ± 1.6 +11.7 -1.3 2.0 ± 0.5 +4.4 -0.8 Photon flux, F (≥0.1 GeV, 10 -9 ph cm -2 s -1 ) 13.3 ± 3.7 +0.2 -0.2 2.8 ± 1.1 +2.5 -1.0 4.3 ± 1.0 +2.5 -1.3 Energy flux, G (≥0.1 GeV, 10 -12 erg cm -2 s -1 ) 10.2 ± 1.5 +0.3 -0.4 3.1 ± 0.8 +0.7 -0.4 3.8 ± 0.7 +0.4 -0.3 Luminosity, L γ /f (≥0.1 GeV, 10 33 erg s -1 ) 1 . 5± 0.3 +0.3 -0.3 0.43 ± 0.12 +0.11 -0.07 0.5 ± 0.2 +0.2 -0.2 Efficiency, η/f (≥0.1 GeV) 0.44 ± 0.10 +0.08 -0.08 0.020 ± 0.005 +0.005 -0.003 0.10 ± 0.03 +0.03 -0.03
a Please see the discussion about this pulsar's efficiency in Section 2.2.1. b These values were fixed in the spectral analysis, and in these cases the systematic errors were calculated in a different manner. Details on the measurement of these parameters are given in Section 2.2.

other four MSPs, the systematic errors caused by fixing the photon index to a nominal value of 1.3 probably dominate the errors because of uncertainties in the IRFs. To estimate the systematic uncertainties resulting from this choice, we fitted the data using photon indices of 0.7 and 2, these values representing the extrema observed for 32 strong gamma-ray MSPs (Abdo et al., in preparation). The best-fitting parameters for = 0.7 and = 2 provide the limits on the actual values listed in Table 3. We note that the fits obtained with = 2 were inconsistent with the data in all cases, which resulted in very large values of E c for PSRs J1600-3053 and J1741+1351 and, therefore, large systematic errors for these objects. Finally, the gamma-ray luminosities L γ = 4πf Gd 2 and the efficiencies η = L γ / Ė for the conversion of spin-down luminosity into gamma-ray radiation were calculated assuming a beaming factor of f = 1 (see Watters et al. 2009, for more details), which is common for outer magnetospheric emission models (e.g. Venter, Harding & Guillemot 2009). These quantities are given in Table 3. Two objects present anomalous efficiencies. We obtain η > 0.8 for PSR J1024-0719, determined by the maximum Ė allowed by the Ṗ value adopted in Section 2.1.1, which is very low. However, if this MSP were at 350 pc, Ė would be larger and thus the efficiency would be smaller, and closer to values commonly observed (see the discussion on this MSP's distance in Section 2.1.1). We note that this problem could also be alleviated if f was smaller than 1. By modelling the gamma-ray emission of a few MSPs, Venter et al. ( 2009) calculated f for different geometric configurations (magnetic inclination and line of sight). Although they have found f values close to 1 for known gamma-ray MSPs, under some specific geometric configurations their calculations have shown that f could also be very small.

The other case is the large efficiency obtained for PSR J0610-2100 (η ∼ 10; see Table 3). An efficiency above 1 is unphysical, and to have η < 1 requires d 2 f < 1 kpc 2 . As noted above, it is possible that the flux correction factor is not exactly f = 1, but smaller. We also consider the possibility that the large DM-based distance (Table 2) is a result of material along the line of sight, not modelled in NE2001. Infrared images acquired by the Wide-field Infrared Survey Explorer (WISE; 4.6, 12 and 22 μm) and by IRAS (25, 60 and 100 μm) show pronounced nebulosity around PSR J0610-2100. A distance of 1 kpc in this direction corresponds to 15 pc cm -3 in NE2001. Typical cloud sizes and overdensities can accommodate this DM discrepancy. For a distance of 1 kpc (taking f = 1), the gamma-ray luminosity would be as low as 0.8 × 10 33 erg s -1 , the Shklovskii-corrected spin-down power would be 6 × 10 33 erg s -1 and the efficiency would be η = 0.1. A smaller distance also improves the transverse velocity estimate. The measured proper motion (Table 2) and a distance of 4 kpc imply a transverse velocity V T = 345 km s -1 , a rather large value, given that the mean value for MSPs has been estimated to be ∼90 km s -1 , with a dispersion of 20 km s -1 (Lyne et al. 1998). With a distance of 1 kpc, the transverse velocity of PSR J0610-2100 would be ∼85 km s -1 , very similar to the mean value for MSPs. However, we note that if f < 1, the distance (and thus the transverse velocity) could be larger. For example, if we impose the maximum value of V T allowed from the dispersion (i.e. V T = 110 km s -1 ; Lyne et al. 1998), then d = 1.3 kpc and f < 0.6.

Search for pulsations

It has been shown that pulsation searches can be made more sensitive by weighting each event by the probability that it originates from the considered gamma-ray point source (Kerr 2011;Guillemot et al. 2012b). These event probabilities depend on the spectra of the sources of interest and of the other sources in the region. Thus, the best-fitting spectral models and the Fermi ST gtsrcprob have been used to compute the event probabilities. The event arrival times were finally phase-folded with the Fermi plug-in distributed in the TEMPO2 pulsar timing package (Hobbs et al. 2006;Ray et al. 2011). Table 3 lists the weighted H-test statistics (Kerr 2011) obtained by reducing the regions of interest to 5 • around the pulsars. These test statistics values all correspond to pulsation significances larger than 5σ . Therefore, we have detected pulsed gamma-ray emission from PSRs J0610-2100, J1024-0719, J1600-3053, J1713+0747 and J1741+1351 for the first time, and we have confirmed the marginal detection of PSR J2051-0827 presented by Wu et al. (2012).

The probability-weighted gamma-ray light curves for the six MSPs are shown in Figs 123456, along with radio and X-ray profiles when available. The gamma-ray background levels in these figures were obtained by summing the probabilities that the events are not a result of the pulsar, as described in Guillemot et al. (2012b). Statistical error bars have been calculated as i w 2 i , where w i is the event probability and i runs over events in a given phase bin (Pletsch et al. 2012). In Table 3, we list the radio-to-gamma-ray lags δ and gamma-ray peak separations for pulsars with multiple gamma-ray peaks. Here, the positions of the gamma-ray peaks were determined by fitting the integrated light curves above 0.1 GeV with Lorentzian functions, and the positions of the radio peaks were defined as the maxima of the radio light curves.

M U LT I WAV E L E N G T H P RO P E RT I E S O F T H E S I X M S P s

Here, we present some basic information on these MSPs and detailed descriptions of their pulse profiles, in gamma-rays and radio. If X-ray observations are available, then descriptions of their main X-ray properties are also given.

PSR J0610-2100

PSR J0610-2100 is a 3.8-ms pulsar in a 6.9-h orbit with a ∼0.02-M white dwarf (Burgay et al. 2006). Given the low mass of the companion and the short orbital period, PSR J0610-2100 is one of the black widow binary systems. However, in this case, there are no radio eclipses or DM variations caused by the ablation of the companion star (Burgay et al. 2006), which are usually associated with these systems (e.g. PSR B1957+20; Fruchter, Stinebring & Taylor 1988). Pulsed gamma-ray emission seems to come mostly from the 1-3 GeV band (Fig. 1).

In X-rays, the position of PSR J0610-2100 has been imaged with the Swift X-ray Telescope (XRT) in six observations performed in 2010 March and 2011 April and May. The effective exposures of those observations were very short, between 0.4 and 4.1 ks, and with the pulsar position offset from the XRT aim point by angles between 0.2 and 4.6 arcmin. PSR J0610-2100 was not detected with Swift. We used the level 2 data taken in the photon counting mode from these six data sets (10.0 ks of the total effective exposure). We applied the approach developed by Weisskopf et al. (2007) for statistical estimates on source detections to put a 3σ limit on the XRT source count rate of 1.3 cnt ks -1 (in a 30-arcsec radius aperture, centred at the pulsar radio position). Assuming, for simplicity, a power-law X-ray spectrum with a photon index X = 2 and absorbing hydrogen column density N H = 0.9 × 10 21 cm -2 (equal to the total Galactic H I column in the pulsar direction5 ), this count rate translates into a 3σ upper limit of 1.1 × 10 32 (d/4.0 kpc) 2 erg s -1 for the X-ray luminosity of the pulsar in the band 0.3-10 keV. This limit would be more stringent if the distance were overestimated (see the end of Section 2.2.1).

PSR J1024-0719

This is a 5.2-ms isolated pulsar located less than 500 pc from the Sun (Bailes et al. 1997). An X-ray counterpart was proposed by Becker & Trümper (1999) based on ROSAT observations. Later, using XMM-Newton data, Zavlin (2006) detected pulsed emission. A candidate optical counterpart was reported by Sutaria et al. (2003), but the association could be a result of positional coincidence and further observations are necessary.

The gamma-ray pulse profile of PSR J1024-0719 exhibits one broad gamma-ray peak, roughly half a rotation wide, with its centre preceding the radio peak by about 0.5 rotations. Pulsed emission seems to be present mostly in the 1-3 GeV band.

We include the X-ray pulse profile of PSR J1024-0719 in Fig. 2, obtained from data collected with the EPIC-pn instrument, operated in the timing mode, in an XMM-Newton observation conducted in 2003 December for a 66-ks effective exposure (Zavlin 2006). We have processed the data as in Zavlin (2006), but using the latest XMM-Newton data reduction software (SAS v. 11.0.0). The light curve was obtained from 815 events extracted from columns 38-39 in the one-dimensional EPIC-pn CCD image and in the 0.3-2 keV energy range, to maximize the signal-to-noise ratio. We used the photons plugin6 for TEMPO2 to assign a phase to each selected photon, together with a timing ephemeris obtained from JBO observations between 2002 November and 2004 October.

The H-test reports a 7.5σ detection of the pulsed emission, confirming the result of Zavlin (2006). The estimated intrinsic pulsed fraction, corrected for the background contribution, is 54 ± 20 per cent. The X-ray pulse peaks at phase 0.1-0.15, as determined by fitting a series of harmonics to the pulse profile, and its shape indicates a predominantly thermal origin of the pulsar X-ray emission (see Zavlin 2006, for more details).

At ∼1.4 GHz, this pulsar presents a complex average radio pulse profile, consisting of about nine components spanning 0.2 rotations, which emerge from a wider base almost half of a rotation wide. The profile at 2 GHz looks very similar (Fig. 2). The emission at 1.4 GHz is almost completely linearly polarized in the leading part of the pulse, involving the three main peaks, but it is not polarized on the trailing part (Yan et al. 2011). There is little variation in the position angle (PA) detected across the profile.

PSR J1600-3053

This is a 3.6-ms pulsar in a 14.4-d orbit with its binary companion (Jacoby et al. 2007). Its gamma-ray peak is ∼0.3 rotations wide and exhibits a sharp leading edge and a slowly decaying trailing edge. Most pulsed emission comes from the two higher energy bands, although some pulsed emission at earlier phases might also come from the 0.1-1 GeV band. PSR J1600-3053 was observed with XMM-Newton in 2008 February with the EPIC-MOS and EPIC-pn instruments, oper- ated in the full window and timing modes, respectively, for ∼30ks effective exposures. Data were reduced with SAS v. 11.0.0. An examination of the EPIC-MOS images around the radio position yielded no detection of the pulsar. The emission detected in the one-dimensional EPIC-pn image was heavily contaminated by enhanced background and by another bright source in the field of view. Therefore, despite the much higher sensitivity of the EPICpn instrument, no timing analysis of these data turned out to be meaningful (see the example of PSR J0034-0534 in Zavlin 2006). Using the approach described in Section 3.1, we put a 3σ upper limit of 1.2 cnt ks -1 on the pulsar EPIC-MOS count rates, in the 0.3-10 keV band (as measured from a 30-arcsec radius aperture centred at the pulsar radio position). For the power-law X-ray model of X = 2 and absorbing hydrogen column density N H = 1.0 × 10 21 cm -2 , these count rates translate into a rather deep 3σ upper limit of 7.1 × 10 30 (d/2 kpc) 2 erg s -1 for the X-ray luminosity, in the 0.3-10 keV band.

The average pulse profile at 1.4 GHz is dominated by a sharp pulse, which is preceded by a broader additional component, less than 0.1 rotations earlier, with roughly half its amplitude. There are no major differences with the pulse profile at 2 GHz (Fig. 3). The emission at 1.4 GHz is 30 per cent linearly polarized and exhibits two orthogonal PA jumps, the second coincident with a sense reversal of the circular polarization (Yan et al. 2011).

PSR J1713+0747

This is a 4.6-ms pulsar in a 68-d orbit with a white dwarf (Foster, Wolszczan & Camilo 1993). The gamma-ray pulse profile of PSR J1713+0747 is wide and roughly triangular but not exactly symmetric, with a shallower trailing edge (Fig. 4).

In the X-ray band, the position of PSR J1713+0747 has only been observed in four short exposures by the Swift XRT in the photon counting mode (3.1 ks of the total effective exposure). Using these XRT level 2 data and using the same approach as in Section 3.1, we obtained a 3σ limit of 1.9 × 10 31 (d/1.1 kpc) 2 erg s -1 on the X-ray luminosity of this MSP in the range 0.3-10 keV (adopting N H = 0.5 × 10 21 cm -2 ).

The average pulse profile of PSR J1713+0747 at 1.5 GHz consists of one sharp pulse with a weak trailing component and two small components on its base, preceding the main pulse by 0.04 and 0.08 rotations (Yan et al. 2011). Beside some broadening at lower frequencies, there are no major differences between the pulse profile at 0.8 and 2 GHz. There is a third, very shallow component visible in the three radio bands trailing the main peak by a little more than 0.1 rotations (see Fig. 4). The emission at 1.4 GHz is almost 100 per cent linearly polarized at the leading and trailing edges of the profile. There are two orthogonal PA jumps at each side of the main peak and a third just before the shallow trailing component. The second jump (and possibly also the first) is coincident with a sense reversal of the circular polarization (Yan et al. 2011).

PSR J1741+1351

This is a 3.8-ms pulsar in a 16-d period binary system (Jacoby et al. 2007). The gamma-ray peak of PSR J1741+1351 spans ∼0.3 rotations and exhibits relatively sharp edges. Most of the emission is in the two upper energy bands and leads the radio peak, an unusual situation among gamma-ray pulsars. Other cases like this are the MSPs PSR J1744-1134 and PSR J2124-3358 (e.g. Abdo et al., in preparation).

PSR J1741+1351 was in the field of view of the Swift XRT in 2012 January for 3.7 ks. However, because it was projected very close to the edge of the XRT image (with a 7.3-arcmin offset angle), no meaningful limit for the pulsar X-ray flux could be derived.

The pulse profile of PSR J1741+1351 at 1.4 GHz exhibits one sharp main pulse and a smaller additional component, 0.1 rotation wide, leading the main pulse by ∼0.35 rotations (see Fig. 5). There is also a small component on the base of the main pulse, which appears more prominent in the 0.3-GHz profile. This component precedes the main pulse by a ∼0.05 rotation.

PSR J2051-0827

This 4.5-ms pulsar is one of the so-called black widow systems. It is in a 2.4-h orbit and exhibits radio eclipses for about 10 per cent of the orbital period, observable at low frequencies (≤0.6 GHz; Stappers et al. 1996).

The gamma-ray profile of PSR J2051-0827 exhibits one steep, ∼0.4 rotation wide, peak preceding the radio peak by half a period. The pulsed emission is similar in the three energy bands (Fig. 6). This profile is compatible with the light curve presented by Wu et al. (2012).

XMM-Newton observed PSR J2051-0827 in 2009 April with the EPIC-MOS and EPIC-pn detectors, all operated in the full window mode, for 44-and 37-ks effective exposures, respectively. Unfortunately, these observations suffered from numerous strong particle flares, which increased the background level by a factor of up to 10-12 (compared to the normal one) during about 80 per cent of the observations. This made the XMM-Newton data practically useless for a purposeful analysis.

Five observations of this object were conducted with the spectroscopic array of the Chandra Advanced CCD Imaging Spectrometer (ACIS-S) in the very faint mode in 2009 March and July for 9-ks effective exposures (see also Wu et al. 2012). The aim of these observations was to search for orbital variability of the pulsar X-ray flux and they did not provide a time resolution suitable to perform a timing analysis at the pulsar spin period. We reduced the data using the CIAO software, version 4.3. Despite the small number of source counts collected in each observation (8-12 counts in a 1-arcsec radius aperture centred at the pulsar radio position), the object was clearly detected. The combined data-set totalled 44 source counts in the 0.3-10 keV range, with a negligible (less than 1 per cent) background contamination. The spatial distribution of these counts was found to be consistent with the ACIS-S point-like source image. Of the extracted source counts, 90 per cent were detected at photon energies below 2 keV, indicating that the spectrum of PSR J2051-0827 is soft and likely of a thermal origin rather than of a non-thermal (magnetospheric) origin. Indeed, fitting a power-law model to the extracted spectrum resulted in a large best-fitting photon index X 4 and a hydrogen column density N H 2.3 × 10 21 cm -2 , which is significantly greater than the value of 0.6 × 10 21 cm -2 expected from the pulsar dispersion measure. Fixing N H = 0.6 × 10 21 cm -2 results in a power-law model with X 2.6 and unabsorbed flux F X 7.9 × 10 -15 erg cm -2 s -1 in the 0.3-10 keV range, in agreement with the result reported by Wu et al. (2012). The thermal blackbody model provided more reasonable parameters: N H < 1 × 10 21 cm -2 , a (redshifted) temperature T ∞ pc 2.7 × 10 6 K and a radius of the emitting area (hot polar caps on the pulsar's surface) R ∞ pc 0.04 (d/1 kpc) km. The corresponding bolometric X-ray luminosity is L ∞ pl 0.6 × 10 30 (d/1 kpc) 2 erg s -1 . Using a non-magnetic hydrogen atmosphere model for the thermal emission of MSPs (Zavlin, Pavlov & Shibanov 1996;Zavlin 2006[START_REF] Zavlin | Astrophysics and Space Science Library[END_REF], we obtained T pc 1.6 × 10 6 K, R pc 0.18 (d/1 kpc) km and L pc 0.8 × 10 30 (d/1 kpc) 2 erg s -1 . These are all unredshifted values (i.e. as measured at the neutron star surface).

Using the photons plugin for TEMPO2, we have searched for orbital variability of the X-ray photons detected in the Chandra observations. With the scant count statistics available, we have found that the significance of a possible modulation does not exceed 2.2σ . Xray data of much better quality are required to elucidate the origin of the emission (whether it is thermal polar cap radiation or nonthermal flux from the pulsar magnetosphere and/or interaction of the pulsar wind with the companion).

The pulse profile of this MSP at ∼1.4 GHz exhibits one main pulse and one overlapping trailing component, about 0.05 rotations apart. Additionally, there are at least three smaller trailing components extending up to 0.2 rotations away, forming a shoulder. The pulse profile at 0.4 GHz keeps the same principal structure but the two main peaks appear broader and slightly more separated. There are no major differences with the profile at 2 GHz (see Fig. 6). The emission at 1.4 GHz is mildly linearly polarized (∼12 per cent) throughout the whole profile (Xilouris et al. 1998).

S E A R C H F O R O R B I TA L M O D U L AT I O N O F T H E G A M M A -R AY E M I S S I O N

PSRs J0610-2100 and J2051-0827 are in close binary orbits, where interactions between the wind of the pulsar and the atmosphere of the companion star are likely to occur. Accordingly, we have searched for orbital gamma-ray flux modulations. For completeness, we have also searched for such modulations in the other three binary systems, PSRs J1600-3053, J1713+0747 and J1741+1351.

The Fermi LAT sensitivity to a source varies on multiple timescales, most notably the spacecraft's orbital period (∼95 min) and period of precession (∼53 d). Beats between the binary orbital frequency and harmonics of these time-scales induce an apparent modulation of photon rate as a function of orbital phase. Therefore, when conducting searches for intrinsic modulation, it is important to carefully correct for the time-dependent sensitivity (e.g. Corbet & Kerr 2010;Guillemot et al. 2012b). Because of the relatively low count rates, we decided to adapt an unbinned pulsation statistic, the H-test (de Jager, Raubenheimer & Swanepoel 1989), for use with uneven exposure. To do this, we computed the exposure to the source with 30-s resolution over broad energy bands (two per decade). We used TEMPO2 and the timing solutions of Section 2.1.1 to compute the orbital phase for each 30-s interval. We have taken the resulting distribution of phases, F , to represent the null hypothesis of no intrinsic modulation. By definition, the quantity F (φ), where φ is the observed orbital phase, is a uniformly distributed random variable in the absence of intrinsic modulation, and thus it is suitable for use in the H-test.

Using this modified H-test, along with the photon weights (see Section 2.2), we searched for a signal in individual energy bands (100-300, 300-1000, 1000-3000, 3000-10 000 and 10 000-30 000 MeV) as well as cumulatively (>100, >300, <300 MeV, etc.). We detected no significant orbital modulation of the gamma-ray signal from any of the sources.

This null result is physically expected for binary systems with long orbits, where no direct interaction between the two bodies is feasible. In the case of PSR J2051-0827, this result is also unsurprising because no evidence for interaction of the pulsar wind with the companion star is seen in the X-ray spectrum (see Section 3.6). Interestingly, the weighted orbital light curve at energies above 3 GeV for PSR J0610-2100 presents two peaks, suggesting the presence of orbital modulation of the emission at these energies (Fig. 7). However, there are very few events above 3 GeV and the peaks observed could be the effect of low statistics; indeed, the H-test for this set of events indicates only a 2σ significance. Furthermore, the light curve of the background photons (obtained by assigning a weight 1w i to every event in the same 5 deg ROI; see Section 2.2.2.) for this source seems to follow the same shape in Fig. 7. Although fairly flat, the exposure folded curve for PSR J0610-2100 also presents local peaks at similar orbital phases (∼0.2 and ∼0.5).

Finally, we performed a Monte Carlo study to determine the sensitivity to a few particular types of modulation. For each pulsar, we have used the observed weights and an assumed morphology for the orbital modulation to generate 100 random realizations of the orbital phases. Then, we have varied the strength of the modulation until 95 of the 100 simulations exceed the statistical threshold, taken to be H = 8, or a ∼2σ detection.

In general, we find that only strong modulations are detectable. Sinusoidal modulation is detectable for PSRs J0610-2100 and J1713+0747 if ∼70 per cent of the total flux is modulated, while not even 100 per cent sinusoidal modulation is reliably detectable for the remaining sources. However, such high levels of modulation, which imply additional, non-magnetospheric sources of gammaray emission from the system, are incompatible with the estimated background levels of the rotational phase light curves (Figs 123456).

For the case of a notch in the otherwise steady emission (an eclipse), we find that the sensitivity depends strongly on the notch width. Any gamma-ray eclipse is expected to be narrow, and unfortunately the LAT is largely insensitive to such small features. Eclipses become detectable for PSRs J0610-2100 and J1713+0747 when they span ∼20 per cent of the orbit, while even greater values (30-50 per cent) are required for the remaining sources. Thus, while we can rule out extreme cases of modu-lation (fully modulated sinusoids, broad eclipses), our sensitivity study shows that the null results discussed above are not highly constraining.

R A D I O A N D G A M M A -R AY P RO P E RT I E S O F T H E G A M M A -R AY M S P P O P U L AT I O N

Our sample of 30 MSPs is constituted by all MSPs known to emit in gamma-rays at the time of writing (published or publicly announced). We have compiled published radio pulse profiles and measurements of the radio spectral index (s), the radio flux densities, the average degree of radio linear polarization through the main radio pulse, the phase lag of the closest gamma-ray peak to the main radio peak (δ), the gamma-ray photon index ( ) and the flux for each of the 30 MSPs in the sample. Table 4 gives some of these quantities and also lists references where gamma-ray and radio pulse profile plots can be found.

Radio and gamma-ray pulse profile properties

Considering properties such as the phase-alignment between radio and gamma-ray peaks (above 0.1 GeV), the radio duty cycle and the number of pulse components, we define three main groups and label them A, N and W (Table 4). Although this analysis is based mainly on the inspection of average pulse profiles at ∼1.4 GHz, we have also inspected pulse profiles at lower and higher radio frequencies. None of the MSPs studied exhibits the appearance or disappearance of pulse components (although, see comments on PSR B1957+20 below and on PSR J0218+4232 in Section 5.1.2). We note that pulse-profile evolution with frequency is observed in some other MSPs.

Type A(ligned)

These are the MSPs that have their main gamma-ray peak aligned with the main radio pulse. There are five objects in this group (Table 4).

All of the A-type MSPs exhibit two main radio peaks, with each of them composed of one or more components. In two cases (PSRs B1937+21 and B1957+20), the second peak appears as an interpulse (i.e. half a rotation away from the main pulse). However, in the case of PSR B1937+21, low-level emission detected between the peaks suggests that both peaks might be a result of emission produced not at opposite magnetic poles but somewhere in the outer magnetosphere (Yan et al. 2011). Only five MSPs in the whole sample present secondary peaks about half a rotation away from the main pulse.

For the five A-type MSPs, all radio peaks have a gamma-ray counterpart. The only exception might be PSR B1957+20, for which the second radio peak visible at 1.4 GHz, leading the main peak by ∼0.26 rotations, appears not to have an obvious gamma-ray counterpart. We note that this peak is not visible at lower radio frequencies (e.g. 0.3 GHz).

Four MSPs in this group exhibit very low or undetectable levels of linear polarization in their radio emission (Table 4). This seems to be a property exclusive to A-type MSPs; all other gamma-ray MSPs in the sample for which polarization data were available exhibit normal levels of linear polarization. The only A-type MSP showing some degree of linear polarization is PSR B1937+21.

Type N(on-aligned)

These are MSPs that have their main gamma-ray peak out of phase with the main radio peak and their radio emission dominated by a single peak. Two MSPs in this group present gammaray pulses almost phase-coincident with their main radio pulses (PSRs J1744-1134 and J2214+3000). However, they are not exactly aligned and their peaks lead the radio peak, which is unusual among gamma-ray pulsars -although it is also seen in the light curves of the MSPs PSR J1741+1351 and PSR J2124-3358. These objects might belong to a different class of MSP, in terms of morphology (cf. Venter et al. 2009;Johnson 2011).

There seems to be some gamma-ray emission, besides the main peak of emission, at the radio pulse phase in at least three objects (PSRs J0340+4130, J0610-2100 and J1024-0719). There are 21 N-type MSPs in total, two of them exhibiting an interpulse (PSRs J0030+0451 and J0101-6422).

Type W(ide)

This is a group of four MSPs that have their main gamma-ray peak out of phase with the main radio peak and whose radio emission consists of multiple peaks of comparable amplitude, covering most of the rotational period. PSR J0218+4232 presents a secondary radio pulse component coincident with the gamma-ray peak of emission. This radio component is stronger at lower radio frequencies (<0.6 GHz), becoming comparable in strength to the main radio peak (Kuiper et al. 2002). There also seems to be alignment between some secondary radio components and the main gamma-ray peak for PSRs J2124-3358 and J2302+4442. Nevertheless, they are not proper A-type pulsars because not all of their radio peaks are aligned with gamma-ray peaks. In addition, it must be noted that given the large number of radio components and their long duty cycles, these could be mere coincidences. Besides these practical reasons, there is a physical reason to keep the W and N groups separated. The W-type pulsars are possibly aligned rotators (i.e. pulsars with their magnetic axes almost coincident with their rotation axes), a situation that could affect their observed properties.

Magnetic field at the light cylinder

The dipole magnetic field strength of a pulsar is normally calculated by assuming that all rotational energy losses come from dipole radiation. Its value at the light cylinder can be estimated by (e.g. [START_REF] Lorimer | Handbook of Pulsar Astronomy[END_REF])

B LC = 9.2 P 1 s -5/2 Ṗ 10 -15 1/2 G. (2)
At the light cylinder, the corotational velocity is equal to the speed of light, c, and its radius is therefore given by R LC = cP /2π. This implies that gamma-ray MSPs have the smallest light cylinder radii, followed by rapidly spinning young pulsars, such as the Crab pulsar. Because of the weak dependence of B LC on Ṗ and the different Ṗ ranges exhibited by MSPs and young pulsars, both populations share similar B LC values, the largest among the whole pulsar population (>10 4 G).

The A-type MSPs tend to have higher B LC values than the rest of the population. Fig. 8 shows B LC as a function of P for all known pulsars having P ≤ 0.1 s. If proper motion measurements were available through the Australia Telescope National Facility (ATNF) Pulsar Catalogue, the B LC values were corrected for the Shklovskii effect (Table 4). Lines of constant spin-down energy rate Ė = 4π 2 I Ṗ /P 3 (where I is the moment of inertia of the star, assumed to be 10 45 g cm 2 ) are plotted with dashed lines. Most MSPs populate the lower branch crossing the centre of the plot and the 30 gamma-ray MSPs in our sample populate its higher end, towards the shortest periods, and larger Ė and B LC values. The MSPs from the three types described in Section 5.1 are plotted using different symbols. The A-type MSPs are near the top-left corner, with the shorter periods and the higher B LC and Ė values (Johnson 2011).

PSR B1823-30A is one of the A-type MSPs in the plot and, considering its period, it shows a relatively large B LC . This MSP is in a globular cluster (GC) and has a very large Ṗ , hence its large B LC . Freire et al. (2011a) have studied the gamma-ray emission of this MSP and have concluded that the observed Ṗ is intrinsic to the pulsar and is not a result of acceleration in the cluster. MSPs in GCs are marked in Fig. 8 A two-dimensional Kolmogorov-Smirnov (KS) test (Fasano & Franceschini 1987;[START_REF] Press | Numerical Recipes: The Art of Scientific Computing[END_REF]) between A-type MSPs and N and W types combined indicates that the probability that they belong to the same P-B LC distribution is less than 1 per cent. However, the significance of the two-dimensional KS test remains valid while N e = N 1 N 2 /(N 1 + N 2 ) ≥ 20 [START_REF] Press | Numerical Recipes: The Art of Scientific Computing[END_REF]) and, in this case, N e = 4.2. To check these results further, we also tried the standard one-dimensional KS test, which can be applied for a lower number of objects (valid for N e ≥ 4; [START_REF] Press | Numerical Recipes: The Art of Scientific Computing[END_REF]. Comparing the B LC values of the five A-type MSPs with the 25 N and W types gives 0.02 per cent probability that they belong to the same distribution.

To quantify how likely it is to obtain such a value by selecting a small group at random, we randomly picked five MSPs from the sample and calculated the KS test against the remaining 25 objects. After repeating this 30 000 times, we found that for more than 99 per cent of the cases, the KS test null hypothesis probability is greater than 0.02 per cent, with 67 per cent of the cases having a probability greater than 40 per cent. Therefore, we conclude that the distribution of the B LC values of A-type MSPs is significantly different to that of the rest of the gamma-ray MSPs.

None the less, we acknowledge that the above results might be driven by the very short periods that most of the A-type MSPs have, compared to those of the N and W types together. Indeed, while the KS test gives a probability of 41 per cent for the Ṗ values being from the same distribution, the same test gives only 0.4 per cent probability that the periods belong to the same distribution. However, this is about 10 times larger than the probability obtained for B LC . If Ė values are considered instead, the probability is 1 per cent, which is 50 times larger than what is obtained for B LC values. Moreover, as discussed in Section 6, a comparison with the Crab pulsar plus other considerations hint at B LC being the relevant parameter, rather than P or Ė.

Radio spectral indices

The radio spectral index s corresponds to the slope of the spectrum in a logarithmic scale, in which the flux density is described by S ∝ ν s . Out of the 30 MSPs in the sample, we have spectral information for 19 of them (Table 4). These data show that A-type MSPs tend to have lower spectral indices than the rest of our sample. Below, we assess the significance of this observation.

The distribution of spectral indices for all the gamma-ray MSPs, together with 33 additional MSPs for which spectral information was available, is shown in Fig. 9, indicating the contribution made by the A-type MSPs. Data were collected from the ATNF Pulsar Catalogue, from the compilations in [START_REF] Lorimer | [END_REF], Kramer et al. (1998) and Toscano et al. (1998), and from flux densities at different frequencies obtained via private communication with M. Keith and S. Ransom. Note that the error bars of the s values sometimes might not reflect the intrinsic difficulties associated with flux density measurements, which directly affect the spectral index estimate. While the s distribution of gamma-ray MSPs seems to follow the general distribution, with a mean spectral index above -2.0 (consistent with Kramer et al. 1998;Toscano et al. 1998), the mean value for the A-type gamma-ray MSPs clearly falls below -2.5 (see Fig. 9).

To test how significant the difference between the s distributions of the different sets is, we again use the KS test. A KS test over the s values of four A-type MSPs against the remaining 15 gamma-ray MSPs (Table 4) indicates a probability of 0.6 per cent that both sets belong to the same distribution. If the test is applied for the A-type MSPs against all other MSPs in Fig. 9, then the probability falls to 0.2 per cent. Therefore, it appears that the A-type MSPs have distinctively low s values. We note that for the first KS test mentioned above, N e = 3.2 (see Section 5.2). However, the test indicates a probability of 56.7 per cent that the s values of all gamma-ray MSPs belong to the overall s distribution. This indicates that, in general, there is no evidence for different radio spectral behaviour between gamma-ray MSPs and the rest of the MSP population.

We have also checked how likely it is to obtain an s distribution, such as the one exhibited by the A-type MSPs, by pure chance. To do so, we randomly picked four values from the s distribution of all gamma-ray MSPs and we calculated the probability of them being from the same distribution as the remaining 15. After repeating this process many times, we have found that more than 99 per cent of the trials give a probability greater than 0.6 per cent, with more than 66 per cent of the trials giving a probability greater than 40 per cent. Thus, the particular s distribution of A-type MSPs appears to be truly different.

We find no correlation between s and the gamma-ray spectral index or photon index, , for the MSPs in the sample. None the less, we note that values of are tightly clustered, with a scatter smaller than unity (Abdo et al. 2010b). Also, no obvious correlation was found between and any other investigated parameter.

Fluxes

There are radio flux densities at 1.4 GHz available for 23 of the 30 MSPs in our sample and gamma-ray fluxes for all of them. These last values are preliminary results of the ongoing effort by the Fermi collaboration to produce the second Fermi LAT Catalogue of gamma-ray pulsars (Abdo et al., in preparation). For many MSPs in the sample, flux values can also be found in the references given for the gamma-ray light curves in Table 4. We find no obvious correlation between radio flux densities and gamma-ray fluxes (Ackermann et al. 2012b). Fluxes are expected to be at least slightly correlated because of their mutual dependence on distance. The fact that this correlation is not observed in the sample suggests that radio and/or gamma-ray fluxes are highly dependent upon geometrical factors and probably other intrinsic properties dependent on P and Ṗ . No correlation between these fluxes, or their ratio, and any other quantity was found.

Linear polarization

We have found published average radio linear polarization L/I data for 17 MSPs in our sample. No correlation is found between this parameter and any other studied quantity. However, we note that the four MSPs with very low (or zero) L/I are all A-type MSPs (PSRs J0034-0534, B1820-30A, J1902-5105 and B1957+20; see Fig. 10).

D I S C U S S I O N

Gamma-ray MSPs at a glance

By studying the radio and gamma-ray pulse profiles of 30 MSPs, we have classified gamma-ray MSPs into three types, labelled A, N and W. The first group is composed of MSPs exhibiting pulsed gamma-ray emission phase-aligned with their radio emission. The N-type and W-type MSPs exhibit misaligned emission, with the N-types presenting radio emission dominated by a single pulse and the W-types presenting multiple and wide radio pulses of comparable amplitude. While the N-types and W-types appear to constitute the norm among gamma-ray MSPs, A-type MSPs appear less frequently and present various notable differences from the rest of the population. We note that these properties are not shared by PSRs J1744-1134 and J2214+3000, the two MSPs that could have possibly been regarded as A-type (see Sections 5.1.1 and 5.1.2).

The W-type MSPs are more likely to be aligned rotators, with their rotation axes almost coincident with their magnetic axes. However, this situation does not appear to drive obvious trends in either their gamma-ray profiles or any of the other studied parameters.

The difference in rotational phase between radio and gammaray pulses has commonly been interpreted as a difference in the location, within the pulsar's magnetosphere, where the emission is generated. Standard models normally place the source of radio emission at lower altitudes compared to the source of gamma-rays, in accordance with the radio/gamma-ray misalignment observed for most gamma-ray pulsars and MSPs. However, in light of the newly detected MSPs exhibiting phase-aligned emission, other possibilities ought to be considered. In the two most likely scenarios, both radio and gamma-ray emission are generated either relatively near the surface of the neutron star (low-altitude Slot Gap (laSG) models) or over an extended region at higher altitudes (altitude-limited Open Gap (alOG) and altitude-limited Two Pole Caustic (alTPC) models; see Venter et al. 2012). Although the alTPC model seems to offer better fits for some A-type MSPs, further studies are necessary to better understand the emission mechanisms operating in these pulsars.

None the less, models involving a caustic origin for the radio emission (at high altitudes) predict low levels of linear polarization and rapid swings of the polarization position angle (Dyks, Harding & Rudak 2004). In such models, rapid position angle swings with phase and depolarization of the total emission would be an effect of mixing signals from different altitudes (Venter, Johnson & Harding 2012). This is observed in some, but not all, of the A-type MSPs (Section 5.5; Table 4).

It has been noted that for high Ė non-MSP pulsars, the degree of linear polarization might be correlated with Ė (Weltevrede & Johnston 2008). Moreover, based on their wide beams, it has been argued that the radio emission of some of these pulsars is produced at a large range of heights. However, their polarization levels appear high, contrary to what is argued above to explain the lack of polarized emission in some A-type MSPs. Also, A-type MSPs have the largest Ė values among MSPs, defying the correlation found by Weltevrede & Johnston (2008). Nevertheless, normal pulsars have larger light cylinder radii than MSPs, and it might be that depolarization is less efficient in their larger magnetospheres.

In the smaller magnetospheres of MSPs, the production of radio emission could naturally extend up to larger fractions of R LC , producing the observed wide beams (Yan et al. 2011) and the phasealigned profiles of A-type MSPs [START_REF] Harding | Proceedings of the 22nd Texas Symposium on Relativistic Astrophysics[END_REF]. In fact, A-type MSPs have shorter rotation periods than other MSPs (Fig. 8). Based on a study of the beaming fractions of the radio and gamma-ray emission, Ravi, Manchester & Hobbs (2010) have concluded that the radio emission of high-Ė pulsars (including MSPs) must originate higher up in the magnetosphere, near to where the gamma-ray emission is produced. This is consistent with the high-Ė values of A-type MSPs. The Crab pulsar is the only non-MSP known to exhibit nearly phase-aligned radio and gamma-ray emission (Kuiper et al. 2003). We note that Ė for the Crab pulsar is almost three orders of magnitude larger than the highest values among MSPs and that many young gamma-ray pulsars have Ė values similar to those of the A-type MSPs, but show misaligned emission. Also, while the Crab pulsar has the second smallest R LC value among normal pulsars, it is more than one order of magnitude larger than those of the A-type MSPs. If R LC were the main factor determining the alignment between radio and gamma-ray emission, we would expect to see this alignment for most MSPs, which is not observed. If Ė were the main factor, we would expect to see more young gamma-ray pulsars exhibiting aligned radio/gamma-ray emission, which is not the case. Considering the KS tests described in Section 5.2 and the position of the Crab pulsar in Fig. 8, B LC appears to be a natural common property among A-type MSPs and the Crab pulsar. We note that the polarization levels of the Crab pulsar's radio emission appear to be well above zero (Gould & Lyne 1998).

The Crab pulsar and A-type MSPs all have similar radio spectra. The radio spectral index of the Crab pulsar is -3.1 [START_REF] Lorimer | [END_REF], substantially steeper than the averages of -1.9 for MSPs and of -1.8 for normal pulsars (Toscano et al. 1998;Maron et al. 2000). Although viewing angles and other geometric factors can bias our measurements, the steep radio spectra of the A-type MSPs and the Crab pulsar, together with their high B LC values, are likely to be related to magnetospheric similitudes and common processes on the generation of their emission.

The emission of giant radio pulses (GPs) is another common feature among A-type MSPs and the Crab pulsar. GPs are sporadic, short and intense bursts of radio emission, following power-law energy statistics. They were originally detected in the emission of the Crab pulsar, and later in PSR B1937+21 [START_REF] Heiles | [END_REF]Cognard et al. 1996). Currently, there are eight pulsars and five MSPs known to emit these type of pulses (Knight et al. 2005;Knight 2006). The GP properties exhibited by the Crab pulsar differ from those of most other pulsars but are very similar to those exhibited by the five MSPs found to emit GPs (Knight 2006). Three of these MSPs are A-type MSPs: PSRs B1937+21, B1957+20 (Knight et al. 2006) and B1820-30A (Knight et al. 2005) and the other two are PSRs J0218+4232 [START_REF] Joshi | Proc. IAU Symp. 218, Young Neutron Stars and Their Environments[END_REF]), a Wtype MSP and B1821-24 (Romani & Johnston 2001), for which no gamma-ray pulses have been detected with a confidence above 5σ (Pellizzoni et al. 2009). It has been proposed that the relatively large B LC values exhibited by these five MSPs and the Crab pulsar could be the main physical factor determining the generation and the main properties of their GPs (Cognard et al. 1996, but see the discussion by Knight et al. (2006)). Future observations and careful analyses should discern whether the emission of GPs is connected to the alignment between radio and gamma-ray emission or is mere coincidence. No GPs were detected for the A-type PSR J0034-0534 (Knight et al. 2005).

The wide pulse profiles that MSPs generally have, compared to normal pulsars (Yan et al. 2011), could be understood as evidence for outer magnetosphere caustic radio emission for all types of MSPs. We could argue that the outer magnetosphere offers enough room for different emission locations for radio and gamma-rays. However, we have shown that A-type MSPs exhibit different emission properties. Could the availability of larger magnetic fields at the light cylinder generate conditions, in relatively smaller magnetospheres, that favour the generation of radio emission at higher altitudes (at least as a fraction of R LC ) and co-located with the production of gamma-rays? Would this different emission mechanism naturally produce a steeper spectra and favour the production of giant pulses? Any model describing phase-aligned radio and gamma-ray emission should take these properties into consideration.

S U M M A RY

We have presented the detection, by the Fermi LAT, of gamma-ray pulsations from six MSPs, five of them detected for the first time and a sixth, PSR J2051-0827, confirmed at the 5σ level. The six MSPs present properties that are common among the gamma-ray MSP population. All of these pulsars are significantly detected in the radio domain, but because of the low conversion of spin-down energy into X-rays, only two are detected significantly in X-rays.

By studying the morphology and phase relationship of radio and gamma-ray pulse profiles of a sample of 30 MSPs, we have grouped gamma-ray MSPs into three types. The most distinctive type of gamma-ray MSPs consists of those exhibiting phase-aligned radio/gamma-ray emission. We find some clear trends in their emission properties, which differ significantly from the rest of the MSP population. We have shown that the MSPs in this group have a radio spectra steeper than the rest of the MSP population and that they also have among the highest inferred magnetic field strengths at the light cylinder. Additionally, some MSPs in this group have distinctively low degrees of radio linear polarization and some MSPs of this type are amongst the handful of MSPs known to emit giant radio pulses. Many of these properties are also observed from the Crab pulsar, the only normal radio pulsar known to emit gamma-ray emission phase-aligned with its radio emission.

The use of combined information obtained through the study of phase-aligned gamma-ray, X-ray and radio emission, along with their intrinsic properties, offers a wide perspective that certainly helps to improve our understanding of the emission mechanism of pulsars.

Figure 1 .

 1 Figure1. Gamma-ray and radio pulse profiles for PSR J0610-2100. The profiles are shown twice, here and elsewhere in the paper, to elucidate their complex nature and to clearly visualize alignment or otherwise. The gammaray light curve for events with energies above 0.1 GeV is shown in the top panel. The following panels show the light curve for different gamma-ray energy bands, as indicated. The histograms were made using photon weights, and the horizontal dashed lines represent the estimated background emission (Section 2). The bottom plot shows the radio pulse profile at the indicated frequencies.

Figure 2 .

 2 Figure 2. Gamma-ray and radio pulse profiles for PSR J1024-0719 (see Fig. 1 for details). The X-ray light curve is presented in the fifth panel from the top.PSR J0610-2100 has a 1.4-GHz pulse profile consisting of three components: a main pulse and precursor separated by 0.1 rotations and another component that lags the main pulse by about 0.3 rotations (see Fig.1).

Figure 3 .

 3 Figure 3. Gamma-ray and radio pulse profiles for PSR J1600-3053 (see Fig. 1 for details).

Figure 4 .

 4 Figure 4. Gamma-ray and radio pulse profiles for PSR J1713+0747 (see Fig. 1 for details).

Figure 5 .

 5 Figure 5. Gamma-ray and radio pulse profiles for PSR J1741+1351 (see Fig. 1 for details).

Figure 6 .

 6 Figure 6. Gamma-ray and radio pulse profiles for PSR J2051-0827 (see Fig. 1 for details).

Figure 7 .

 7 Figure7. Gamma-ray emission above 3 GeV folded on the orbital period for PSR J0610-2100. The light curve was made using photon weights and it is shown twice for clarity. Vertical dotted lines indicate the phase at which eclipses would be expected and the horizontal dashed line represents the photon background. The H-test for this light curve is low, suggesting that the modulation is not significant.

  and, because of possible contamination of their Ṗ values produced by movement in the gravitational potential of a GC, their B LC values are uncertain.

Figure 8 .

 8 Figure 8. Magnetic field at the light cylinder, corrected for the Shklovskii effect if possible, versus the rotational period for all known pulsars with P ≤ 0.1 s. The MSPs in GCs are shown using thick black dots to indicate that their B LC values are uncertain (see Section 5.2). Lines of constant spin-down energy rate are dashed.

Figure 9 .

 9 Figure 9. Distribution of radio spectral indices for 52 MSPs highlighting the different values for gamma-ray MSPs and A-type MSPs.

Figure 10 .

 10 Figure 10. Distribution of the average levels of radio linear polarization for 18 MSPs in the sample.

Table 2 .

 2 Main pulsar parameters for the six detected MSPs. The first eight columns contain the pulsar names, positions in galactic coordinates, periods, period derivatives, orbital periods, proper motions and distances. The last two columns give the values for the period derivative and spin-down energy rate, both corrected for the Shklovskii effect (P corr and Ėcorr ). Uncertainties in the last quoted digits are in parentheses.

	Pulsar		b	P	Ṗ	P B	μ	d	Ṗ corr	Ėcorr	Refs. a
		(deg)	(deg)	(ms)	(10 -20 )	(d)	(mas yr -1 )	(kpc)	(10 -20 )	(×10 33 erg s -1 )
	J0610-2100	227.8 -18.2	3.86	1.24	0.3	18.2(2)	4(1) b	0.1(3)	1(2)	1, 2
	J1024-0719 c	251.7	40.5	5.16	1.85	-	59.9(2)	0.5(1)	<0.05	<0.1	1, 3
	J1600-3053	344.1	16.5	3.60	0.95	14.4	7.2(3)	2(1)	0.84(4)	7.1(4)	4, 3
	J1713+0747	28.75	25.2	4.57	0.85	67.8	6.30(1)	1.1(1)	0.805(4)	3.33(2)	5
	J1741+1351	37.89	21.6	3.75	3.02	16.3	11.71(1)	1.08(5)	2.89(1)	21.68(4)	6
	J2051-0827	39.19 -30.4	4.51	1.27	0.1	7.3(4)	1.0(2)	1.21(1)	5.21(5)	7, 2
	a The references in the last column are for proper motion and distance: (1) this work; (2) NE2001 model (Cordes & Lazio 2002);
	(3) Verbiest et al. (2012); (4) Verbiest et al. (2009); (5) Splaver et al. (2005); (6) preliminary, Freire et al. (in preparation);
	(7) Lazaridis et al. (2011).
	b See a discussion about this value at the end of Section 2.2.1.
	c See a discussion about the distance to this MSP and the Shklovskii corrections in Sections 2.1.1 and 2.2.1.

Table 3 .

 3 Gamma-ray light curve and spectral properties for the six detected MSPs. The first quoted uncertainties are statistical, while the second are systematic. For PSR J1024-0719, the reported gamma-ray efficiency η is a lower limit, because only an upper limit on the Shklovskii-corrected period derivative, Ṗ , is known (see Sections 2.1.1 and 2.2.1).

	J0610-2100 a	J1024-0719	J1600-3053

https://confluence.slac.stanford.edu/display/GLAMCOG/Public+List+ of+LAT-Detected+Gamma-Ray+Pulsars

http://fermi.gsfc.nasa.gov/ssc/data/access/lat/ephems/

http://fermi.gsfc.nasa.gov/ssc/data/analysis/scitools/overview.html

These diffuse models are available for download from the Fermi Science Support Center; see http://fermi.gsfc.nasa.gov/ssc/data/access/lat/ BackgroundModels.html.

http://heasarc.gsfc.nasa.gov/cgi-bin/Tools/w3nh/w3nh.pl

http://www.physics.mcgill.ca/˜aarchiba/photons_plug.html

AC K N OW L E D G M E N T S

The Fermi LAT Collaboration acknowledges generous ongoing support from a number of agencies and institutes that have supported both the development and the operation of the LAT, as well as scientific data analysis. These include the National Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l'Energie Atomique and the Centre National de la Recherche Scientifique/Institut National de Physique Nucléaire et de Physique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di Fisica Nucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), the High Energy Accelerator Research Organization (KEK) and the Japan Aerospace Exploration Agency (JAXA) in Japan, and the K. A. Wallenberg Foundation, the Swedish Research Council and the Swedish National Space Board in Sweden. We gratefully acknowledge additional support for science analysis during the operations phase from the Istituto Nazionale di Astrofisica in Italy and the Centre National d' Spatiales in France.

The Nanc ¸ay Radio Observatory is operated by the Paris Observatory, associated with the French Centre National de la Recherche Scientifique (CNRS). The Lovell Telescope is owned and operated by the University of Manchester as part of the Jodrell Bank Centre for Astrophysics, with support from the Science and Technology Facilities Council of the United Kingdom. The Westerbork Synthesis Radio Telescope is operated by the Netherlands Foundation for Radio Astronomy, ASTRON. The Parkes Radio Telescope is part of the Australia Telescope, which is funded by the Commonwealth Government for operation as a National Facility managed by CSIRO. We thank our colleagues for their assistance with the radio timing observations. The Arecibo Observatory is operated by SRI International under a cooperative agreement with the National Science Foundation (AST-1100968), and in alliance with Ana G. Méndez-Universidad Metropolitana and the Universities Space Research Association.

Table 4. Radio and gamma-ray properties of 30 MSPs detected by the LAT. The columns show the pulsar name, magnetic field at the light cylinder (B LC ), mean degree of linear polarization L/I , radio spectral index s, radio lag (δ) and type of MSP, according to the classification in Section 5. Uncertainties on the last quoted digit are indicated in parentheses. For the average linear polarization, the errors quoted correspond to a rough and conservative estimate of the standard deviation from several measurements made by different authors, preferably around 1.4 GHz. b MSPs with radio and gamma-ray pulses overlapping in phase but not considered as aligned (Section 5.1.1).

Pulsar