Individual subpulses of PSR B1916+14 and their polarization properties

Individual subpulses of pulsars are regarded as the basic emission components, providing invaluable information to understand the radio emission process in the pulsar magnetosphere. Nevertheless, subpulses are overlapped with each other along the rotation phase for most pulsars, making it difficult to study the statistical properties of subpulses. Among the pulsars observed by the Five-hundred-meter Aperture Spherical radio Telescope, PSR B1916+14 has a large number of isolated well-resolved subpulses in the high time resolution observations, having a typical width of 0.15 ms and a high linear polarization. We find that the number distribution of subpulses contributes dominantly to the mean profile. According to the emission geometry, these emission units come from a region roughly 155 km above the polar cap in the pulsar magnetosphere, and the length scale of basic emission units is approximately 120 m. The deviations of polarization position angles for these single subpulses from the standard S-shaped curve are closely related to their fractional linear and circular polarization, and the large deviations tend to come from drifting subpulses.


INTRODUCTION
The pulsar radio emission is polarized.The mean pulse profile delineates the emission window in the pulsar magnetosphere (Lyne & Manchester 1988), and the polarization position angle (hereafter PA) at each rotation phase corresponds to the plane of curved magnetic field lines according to the Rotation Vector Model (hereafter RVM, Radhakrishnan & Cooke 1969;Komesaroff 1970).An individual pulse consists of many subpulses emerging at different phases with various strengths.Stinebring et al. (1984a,b) studied the statistical properties of single pulses and found the orthogonal polarization modes.Mode-separated pulse profiles are often highly polarized (McKinnon & Stinebring 2000), which implies a depolarization caused by overlapped subpulses with different polarization modes.Distinguishable subpulses (or single subpulses), such as micro-pulses and giant pulses, have been detected from a few pulsars, and these reflect the characters of the basic emission units coming from charged particle bunches.Both micro-pulses and giant pulses have a narrow width smaller than hundreds of microseconds, or even nanoseconds (Hankins & Boriakoff 1978;Hankins et al. 2003;Soglasnov et al. 2004;Jessner et al. 2010;Mitra et al. 2015;Kazantsev & Potapov 2018), and are usually highly polarized.
The accurate polarization measurements of narrow single subpulses give first-hand information on emission processes, as demonstrated by the detected dwarf pulses from PSR B2111+46 in nulling mode (Chen et al. 2023), providing the following observation conditions: 1) a high time resolution enabled by a good digital backend; 2) a good quality of polarization calibration; 3) a high signal-to-★ E-mail: hjl@nao.cas.cnnoise-ratio of pulse detection; 4) a pulsar with a substantial number of well-resolved subpulses.The first three conditions can be satisfied by a large sensitive telescope, but the last one is hard to satisfy because subpulses are often overlapped for most pulsars.
The Five-hundred-metre Aperture Spherical radio Telescope (FAST, Nan 2006) is currently the most sensitive radio telescope for pulsar observations.It has been used for the FAST Galactic Plane Pulsar Snapshot (GPPS) survey (Han et al. 2021), the most efficient pulsar-search survey on the galactic plane currently in the world.Up to now, the survey has detected 637 new pulsars.During test observations for the survey, a large number of known pulsars have been detected.We examined high-time-resolution profiles of individual pulses and found that PSR B1916+14 is the best to have resolved single subpulses.It is a bright pulsar discovered by Hulse & Taylor (1974), and has a period of P 0 = 1.181s and DM = 27.202pc cm −3 (Hobbs et al. 2004b) and RM = −41.7 ± 3.5 rad m −2 (Han et al. 2018).The mean pulse profile is very narrow with a full width of 28 ms at the half maximum at 1.3 GHz, occupying only 2.4 percent of its period (Hobbs et al. 2004a).This narrow profile has three predominant components (Blaskiewicz et al. 1991;Rankin 1993;Weisberg et al. 1999).The emission geometry has been determined by Olszanski et al. (2019) by fitting the observed PA with a standard S-shape PA curve according to the RVM.The inclination angle of the magnetic field axis to the rotation axis is  = 79 • and the impact angle, the closet approach between the line of sight and the magnetic axis, is  = 1.2 • (Rankin 1993;Olszanski et al. 2019).Some subpulses are found to drift, with a drifting rate of  3 = 59 periods and the phase-shift in a period being  2 = −3.9• (Song et al. 2023).
In this paper, we present the sensitive polarization observations for individual subpulses of PSR B1916+14 observed by the FAST during  the GPPS survey.In Section 2, the FAST observations and processing procedures are briefly introduced.In Section 3, individual subpulses and their polarization properties are statistically studied.The related physics of emission and propagation processes are discussed, and conclusions are given in Section 4.

FAST DATA AND PROCESSING
The FAST observations for PSR B1916+14 were made on May 28th, 2020 (hereafter 2020-05-28) during the regular survey observations of the GPPS survey (Han et al. 2021) by using the 19-beam L-band receivers (Jiang et al. 2020).The pulsar is located in the 11th beam but 1.1 ′ away from the center, see Table 1.To verify the results on single subpulses, one more FAST session was carried out for 15 minutes on March 8th, 2023 (hereafter 2023-03-08) by using the central beam of the 19-beam receiver.The FAST has a full gain of G = 16 K/Jy when it points to any object with a zenith angle smaller than 28.5 • .The receivers cover the frequency range from 1.0 to 1.5 GHz, and the data from all 19 beam receivers are recorded for four polarization channels of  ,  , ( *  ) and ( * ) with a time resolution of 49.152 μs.In case there is any known pulsar in the targeted region for a GPPS survey observation, data are recorded for two extra minutes but with the injected calibration signals of 1 K On-Off noise with a period of 2 s, which are used to calibrate the band-pass and polarization performance of the receiver.
The raw data are saved in a search mode with 2048 channels over the observation band.Using the ephemeris from the Australia Telescope National Facility Pulsar Catalogue (Manchester et al. 2005), the data are de-dispersed and folded by using the open source package dspsr (van Straten & Bailes 2011).To uncover the details of individual pulses in the 261/771 periods of the observations on 2020-05-08/2023-03-08, we fold data with the original sampling time of 49.152 μs (24030 bins).To get the mean profile, observational data are also folded with 2048 bins per period, as shown in Figure 1.Careful data analysis processes were accomplished by the pulsar toolkit, psrchive (van Straten et al. 2012).For example, by using the command psrzap, the RFIs in folded data are cleaned in the two-dimensional frequency-time images interactively and the relative channels are weighted to zero.With pac, the polarization leakage is corrected and the frequency bands are calibrated by using the solutions derived from the 2-minute calibration data.We get the observed rotation measurements from two FAST observation sessions via rmfit,  = −36.15± 0.14 rad m −2 on 2020-05-08 and  = −34.56± 0.10 rad m −2 on 2023-03-08.After the RM contribution from the earth ionosphere is estimated via IonFR (Sotomayor-Beltran et al. 2013), the true  on 2020-05-08 and 2023-03-08 are −36.57± 0.17 and −37.04 ± 0.21 rad m −2 respectively, consistent with but much improved from the previous value (Han et al. 2018).

FAST OBSERVATION RESULTS
All polarization profiles, both the mean profile and individual pulses, are obtained after removing the effect of RMs and averaging over all channels by using pam.

The mean polarization profile and emission geometry
The mean polarization profile of PSR B1916+14 observed by FAST on 2020-05-28 and 2023-03-08 are shown in Figure 1, which are consistent with that in Weisberg et al. (1999) observed by using the Arecibo telescope.Three components are predominant in the mean profile, marked as C 1 , C 2 and C 3 for the leading, central, and trailing components in the phase range of [-5.5 respectively.The fractional linear/circular polarization, / and / of these three components are 73% and 6%, 53% and 15%, and 48% and 5%, respectively  is plotted with parameters of =79 • and =1.2 • given by Rankin (1993) and Olszanski et al. (2019).The fitted RVM curve can be taken as a standard PA curve.A small number of bins around the longitude of −4 • in the panel (b) of Figure 1 show the orthogonal mode.
The emission beam of PSR B1916+14 complies with the cone+cone model for three reasons: (1) the width evolution of the central profile component with frequency is consistent with that of the cone (see Fig. A18 of Rankin et al. 2023), indicating that the central profile component comes from conal emission; (2) the subpulse drifting (see Section 3.4) detected from the central component is a common property of conal emission (Rankin 1986;Rankin et al. 2023); (3) no sign reversal of circular polarization is detected for the central component, indicating that the central profile id not the core emission from the beam center Lyne & Manchester (1988).Therefore, the three predominant components of the mean profile of PSR B1916+14 should be produced by the line of sight grazing the inner and outer cones of the emission beams, as discussed by Olszanski et al. (2019).
The emission height of pulsars,  em , could be estimated by (Lyne et al. 2022): Here  b is the angular radius of the outer cone of the emission beam, and  10 represents the pulse width at the level of 10% of the peak.
The coefficient 2 describes the farthest point on the last open field line in the unit of the light cylinder radius.For PSR B1916+14,  10 = 7.57 • ,  lol = 1.33 for =79 • .Taking these geometric parameters into the above equation, we get the size of the emission beam,  b ≃ 3.9 • for the outer cone, and then the emission height  em = 155 km.

Morphological parameters of subpulses
We obtained high time-resolution profiles of individual pulses for 261 and 771 periods of PSRB1916+14 (see Appendix A).As shown in Figure 2 for two examples, an individual pulse consists of many individual subpulses, which are well-separated and act as the isolated emission units.We checked individual pulse profiles for many bright pulsars observed by the FAST and found that PSR B1916+14 is the only pulsar having so large number of well-distinguished subpulses.Therefore this is an excellent chance to study the properties of individual subpulses radiated by the individual clusters of particles in the pulsar magnetosphere.
For this purpose, we fit individual sub-pulses with the multi-Gaussian function and finally get 23321 Gaussian fittings to subpulses, among which there are 955 single subpulses (see Sec. 3.3).The morphological parameters, such as the peak phase in the longitude, the pulse energy (), and the full-width at half maximum w sp , are obtained and published online (see Table A1 in Appendix B).
The distributions of these morphological parameters for the 23321 subpulses are given in Figure 3.These subpulses have a typical width w sp in the range from 0.05 ms to 1.0 ms and peak at about 0.15 ms (solid line).For a typical width  sp = 0.15 +0.85  −0.10 ms, the typical length scale of such basic emission units in the magnetosphere is estimated as  em sin  2 sp / 0 = 0.12 +0.68 −0.08 km.The scattering time of this low DM pulsar is 0.03 μs, and the dispersion smearing time of a channel at the central frequency of 1.25 GHz is about 27 s, both are negligible.
The distribution of these subpulses along the rotation phase (the solid) is shown in Figure 3b, roughly following the same function as the mean profile.More subpulses are detected within the  2 component near the profile center than others.The histograms of the subpulses energy () in the logarithm format with a base of 10 for three components are shown in Figure 3c, and they are fitted well by a log-normal distribution function, which suggests a linear wave growth process (Cairns et al. 2004): Finally, the fitting parameters are shown in   b) and their number distributions (subpanels a and c).In subpanel (b), the symbol size is proportional to the PA deviation of a single subpulse from the fitted mean PA curve, "×" for positive offset and the circles for negative deviations.
pulse energy is only a factor of two compared to those of the other two components.However, the number of subpulses represented by a 0 in Table 2 dominates more effectively to the strength of this mean profile component.

Polarization properties of single subpulses
We pick up well-separated single-subpulses to avoid any confusion of overlapped subpulses for better statistics, as demonstrated in Figure 2 for the 110th period observed on 2020-05-28.The single-subpulses of PSR B1916+14 are selected with the following conditions: (1) peak intensity greater than 20 I , so that the subpulse is strong enough to ensure a high confidence for the polarization parameter determination; (2) subpulses nearly isolated as chosen by the algorithm introduced in Appendix B. Those with obvious shoulder components such as the No. 110-5 in Figure 2 are discarded, because the shoulder components may affect the width measurements and polarization measurements for the concerned peak.In the end, from individual pulses of 1032 periods of PSR B1916+14, we get 955 single subpulses with high-quality polarization measurements.
The single-subpulses, as basic emission units, often have strong linear polarization (see Figure 4).More intriguing is that the large PA deviations from the standard PA curve are seen from some single-  2. To explore the polarization properties and their interplay for these individual subpulses as basic emission units, we measure the polarization parameters for all singlesubpulses, as listed in Table A1 in the Appendix B, including the peak flux density  p , linear polarization  p , circular polarization intensity  p (positive for the left-hand sense) and the absolute value | p |.For each of single-subpulses, data bins in the longitude range of w sp are summed to get the integrated flux density for a subpulse on , ,  and | |.Other derived parameters, linear polarization percentage (hereafter, /), and circular polarization percentage (hereafter, /) are all derived for the integrated flux density for every single subpulse.

subpulses, see examples in Figure
As seen in Figure 4, most single subpulses are highly linearly polarized with weak circular polarization, even though a small number of them are highly circularly polarized, as indicated by the histograms shown in the subpanels (a) and (c).Their histograms peak at /=90% and /=10%.The PA deviations are somehow related to / and also /.We then check observation data in Figure 6, see if there is any dependence of Δ on , , and .
The PA deviation (Δ) of a single subpulse from the standard S-shape curve is the weighted average for values in these bins, calculated by considered.Such a weighted average naturally rests on the values in better-measured bins inside a subpulse.We analyzed the Δ distribution of all the single-subpulses and found that it fits an exponential distribution quite well as for Δ > 0, Δ < 0 or |Δ| as shown in Figure 5.Here the scale of the exponential distribution Δ scale reflects the dispersion degree of PA deviations.The single-subpulses with Δ < 0 and > 0 have Δ scale = 4.9 • ± 0.2 • and 4.8 • ± 0.4 • respectively.For all the singlesubpulses, Δ scale = 4.8 • ± 0.3 • .The PA deviations spread over a wide range in ±90, and show non-dependence on  and  (see Figure 6a and c), but are bigger for these single subpulses with weaker linearly polarized intensity (see Figure 6b).We further check the dependence of Δ on / and /.The PA deviations become more dispersed as / becomes smaller (see Figure 6d), in a similar style for positive and negative Δ.But just in contrast, the more dispersed Δ data are found for stronger circular polarization (see Figure 6e and f).We find that each of the Δ distributions is symmetrical with Δ = 0 • , they are folded along Δ = 0 • for simplicity (see Figure 6g to i).Quantitatively, we divide all the points into 4 equal portions in order of / and/or |/|, and get their Δ scale (the red boxes in Figure 6g and i).The scales of PA deviation distribution are linearly related to fractional linear polarization, Δ scale = (−13.3± 2.0) × / + (16.5 ± 1.5), and to fractional circular polarization Δ scale = (18.0±1.4) × |/ | + (2.5± 0.3).

Subpulse drifting and PA deviation
From the FAST data of PSR B1916+14, we find that subpulse-drifting in the components  2 and  3 , as shown in Figure 7. Through the twodimensional Fluctuation Spectra method (2DFS, Edwards & Stappers 2002), we obtained a drift periodicity  3 =58.30 and a longitude spacing of  2 =−3.7 • , consistent with the results given by Song et al. (2023).The large  3 among the known driftings may imply a low acceleration potential drop of this pulsar according to the partially screened gap model (Gil et al. 2003;Basu et al. 2016).
Interestingly, the single subpulses with a large Δ appear mainly in the drifting  2 and  3 components (see the lower panel of Figure 7).The corresponding emission in these two components is relatively strong.In physics, the plasma environment in the emission region of the  2 and  3 components should be denser and more complicated compared to that of the  1 component.We suspect the large PA deviation from the PA curve of individual subpulses may relate to the complicated plasma environment.

DISCUSSION AND CONCLUSION
Based on FAST observations, PSR B1916+14 is found to have a very well S-shape PA curve of the mean profile, consistent highly with the standard PA curve predicted by the RVM model with  = 79 • and  = 1.2 • (Olszanski et al. 2019).Its emission beams comply with the cone+cone model and the emission height of the outer cone is estimated to be 155 km and the size of the outer cone  b ≃ 3.9 • .In addition, the orthogonal modes are detected for a only very small number of subpulses.This pulsar radiates almost in a single natural wave mode (ordinary or extraordinary mode).
The basis properties of individual subpulses have been investigated by a large number of distinguishable single subpulses from the FAST super-sensitive high-time resolution observations.They have a typical width of W sp ∼ 0.15 ms, corresponding to the length scale of 0.12 km.Through investigating the distribution of phase locations and pulse energy, we found that both contribute to the mean profile, but the number of subpulses is the dominant factor.
We find that the PA deviations of single subpulses from the standard PA curve are related to the fractional linear and circular polarization.The PA deviations are more dispersed when the fractional circular polarization becomes larger and the fractional linear polarization becomes smaller.Interestingly, the large PA deviation usually is detected from single subpulses in drifting  2 and  3 components.There are four possible explanations for the PA deviations: (1) Since the subpulses may produce at different emission heights, their PA deviations from the standard S-shape curve can be caused by the aberration and retardation effects and the possible plasma corotation effect (Blaskiewicz et al. 1991).Taking into account these effects, the PA dispersion around the S-shape curve may be explained.However, it fails to relate the PA deviation with the circular polarization.
(2) Incoherent mixing of orthogonal modes may also cause the PA deviation accompanied by the strong de-polarization for both linear and circular polarization (Mitra et al. 2023).This is compatible with the weakly-polarized subpulses with a large PA deviation (see Figure 4).However, it can not explain these with strong circular polarization.
(3) According to Wang et al. (2022), once the departure of the line of sight from the radiation beam center reaches 1/ or more where  is the Lorentz factor of plasma, the significant circular polarization would be created and thus the calculated PA deviates from the value given by RVM curve (see Equations ( 29), ( 30) and (31) in Wang et al. 2022).The fractional circular polarization would be high when the departure gets large enough, however, the intensity of an emission unit would decrease dramatically at the same time.That is not consistent with the finding of strong single subpulses highly-circularly polarized for example No.110 in Figure 2. High fractional circular polarization can achieve hardly in the usual views like Wang et al. (2012).
(4) Propagation effects could be a good scenario for explaining the large PA deviation and high circular polarization for the single subpulses.According to Wang et al. 2010, if the impact angle  is enough low, 1.2 • for PSR B1916+14, a large PA deviation would be generated accompanied by high circular polarization with a single-handedness through wave mode coupling effect by which a natural wave evolves from adiabatic to non-adiabatic in the streaming plasma along its path.As discussed in Section 3.4, strong subpulses tend to appear in the drifting  2 and  3 components.That indicates a circumstance with many high-density plasma bunches, which are the basic units producing single subpulses, in the magnetosphere of the  2 and  3 components.When a natural wave propagates through these plasma bunches, evolving from non-adiabatic to adiabatic and finally to nonadiabatic, its large PA deviation and large circular polarization are produced.On the contrary, because there is a low-density plasma situation in the magnetosphere of the weak  1 component, which leads to a weak propagation effect, the single subpulses are weaklycircularly polarized and their PAs follow the standard S-shape curve.It should be noted that a sense reversal of circular polarization could happen through the quasi-tangential effect when the radiation points are close to the conical surface scanned by the magnetic axis around the spinning axis (Wang & Lai 2009).Nevertheless, the propagation processes have been discussed traditionally based on the plasma background before.Comprehensive discussions are needed for the question of how radio emission passes through many high-density plasma bunches.
In short, large PA deviation of individual subpulses can be generated by four physical mechanisms.All of them could contribute.However, the propagation effect in the pulsar magnetosphere has an advantage in explaining the dependence of PA deviation on circular polarization.

Figure 1 .
Figure 1.The mean polarization profiles and the data distribution of polarized emission in the pulse window are obtained from the two observations on 2020-05-28 and 2023-03-08.In panel (a), the mean profile for the total intensity , linear intensity , and circular polarization  are presented.The longitude ranges for the 3 mean profile components ( 1 ,  2 , and  3 ) are marked.In panel (b) the data distribution of PAs for each bin of individual pulse (2048 bin per period) is plotted only when the PA has an uncertainty less than 5 • (calculated via  Ψ = ( 2  2  +  2  2  ) 1/2 /(2 2 ).The PA curve of the mean profile (dotted line) is plotted together with the RVM fitting (dashed line).The panels (c) and (d) show the distribution of / and / for individual bins with  > 20 I .

Figure 2 .
Figure 2. The polarization profiles of two example individual pulses (the No.110 and 111 periods observed on 2020-05-28) with a time resolution of about 49 μs from the super-sensitive FAST observations, with well-separated and resolved "subpulses" as exampled in the 6 lower panels.The total power , linear polarization intensity , and circular polarization intensity  are shown and compared with the enlarged version (x12) of the mean polarization profiles plotted in a fainter version.The dotted line is plotted for 20 I that is obtained from off-pulse data.The PA values of the individual pulses are shown in the upper subpanels if the uncertainty is less than 5 • .The PA-fitted S-shape curve is plotted for comparison.The polarization profiles for these strong individual subpulses are highlighted for  > 3 I .All data for the high-time-resolution profiles for individual pulses are presented in the Appendix.

Figure 3 .
Figure 3. Statistical properties of subpulses obtained by the FAST.(a) The histogram of subpulse width W sp for all (solid boxes) and recognized single (dashed boxes) subpules.(b) The phase distribution for all (solid boxes) and single (dashed boxes) subpulses, compared with the mean pulse profile (solid grey line) in Fig. 1.(c) The histograms of subpulse energy.The solid, dotted, and dashed histograms represent those of C 1 , C 2 , and C 3 components with their fitting lines, respectively.

Figure 4 .
Figure 4.The distribution of / versus / for single subpulses (subpanelb) and their number distributions (subpanels a and c).In subpanel (b), the symbol size is proportional to the PA deviation of a single subpulse from the fitted mean PA curve, "×" for positive offset and the circles for negative deviations.

Figure 6 .
Figure 6.Distributions of PA deviations (Δ) against (a) , (b) , (c) , (d) /, (e) /, and (f) | |/.The sub-panels (g), (h), and (i) show the distribution of absolute PA deviation against /, /, and | |/ respectively.The symbols "×" and the circles represent the positive and negative PA deviations.The mean of PA deviations (Δ scale ) are obtained for 4 equal data numbers (the red boxes) for fitting them separately, which are shown in the subpanel (g) and (i).

Figure 7 .
Figure 7. Observational pulse sequences on 2020-05-28 (in the panel a) and 2023-03-08 (in the panel b).The single subpulses are marked with '×' (Δ > 0) and '⃝' (Δ > 0) in both two panels for PA deviations, and symbol size is proportional to the magnitude of their |Δ|.Furthermore, the Δ distribution is shown in panel (c) with the mean pulse profile in a grey line.

Figure B1 .Figure B2 .Figure B3 .
Figure B1.The multi-Gaussian fittings to the individual pulse of pulsar period No.66 on 2020-05-28, as an example.The top panel shows the total intensity (solid line) with its 22 Gaussians (dashed lines), and the bottom panel exhibits the fitting residual.The dotted lines represent 5  .

Table 1 .
FAST observation sessions for PSR B1916+14.The offset is the pulsar position from a beam center.

Table 2 .
The Gaussian fitting to the log-normal distribution of pulse energy for the three components: the number , and the mean energy, ε, for the three components.

Table 2
. The  2 component has the biggest values of  and , which implies that subpulses tend to get stronger in the central component than in others.The sub-Num.