A&A 402, 827-836 (2003)
DOI: 10.1051/0004-6361:20030279

Signatures of high energy protons in pulsar winds

E. Amato1 - D. Guetta1 - P. Blasi1

INAF/Osservatorio Astrofisico di Arcetri, Largo E. Fermi 5, 50125 Firenze, Italy

Received 4 July 2002 / Accepted 4 February 2003

The resonant cyclotron absorption model (Hoshino et al. 1992; Gallant & Arons 1994) is very successful in describing particle acceleration in plerions. A prediction of this model is the presence of a substantial amount of relativistic protons in pulsar winds. Although difficult to detect, these protons may show up through their interactions either with the photons in the plerion environment or with the thermal gas in the supernova ejecta. Inelastic proton-proton (p-p) collisions are expected to be very effective in young objects, resulting in a copious production of neutral and charged pions. Charged pions produced during the first few hundred years after the supernova explosion may have time to decay into muons, whose subsequent decay may provide an additional source of electrons and positrons in these nebulae, that adds up to the pulsar input. These secondary leptons evolve just as the pairs in the pulsar wind, and signatures of their presence could be found, in principle, even in the synchrotron spectrum of older objects. p-p collisions may remain fairly efficient even in moderately old objects resulting in the production of TeV $\gamma$-rays and neutrinos. We apply our calculations to the case of the Crab Nebula, the best studied plerion thus far, and find that existing data already allow us to infer interesting constraints on the physical properties of the Crab pulsar wind.

Key words: acceleration of particles - stars: pulsars: general - stars: winds, outflows

1 Introduction

The presence in pulsar winds of a baryonic component coexisting with the leptonic one is still an open question. According to the most successful model so far for particle acceleration at the pulsar wind termination shock (Hoshino et al. 1992; Gallant & Arons 1994), baryons should not only be present in the wind, but even carry most of its energy. This makes it mandatory to investigate all possible ways in which signatures of their presence can be found.

Plerions result from the interaction of a pulsar with its surrounding medium. According to the standard picture, pulsars lose most of their energy in the form of a magnetized relativistic wind, whose material component is made up mostly of electron-positron pairs. The wind starts out cold and with a highly relativistic expansion velocity. At some distance from the pulsar the wind then needs to be slowed down in order to match the boundary conditions imposed by the surrounding supernova ejecta expanding at subrelativistic speed. This happens at a termination shock where most of the outflow bulk energy is randomized. Here is where particle acceleration presumably occurs, giving rise to the non-thermal particle spectrum that is the source of the plerion synchrotron and Inverse Compton radiation.

According to both pulsar magnetospheric theories and the available data, the termination shock is most likely a relativistic transverse shock, meaning that the magnetic field is almost perpendicular to the flow velocity. Standard shock acceleration mechanisms (i.e. diffusive shock acceleration or Fermi I process, and shock drift acceleration) are generally believed to be ineffective in this kind of shock (e.g. Hoshino et al. 1992).

An acceleration mechanism that is able to overcome the difficulties of the former two is based on resonant cyclotron absorption (RCA) in a baryon-loaded plasma. Here we briefly review it in order to make clear the assumptions we make on the proton spectrum and its relation to the basic wind parameters.

Within the framework of the RCA model (Hoshino et al. 1992; Gallant & Arons 1994), the ions present in the wind have a Larmor radius $r_{{\rm L}_{\rm i}}$ which is $A m_{\rm p}/Z m_{\rm e}$ times larger than that of the pairs, both species being cold in the wind rest frame. When the heavier particles meet the shock in the pairs, this appears to them as an infinitesimally thin discontinuity in the magnetic field. The sudden enhancement of the magnetic field is expected to cause them to start a bunched gyration, accompanied by the collective emission of cyclotron waves. The fundamental frequency of these waves is the ion Larmor frequency, but also higher harmonics are excited up to the pair gyration frequency. The latter can be resonantly absorbed by the pairs and cause their acceleration, while the ions slowly thermalize after crossing the shock. Efficient acceleration requires the ions to be energetically dominant in the wind, although the pairs are dominant by number.

After crossing the shock the ions drift towards the outer parts of the nebula with a velocity of the order of $c\ r_{{\rm L}_{\rm i}}/R_{\rm N}$ and slowly thermalize to a temperature that is equivalent to a Lorentz factor $\Gamma $ of the order of the initial wind Lorentz factor. The full thermalization takes several gyration periods, so, depending on the initial wind Lorentz factor, it may never be achieved within the nebula.

While drifting towards the edge of the nebula the ions may interact with the plerion radiation field and with the thermal material of the ejecta, which partly penetrates the synchrotron bubble due to Rayleigh-Taylor instabilities. Both these interactions would lead to pion production. The products of the subsequent pion decays could lead to observable signatures if the interaction rate were high enough. The most direct signature would be the instantaneous emission of high energy photons and neutrinos.

The latter are potentially the preferential channel to look at, since their detection would be unequivocally associated with the presence of protons. In fact, a few recent papers have considered the question of neutrino production from protons accelerated in pulsar magnetospheres. Bednarek & Protheroe (1997) computed the flux of high energy neutrinos and $\gamma$-rays that would be expected from nuclear collisions of protons derived from the decay of neutrons produced by the photodisintegration of accelerated nuclei in a pulsar magnetosphere. Beall & Bednarek (2002) computed the neutrino flux that would be produced by high energy protons accelerated in the magnetosphere of a young fast rotating pulsar surrounded by a supernova envelope, and Bednarek (2001) considered the effect of this contribution on the overall neutrino background.

None of these previous investigations explicitly considered the protons that should be part of the pulsar wind if the RCA model for shock acceleration is actually at work in plerions. These are exactly the protons we consider in this paper, with the aim of finding constraints from observations on the characteristics of the pulsar wind. In addition to computing the multi-TeV neutrinos and photon fluxes we also worry about the secondary electrons and positrons produced by $\pi^\pm$ decay. If protons are present, this process acts as an extra source of pairs in the plerion which is expected to be especially active when the plerion is young and compact, but could in principle leave relatively long-lasting traces in the nebular radiation spectrum.

In the following we study the efficiency of the pion production processes as a function of time during the plerion evolution. We approximate the plerion as spherical and expanding at a constant velocity. For the evolution of the nebular non-thermal emission we use the Pacini & Salvati (1973) model (PS hereafter) for a homogeneous plerion.

2 Pion production and energy losses

Charged and neutral pions can be generated either by photo-pion production off the photons in the nebula or by inelastic p-p scattering against the thermal gas in the ejecta. While the decay of neutral pions results in $\gamma$-ray production, the decay of charged particles results in neutrinos and electron-positron pairs (we will refer to both as electrons).

We describe the efficiency of both photomeson production and p-p collisions through a parameter, $f^\pi $, defined as the probability for a single relativistic proton to create a pion before losing energy through other kinds of interactions:

f^\pi={\rm min} \left[1, \frac{t_{\rm life}}{t^\pi} \right] \ ,
\end{displaymath} (1)

where $t_{\rm life}$ is the proton life-time against losses different from p-$\gamma$ and p-p, and $t^\pi$ is the timescale for the process of pion production. The luminosity output in pions will be:

L_\pi=f^\pi L_{\rm p}\ ,
\end{displaymath} (2)

with $L_{\rm p}$ the wind luminosity in protons.

In the presence of an intense magnetic or radiation field, charged pions may lose a considerable fraction of their energy through synchrotron or Inverse Compton emission before decaying into muons. The same could happen to muons before decaying into electrons and neutrinos. If any of these effects were relevant, this would lead to a suppression of the neutrino and electron fluxes at high energies. We shall define the quantities $f^\mu $ and $f^\nu $ as:

f^\mu={\rm min} \left[1,{{t^\pi}_{\rm losses} \over
..., {{t^\mu}_{\rm losses} \over
{t^\mu}_{\rm decay}} \right]\ ,
\end{displaymath} (3)

where $t^\pi_{\rm decay}=2.6 \times 10^{-8}\ \gamma_\pi\ {\rm s}$ and $t^\mu_{\rm decay}=2.2 \times 10^{-6}\ \gamma_\mu \ {\rm s}$. Since the pion energy is roughly evenly shared between the final products of its decay, we can write the number of neutrinos, electrons and $\gamma$-rays injected in the nebula per unit time and unit energy interval as:
                                    $\displaystyle J_\nu(\epsilon_\nu)$ $\textstyle \approx$ $\displaystyle 4\ \left(1+f^\nu[3\epsilon_\nu]\right)\ f^\mu[4\epsilon_\nu]\ \le...
... \over {\rm d}E_{\pi^\pm} {\rm d}t} \right\vert _{E_{\pi^\pm}=4\epsilon_\nu}\ ,$ (4)
$\displaystyle J_{\rm e}(\epsilon_{\rm e})$ $\textstyle \approx$ $\displaystyle 4\ f^\mu[4\epsilon_{\rm e}]\ f^\nu[3\epsilon_{\rm e}]\ \left. {{\...
... d}E_{\pi^\pm} {\rm d}t} \right\vert _{\epsilon_{\pi^\pm}=4\epsilon_{\rm e}}\ ,$ (5)
$\displaystyle J_\gamma(\epsilon_\gamma)$ $\textstyle \approx$ $\displaystyle 4\ \left. {{\rm d}N_{\pi^0} \over {\rm d}E_{\pi^0} {\rm d}t} \right\vert _{E_{\pi^0}=2\epsilon_\gamma} \cdot$ (6)

An important point to be noticed in the above set of equations is that in writing Eq. (4) we take into account only muon neutrinos, since these are the particles that neutrino telescopes are able to detect. The possible effects of neutrino oscillations on the results of our calculations will be discussed in Sect. 3.3.

2.1 Photo-meson production

In the case of photo-meson production the target for high energy protons is the plerion emission. The fractional energy loss rate of a proton with energy $E_{\rm p}$ (= $ \Gamma~ m_{\rm p}~ c^2$) due to pion production results in (Waxman & Bahcall 1997):
                                 $\displaystyle t_{{\rm p} \gamma}^{-1}(E_{\rm p})$ $\textstyle \simeq$ $\displaystyle \frac{2^{p+1}}{p+2}\
\sigma_{\rm peak}\ \xi_{\rm peak}\ \frac{\De...
...}} \left(\frac{d}{{R}_{\rm N}} \right)^2\
\frac{ F_\nu(\nu_{\rm p})}{4\ \pi\ h}$  
  $\textstyle \simeq$ $\displaystyle 2.5 \times 10^{-16} \frac{2^{p+1}}{p+2}
\left( \frac{d_{\rm kpc}} {R_{\rm pc}} \right)^2
F_{\nu_{\rm p}}[{\rm mJy}]\ {\rm yr}^{-1}.$ (7)

Here we have treated the plerion as homogeneous and used the fact that the photon spectrum is a power law, $F_\nu \propto \nu^{-p}$. We have also made the approximation that the main contribution to pion production comes from photon energies $\epsilon_\gamma \approx \epsilon_{\rm peak}=0.3$ GeV, where the p-$\gamma$ cross section peaks due to the $\Delta$ resonance. The numerical values are obtained using: $\sigma_{\rm peak} = 5 \times 10^{-28}\ {\rm cm}^2$, $\xi_{\rm peak} =0.2$, $\Delta \epsilon=0.2$ GeV, $\nu_{\rm p}=\epsilon_{\rm peak}/ (\Gamma\ h)$ and $\beta_p \simeq 1$.

2.2 Inelastic nuclear collisions

The energy loss-rate of a relativistic proton due to inelastic nuclear collisions can be estimated as

t_{\rm pp}^{-1} \approx \zeta\ n_{\rm t}\ \sigma_0\ c\ ,
\end{displaymath} (8)

where $n_{\rm t}$ is the target density, $\sigma_0=5 \times 10^{-26}~ {\rm cm}^2$and $\zeta\simeq 20$% is the average fraction of energy lost by the proton.

As we mentioned in the introduction, the estimate of the effective target density, $n_{\rm t}$, is the main uncertainty in our calculation. The distribution of the thermal material, in all cases when it is observed, is actually far from uniform. The matter of the ejecta penetrates the synchrotron nebulae due to Rayleigh-Taylor instabilities, giving rise to concentrations in the form of "filaments''. The total amount of material contained in these filaments is not easy to estimate with confidence. In fact the possibility that the gas detected through thermal emission hides some colder and denser concentrations cannot be excluded. A comparably difficult task is to estimate, from the projection on the plane of the sky, the fraction of the nebular volume that these filaments occupy, their so-called "filling factor''.

Moreover the effective target density for relativistic protons, which is our reason for worrying about the distribution of thermal matter in plerions, is also affected by the intensity and structure of the magnetic field in and around the "filaments'', which are not known. The magnetic field affects in fact the physics of proton propagation and could in principle be an obstacle to the penetration of protons in the thermal gas, leading to a decrease of the value of $n_{\rm t}$, or conversely be the cause of proton trapping in the filaments, resulting in the opposite effect. We shall summarize all these unknowns in a single parameter $\mu $, which we also take to be time-constant. We define the effective density of material in terms of the average gas density in the plerion as:

n_{\rm t}=\mu\ \bar{n}_t=\mu\ \frac{M_{\rm N}}{m_p}\ \frac{3...\ \frac{M_{{\rm N} \odot}}{{R}_{\rm pc}^3}\ {\rm cm}^{-3}\ ,
\end{displaymath} (9)

where $M_{\rm N}$ is the mass of the ejecta contained in the filaments, which we also take to be time-constant. Using Eq. (9) in Eq. (8) we obtain:

t_{\rm pp}^{-1} \approx 10^{-7}\ \mu\ \frac{M_{{\rm N} \odot}}{{R}_{\rm pc}^3}\ {
\rm yr}^{-1}\ ,
\end{displaymath} (10)

which, when compared with Eq. (7), indicates nuclear collisions as the main mechanism for pion production in plerions, if $\mu $ is not very much less than 1.

2.3 Energy losses

The energy distribution of particles after injection is determined, in a plerion, mainly by the following loss-mechanisms: escape and expansion losses, and radiative cooling through synchrotron and Inverse Compton emission. We approximate the plerion as a homogeneous spherical bubble filled with a uniform magnetic field and expanding at constant velocity v, and estimate the time-scales of the different processes under these assumptions. The rate of adiabatic cooling is simply:

t_{\rm exp}^{-1}= \left(\frac {v} {R_{\rm N}}\right) \simeq 10^{-3} \frac{{V}_3}{{R}_{\rm pc}}\ {\rm yr}^{-1}~
\end{displaymath} (11)

where ${R}_{\rm pc}$ is the nebular radius in parsecs and V3is the plerion expansion velocity in units of 103 km s-1.

The rate of escape of a particle of mass $m_{\rm s}$ and Lorentz factor $\Gamma^s$, whose Larmor radius in the nebular magnetic field is $r_{\rm L}$, can be estimated, based on $\vec \nabla B$ drift, as:

{t^s_{\rm esc}}^{-1}= \frac{c\ {r}_{\rm L}}{{R}_{\rm N}^2}\s...
...\ {R}_{\rm pc}^2}\ \frac{m_{\rm s}}{m_{\rm p}}\ {\rm yr}^{-1},
\end{displaymath} (12)

where B-4 is the nebular magnetic field in units of 10-4 G and  $\Gamma^s_6$ is the particle Lorentz factor in units of 106.

The rate of synchrotron cooling is:

                               $\displaystyle {t^s_{\rm sync}}^{-1}$ = $\displaystyle {4 \over 3}\ \sigma_{\rm T}\ c\ \left({m_{\rm e} \over m_{\rm s}}\right)^2\ {\Gamma^s \over m_{\rm s} c^2}\ U_B$  
  $\textstyle \simeq$ $\displaystyle 7 \times 10^{-14}\ \Gamma^s_{6}\ B_{-4}^2 \left( {m_{\rm p} \over m_{\rm s}} \right)^3 {\rm yr}^{-1}\ .$ (13)

Finally, the cooling rate due to Inverse Compton Scattering (ICS hereafter) is
$\displaystyle {t^s_{\rm ICS}}^{-1}= { 4 \over 3} \sigma_T c\ \left({m_{\rm e} \...
...\nu_{\rm KN}}^{\nu_{\rm max}}
{S_\nu \over \nu^2} {\rm d}\nu \right] \Bigg\}\ ,$     (14)

where $\nu_{\rm KN}=m_{\rm s}c^2/\Gamma^s/h$ is the transition frequency between interaction in the Thompson and Klein-Nishima regimes, $U_{\rm cmb}$ is the energy density of the cosmic microwave-background photons, and $S_\nu$ is the volume emissivity of the plerion.

3 Proton signatures during the plerion history

3.1 Plerion evolution

All the rates in the previous section depend on quantities that change with time during the plerion evolution. We determine the dependency on time of the different quantities within the framework of the PS model. We consider the spin-down power of the pulsar to be continuously converted into injection of fresh particles and magnetic field energy. The spin-down luminosity of the pulsar can be parametrized as:

L(t)={L_0 \over (1+t/\tau)^\alpha}\ ,
\end{displaymath} (15)

where the characteristic spin-down time $\tau$ is a few hundred years, and the index $\alpha$ is expected to be around 2 if the breaking is purely electromagnetic.

Treating the plerion as homogeneous and expanding at constant velocity, the magnetic field evolution, resulting from the competition between the continuous pulsar supply and adiabatic losses, can be determined analytically. For times shorter or comparable with $\tau$, one can make the approximation:

B(t)=\sqrt{3 L_B \over v^3 t^2}\ ,
\end{displaymath} (16)

where LB is the fraction of the pulsar luminosity that goes into magnetic energy.

The other quantity we need to know as a function of time in order to estimate ICS losses and the rate of photo-meson production is the plerion volume emissivity $S_\nu$, due to non-thermal radiation. To this purpose, we need to determine the time-dependences of the spectral density N(E,t) of the emitting particles (pairs). For the wind particles we assume that they are injected with an initial spectrum in the form of a broken power-law:

J_w(E)=K_w \left\{
\left(E /E_1 \right)^{-\...
...~ ~ ~ E_1 \leq E \leq E_{\rm max} \ \\
\end{array}\right. \ ,
\end{displaymath} (17)

where Kw is related to the fraction of pulsar luminosity that is converted into pairs, $L_{\rm pairs}$ through:
$\displaystyle K_w(t)={L_{\rm pairs}(t) \over E_1^2} \left\{
{1 \over 2-\gamma_1...
...[{\left(E_{\rm max} \over E_1\right)}^{2-\gamma_2}-1\right]
\right\}^{-1} \cdot$     (18)

In the following we shall assume E1 to be constant with time and independent of the wind Lorentz factor: we take $E_1=1.5\ {\rm TeV}$, roughly corresponding to the energy of the particles radiating at the frequency at which the high energy break of the spectrum of the Crab Nebula is observed. As for $E_{\rm max}$, this will depend on both time and the wind Lorentz factor, due to the dependence on the latter of the proton energy and on the dependence on the former of the nebular magnetic field strength. In fact we assume the maximum energy to which the pairs can be accelerated to be given by $E_{\rm max}={\rm min}(E_{\rm p},E_{\rm cut})$, with $E_{\rm cut}$determined by the condition $T_{\rm acc}(E_{\rm cut}) = T_{\rm losses}(E_{\rm cut})$, i.e. the acceleration time cannot exceed the time-scale for energy losses. Assuming $T_{\rm acc} \approx \Omega_{\rm L}^{-1}$ (with $\Omega_{\rm L}$ the particle Larmor frequency) and the losses as due to synchrotron emission in the average nebular magnetic field, one finds $E_{\rm cut}(t) = m_{\rm e} c^2 \sqrt{6 \pi e/(\sigma_{\rm T} B(t))}$.

We mentioned that charged pion decay provides an additional source of pairs in the nebula. Their evolution is completely analogous to that of the wind particles, so that apart from the different initial spectrum, the treatment below applies to both components.

The spectral distribution at a time t, N(E,t), is found as a function of the injection spectrum J=Jw+Je (see Eq. (5)), as:

N(E,t)=\int_E^{E_{\rm max}}J(E_i,t_i)\ {\partial t_i \over \partial E}(E,E_i,t)\
{\rm d}E_i\ ,
\end{displaymath} (19)

where ti is the time at which the particle was born and Ei its initial energy. N(E,t) is readily computed once $\partial t_i / \partial E$ is known as a function of E and t.

The needed relation between E, Ei and ti may be found from the equation describing the single particle energy evolution due to synchrotron and adiabatic losses:

{1 \over E\ t} -{1 \over E_i\ t_i}=\int_{t_i}^t c_1 {B^2(s) \over s} {\rm d}s\ ,
\end{displaymath} (20)

where $c_1=(1/6\pi)\ \sigma_{\rm T}\ c/(m_{\rm e} c^2)^2$ is the constant relating the synchrotron power to magnetic field strength and particle energy. For $\partial t_i / \partial E$ we find:

{\partial t_i \over \partial E}={E_i\ t_i^2(E,t) \over E^2\ t}\
\left[ 1 + {2\ E_i \over E_b(t_i)} \right]^{-1}\ ,
\end{displaymath} (21)

where we have used the definition of the so-called break energy:

E_{\rm b}(t)={2 \over c_1\ t\ B^2(t)}\cdot
\end{displaymath} (22)

For times shorter than $\tau$, when $B(t) \propto t^{-1}$ (Eq. (16)), so that $E_{\rm b}(t)\propto t$, one finds:
$\displaystyle t_i={t \over 2 (1 +E/E_{\rm b}(t))} {E \over E_i} \left\{ 1\!+\!
...E_{\rm b}(t)}
\left(1+ {E \over E_{\rm b}(t) } \right) \right]^{1/2} \right\} ,$     (23)

and the solution of Eq. (19) can be approximated as:

N(E,t) \approx {2\ t \over \gamma_1+\gamma_2}
...E) J(E,t_0) \ , ~ ~ ~ ~ E_{\rm b}(t) \leq E
\end{displaymath} (24)

with $t_0\simeq 0$ for $t<\tau$ and $t_0 \simeq \tau$ afterwards.

The evolution of the wind particle spectrum as a function of time is plotted in Fig. 1 for the case $\Gamma=10^6$, $\gamma _1=1.5$, $\gamma _2=2.2$ and $L_0=3 \times 10^{39}\ {\rm erg}\ {\rm s}^{-1}$, with $L_B=L_{\rm pairs}=L_{\rm p}=(1/3) L$. We plot the full solution of the problem taking into account the real time-dependences arising from the varying pulsar input (Eq. (15) with $\tau=700$ yr) together with the approximate expression in Eq. (24).

\par\includegraphics[width=8.8cm,clip]{2873.f1}\end{figure} Figure 1: The time evolution of the wind particle spectrum in the case $\gamma _1=1.5$, $\gamma _2=2.2$, $L_{\rm pairs}=L_B=10^{39}\ {\rm erg}\ {\rm s}^{-1}$. The different curves refer to different epochs after the supernova explosion: $t=1~ {\rm yr}$ (solid), $t=10~ {\rm yr}$ (dotted), $t=100 ~{\rm yr}$ (dashed), $t=500~ {\rm yr}$ (dot-dashed). The points trace, at each time, the approximation given in the text.

Once N(E,t) and B(t) are known, the volume emissivity $S_\nu$ that enters Eq. (14) is computed from synchrotron theory as:

S^{\rm sync}_\nu(\nu,t)={c_1 B \over 2 c_2} {3 \over 4\ \pi\... c_2\ B}}\ N\left[\sqrt{{\nu \over c_2\ B} }\ , t \right]\ ,
\end{displaymath} (25)

where $c_2=0.3\ e\ c/(2\ \pi\ (m_e\ c^2)^3)$. Knowing $S^{\rm sync}_\nu$one can also compute the ICS emissivity of the plerion, $S^{\rm ICS}_\nu$, following, for instance, Blumenthal & Gould (1970). From the total emissivity $S_\nu=S^{\rm sync}_\nu+S^{\rm ICS}_\nu$ we also obtain the flux $F_\nu$ that appears in Eq. (7):

F_\nu(\nu)={R_N^3 \over 3\ d^2}\ S_\nu(\nu)\ .
\end{displaymath} (26)

3.2 Pion production efficiency

The time-evolution of the magnetic and radiation fields in the nebula as derived in the previous section was used to compute the time-history of pion production efficiency.

In the plots below we show the behaviour of $f^\mu $, $f^\nu $, $f^\pi_{{\rm p} \gamma}$ and $f^\pi_{\rm pp}$ as a function of time after the supernova explosion for various values of the pulsar luminosity L, of the wind Lorentz factor $\Gamma $ and of the remnant expansion velocity v. In all the plots in this section we assumed the pulsar energy input to be evenly shared between magnetic field, protons and pairs: $L_B=L_{\rm p}=L_{\rm pairs}=(1/3)\ L_0$. To compute the radiation field, we assumed the injection spectrum of the pairs as in Eq. (17) with $\gamma _1=1.5$ and $\gamma _2=2.2$. We consider a monoenergetic distribution of protons corresponding to the Lorentz factor of the wind.

In Fig. 2 we plot $f^\mu(E_\pi,t)$ for the value of $E_\pi$ corresponding to the highest energy pions produced for each of the values of $\Gamma $ considered.

\par\includegraphics[width=8.8cm,clip]{2873.f2}\end{figure} Figure 2: The time evolution of $f^\mu $ for different values of the pulsar luminosity, of the remnant expansion velocity and of the wind Lorentz factor. The curves plotted in different panels are computed for different values of the wind Lorentz factor. The different curves in each panel are marked as explained below. Empty triangles: $L[{\rm erg}/{\rm s}]=3 \times 10^{39}$, $v[\rm km~s^{-1}]=10^3$. Empty squares: $L[{\rm erg}/{\rm s}]=3 \times 10^{39}$, $v[\rm km~s^{-1}]=10^4$. Filled triangles: $L[{\rm erg}/{\rm s}] =3 \times 10^{40}$, $v[\rm km~s^{-1}]=10^3$. Filled squares: $L[{\rm erg}/{\rm s}] =3 \times 10^{40}$, $v[\rm km~s^{-1}]=10^4$. The vertical lines represent the time at which the optical depth for neutrinos drops below 1 in the case when $v=10^3\ \rm km~s^{-1}$ (solid) and when $v=10^4\ \rm km~s^{-1}$ (dashed).

\par\includegraphics[width=8.8cm,clip]{2873.f3}\end{figure} Figure 3: The time evolution of $f^\nu $ for the highest energy pions corresponding to different values of the pulsar luminosity, of the remnant expansion velocity and of the wind Lorentz factor. Same notation as in Fig. 2.

It is apparent from the figure that for Lorentz factors up to 105 pions never undergo any severe energy losses, not even at early times, when the magnetic and radiation field in the plerion are most intense. For energetic enough pions, instead, radiation losses at early times do become important, decreasing the decay probability by a factor that can be of the order of 102 for the lowest considered value of the expansion velocity and highest considered pulsar luminosity (filled triangles).

Losses at early times are found to be due mainly to synchrotron emission and their dependence on v and L arises from the relation between the magnetic field intensity and these two quantities (Eq. (16)). It is interesting to note that in all the cases considered, $f^\mu $ becomes 1 at very early times, and actually earlier than the time at which the remnant stops being optically thick to neutrinos, which we mark as a reference time in Fig. 2. This reference time is computed considering the fact that the main source of optical depth for high energy neutrinos are nuclear collisions that give rise to muons: $\nu + X \rightarrow \mu + Y$. The cross section for this process (e.g. Learned & Mannheim 2000) is a function of the neutrino energy alone at the energies we consider. Hence, the time at which $\tau_\nu=R_{\rm N}\ n_{\rm t}\ \sigma_{\nu X}$ drops below unity for the highest energy neutrinos corresponding to each value of the Lorentz factor we consider depends only on the Lorentz factor itself and on the density evolution of the remnant (i.e. on its expansion velocity, once the total mass is fixed and assumed as uniformily distributed). One has: $\tau_\nu \propto \Gamma^\beta/v$, with $\beta=0.5$ if $\Gamma<4 \times 10^5$ and $\beta=0.2$ otherwise.

The same legenda as for Fig. 2 applies to Fig. 3, where the evolution of $f^\nu $ is shown (again in the case when it is minimum, i.e. for the maximum muon energies).

As we could expect, radiation losses for muons are more severe than for pions, and for the most energetic we are considering (lower right panel) $f^\nu $actually becomes 1 later than $f^\mu $.

It is important to notice that as far as only the non-thermal radiation field is taken into account, synchrotron losses are always dominant with respect to losses due to ICS, according to our calculation. The latter can become important, though, if also the thermal radiation field is taken into account, but this effect should be limited to earlier times than that at which $\tau_\nu$ drops below 1 (Beall & Bednarek 2002).

As far as $f^\pi $ is concerned, we computed the efficiency of both photo-meson production and p-p collisions. For the latter process we considered $1\ M_\odot$ of material uniformily distributed within the remnant and used $\mu =1$ in the calculations. The results are shown in Fig. 4.

\par\includegraphics[width=8.8cm,clip]{2873.f4}\end{figure} Figure 4: The time evolution of $f^\pi $ for different values of the pulsar luminosity, of the remnant expansion velocity and of the wind Lorentz factor. The solid curves refer to the efficiency of nuclear collisions while the dashed ones show the efficiency of photo-meson production. Curves referring to different values of the pulsar luminosity and remnant expansion velocity are marked according to the same notation as for Figs. 2 and 3.

We find that radiation losses are never effective in limiting the efficiency of p-p collisions, not even in the very intense initial magnetic field of a bright, slowly expanding plerion. Conversely, for high enough values of the wind Lorentz factor $\Gamma $ and low enough values of the wind expansion velocity, synchrotron losses are the main limit to photo-meson production. However, we can see from Fig. 4 that photomeson production never reaches efficiencies higher than about 30%, at least for the prototypical spectrum that has been assumed here.

It should be noticed, finally, that in the plots shown above for $f^\pi $, $f^\mu $ and $f^\nu $, when computing both ICS losses and photo-meson production we have only included the photon spectrum deriving from synchrotron and ICS emission of the wind particles. In doing that, we have neglected the thermal emission from the gas and dust in the filaments and the radiation derived from pion decay. The consistency of the latter approximation was verified a posteriori: after calculating the radiation flux from the products of pion decay we checked that adding this contribution to the target radiation for ICS losses did not change the results plotted in Figs. 2-4.

3.3 Photon and neutrino fluxes from a "young Crab''

In this section we use the theory introduced above to compute both the electromagnetic and neutrino emission from a young plerion with Crab-like parameters, as a function of time after the supernova explosion. We consider a pulsar whose initial spin-down power is $L_0= 3 \times 10^{39}$ erg/s. We assume, as in Sect. 3.2, that L0 is evenly converted into electron, proton and magnetic field energy: $L_B=L_{\rm p}=L_{\rm pairs}=(1/3)\ L_0$. The distance of the pulsar and the plerion expansion velocity are taken to be the same as for the Crab, and 1 $M_\odot$ of thermal material uniformily distributed in the plerion volume is considered as the target for p-p collisions.

The time evolution of the pion production efficiency in such a plerion was computed in Sect. 3.2 (curves marked with empty squares in Figs. 2-4). Taking the results from that section we estimate the global photon spectrum of the plerion. The high energy part is entirely due to $\pi ^0$ decay, while at low energies we estimate the synchrotron emission of both wind and secondary pairs together with the ICS emission of the former. The overall spectrum at different times during the plerion evolution is reported in Fig. 5.

\par\includegraphics[width=8.8cm,clip]{2873.f5}\end{figure} Figure 5: Electromagnetic spectrum of a "young Crab'' at different times after the supernova explosion. The thick curves correspond to the synchrotron (at low energies) and ICS (at higher energies) emission of the wind pairs, while the lighter curves refer to the emission of the pion decay products. The different panels correspond to different values of the wind Lorentz factor $\Gamma $ and the different line types refer to different times.

It is apparent from the figure that, for the set of parameters used, the radiation due to the secondary pairs could only be detectable over the emission due to the wind particles, at early times (up to few tens of years after the supernova explosion) and at photon energies in the range 100 MeV-1 TeV.

We notice that $\gamma$-rays from $\pi ^0$ decay should be detected even by present TeV detectors during the first few tens of years after the supernova explosion, if a similar plerion was born within the galaxy. The $\pi ^0$ decay $\gamma$-ray flux always exceeds the ICS emission at high enough energies for times up to several tens of years at least.

\par\includegraphics[width=8.8cm,clip]{2873.f6}\end{figure} Figure 6: Neutrino spectrum of a "young Crab''. The notation is the same as for Fig. 5. The solid thick line is the detection threshold for ICE-cube (Hill 2001), while the dot-dashed line is the atmospheric neutrino background in a $2^\circ \times 2^\circ $ bin.

In Fig. 6 we show the expected neutrino flux from the same plerion. We find that this is generally well above the atmospheric background of a detector with an angular resolution of $2^\circ \times 2^\circ $. Moreover a galactic object described by our set of parameters could give a flux of neutrinos already detectable by telescopes like AMANDA II for the first ten years. It is clear from the figure that the probability of detecting such objects will increase considerably with the upcoming ${\rm km}^2$ neutrino telescopes, since their lower threshold allows neutrino detection even from older plerions, if the pulsar wind contains a non-negligible fraction of protons.

4 Constraining the parameters of the Crab pulsar wind

In this section, we use the Crab Nebula, the prototype plerion, to discuss how the information coming from high energy photon and neutrino detection could be used to constrain the physics of plerions.

The parameters playing a role in determining the flux of high energy photons and neutrinos are all quite well known for the Crab Nebula, except for the parameters that enter in the determination of $\mu $, the wind Lorentz factor $\Gamma $, and the fraction of the pulsar spin-down power that goes into acceleration of protons. The relevant parameters for our purposes are, from the equations in Sect. 2: the plerion age $t =950 ~{\rm yr}$, its distance $d \simeq 2\ {\rm kpc}$, its expansion velocity, $V_{\rm exp} = 1.5 \times 10^3\ {\rm km}\ {\rm s}^{-1}$, the average magnetic field intensity $B \simeq 3 \times 10^{-4}$ G, and finally the amount of mass contained in the filamentary shell. For the latter, the most recent estimates give $M_{\rm N}=(4.6\pm 1.8)~ M_\odot$(Fesen et al. 1997). We shall consider a value $M_{\rm N}=1~ M_\odot$and values of $\mu $ in the range $1<\mu<20$, according to the upper limit on $n_{\rm t}$ derived by Atoyan & Aharonian (1996), based on the compatibility between the bremsstrahlung $\gamma$-rays and the observed flux from the Crab Nebula above 1 GeV. We notice that if the mass estimated by Fesen et al. (1997) was uniformally distributed in the nebula, this would correspond to $\mu =5$.

Inserting the numbers just quoted into the equations in Sect. 2 we can estimate all the relevant time-scales. One can readily see that what is important for computing $f^\pi $ is the comparison between $t_{\rm esc}$ and $t_{\rm exp}$ (whose relative importance depends on $\Gamma $) and $t_{\rm pp}$, since the other three time-scales are much longer than these. From Eqs. (1), (10)-(12), we have:

                                    $\displaystyle f^\pi$ = $\displaystyle {\rm min} \left[1, {t_{\rm exp} \over t_{{\rm pp}}} {\rm min}
\left(1,{t_{\rm esc} \over t_{\rm exp}} \right) \right]$  
  = $\displaystyle {\rm min} \left[ 1, 2.8 \times 10^{-5} \mu\ {\rm min}
\left(1, {2.2 \over \Gamma_6} \right) \right],$ (27)

where $\Gamma_6$ is the wind Lorentz factor in units of 106.

For the same values of the parameters the muon losses due to synchroton emission or ICS (Eqs. (13) and (14)) are completely negligible on the time-scale of their decay, even for the highest energy we consider. The same is true for pions, whose proper decay times are about a hundred times shorter.

In Figs. 7 and 8 we plot the expected neutrino and photon fluxes from the Crab Nebula as a function of energy, for different assumptions on the wind Lorentz factor and different values of the parameter $\mu $.

\par\includegraphics[width=8.8cm,clip]{2873.f7}\end{figure} Figure 7: The flux of neutrinos expected from the Crab Nebula as a function of energy, for different values of the parameter $\mu $ and for different values of the wind Lorentz factor $\Gamma $. The different curves in each panel correspond to $\mu =1$ (solid), $\mu =5$ (dotted), $\mu =10$ (dashed) and $\mu =20$ (long-dashed). The considered values of $\Gamma $ range from 104 to 107, as specified in each panel. The solid thick line is the detection threshold for IceCube (Hill 2001), while the dot-dashed line is the atmospheric neutrino background in a $2^\circ \times 2^\circ $ bin. We assume 60% of the pulsar wind energy to be carried by protons.

We have assumed that protons are the main energy carriers in the wind so that we put in them about 60% of the pulsar spin-down luminosity $L_{\rm p} = 0.6\ L_{\rm Tot} \simeq 3 \times 10^{38}\ {\rm erg}\ {\rm s}^{-1}$, being $L_{\rm Tot}=5 \times 10^{38}\ {\rm erg}\ {\rm s}^{-1}$. Our choice of the energy fraction into baryons comes from the fact that this is the percentage of the pulsar power that is required to be in baryons if the luminosity contrast of the Crab Nebula wisps is to be explained through magnetic field compression due to the ion gyration in the shocked background pair plasma (Spitkovsky & Arons 1999).

The neutrino and photon fluxes were computed from Eqs. (4) and (6) respectively, with the pion differential spectrum per unit time computed in the scaling approximation for the p-p scattering cross-section (Blasi 2001). In the scaling regime a monoenergetic proton distribution with energy $E_{\rm p}$ leads to the injection of pions described as:

 \begin{displaymath}{{\rm d}N_\pi \over {\rm d}E_\pi {\rm d}t}=K_\pi {E_\pi}^{-1}\ g_\pi(E_\pi/E_{\rm p})\ ,
\end{displaymath} (28)

where the function $g_\pi$ is

 \begin{displaymath}g_\pi(x)=(1-x)^{3.5}+{\rm e}^{-18x}/1.34\ ,
\end{displaymath} (29)

and $K_\pi$ is found from the condition

 \begin{displaymath}f_\pi L_{\rm p} = \int_{m_\pi c^2}^{E_{\rm p}} {{\rm d}N_\pi \over {\rm d}E_\pi {\rm d}t} E_\pi {\rm d}E_\pi\ ,
\end{displaymath} (30)

leading to:

 \begin{displaymath}K_\pi=f^\pi {L_{\rm p} \over E_{\rm p}} \left( \int_0^1 g_\pi(x) {\rm d}x \right)^{-1}\ .
\end{displaymath} (31)

It is apparent from Fig. 7 that, for the generally favoured value of $\Gamma \simeq 10^6$ (Kennel & Coroniti 1984a and 1984b; Gallant & Arons 1994), and for the value of $L_{\rm p} \simeq 0.6 L_{\rm tot}$ favoured by Gallant & Arons (1994) and Spitkovsky & Arons (1999), the flux of neutrinos from the Crab Nebula would be well above the atmospheric neutrino background already for $\mu =1$. The chances of detecting multi-TeV neutrinos from the Crab Nebula with IceCube look promising, based on this figure.
\par\includegraphics[width=8.8cm,clip]{2873.f8}\end{figure} Figure 8: The flux of high energy photons expected from the Crab Nebula, for different values of the parameter $\mu $ and of the wind Lorentz factor $\Gamma $. The notation is the same as in Fig. 7. The expected fluxes are compared with the available data and upper limits (de Jager et al. 1996; Aharonian et al. 2000 and references therein). Also shown as a thick curve in each panel is the ICS flux computed based on this model.

However, part of the section of the parameter space that would lead to a detectable flux of neutrinos can already be excluded based on the available high energy $\gamma$-ray data. In Fig. 8 we plot the flux of high energy photons expected from $\pi ^0$ decay. It is clear that in some cases the computed fluxes at multi-TeV energies are above the measurements and upper limits. The latter, on the other hand, can be comfortably explained, within the error bars, as due to ICS (de Jager et al. 1996; Atoyan & Aharonian 1996; Aharonian et al. 2000). Even within the framework of our crude model, which not only treats the plerion as homogeneous but also neglects the contribution to the target photon field of the ambient thermal emission, the computed ICS flux is only slightly below the TeV data points.

In the following, we review our expectation for the flux of neutrinos from the Crab Nebula in the light of the constraints arising from $\gamma$-ray observations.

Before discussing these limits, however, we should also consider the other possibility we mentioned to put constraints on the three unknowns $\mu $, $L_{\rm p}$ and $\Gamma $. This is offered, at least in principle, by the signatures that may be left by the secondary pairs in the synchrotron radiation spectrum of the nebula. These pairs would be injected in the Nebula according to Eq. (5), but the injection rate (described by Eq. (28) and depending on $f^\pi $ and $L_{\rm p}$ through Eq. (31)) has to be considered now as a function of time. Since we are mainly interested in low energy pairs, with long synchrotron life-times, we need to take into account the whole history of pion production in the nebula. This is done as discussed in Sects. 3.1 and 3.2.

\par\includegraphics[width=8.8cm,clip]{2873.f9}\end{figure} Figure 9: Consequences of pion production on the electromagnetic spectrum of the Crab Nebula. The thick curve and the points represent the observed emission, while the thin curves are fluxes from our model. At low frequencies the model flux is due to synchrotron emission associated to secondary electrons and positrons, while the contribution at high frequencies is the same shown in Fig. 8 and comes from $\pi ^0$ decay. The different panels refer to the emission computed for different values of the wind Lorentz factor while the different line-types are associated to different values of the parameter $\mu $ as specified in the each panel.

We show in Fig. 9 the comparison between the Crab Nebula observed radiation spectrum and the synchrotron flux from secondary pairs computed in the case when $L_{\rm p}(t)=0.6\ L_{\rm Tot} (t)$, with $L_{\rm Tot} (t)$ described by Eq. (15), where $\tau=730\ {\rm yr}$ and $\alpha=2.3$ (Groth 1975). In the same figure, we also show the contribution arising from $\pi ^0$ decay, in the $\gamma$-ray part of the spectrum. It is apparent that the emission produced by the secondary pairs, according to our model, is a completely negligible fraction of the present day synchrotron emission of the Crab Nebula and cannot provide any additional constraint on the parameters we are interested in.

\par\includegraphics[width=6.9cm,clip]{2873.f10}\end{figure} Figure 10: The allowed parameter space for $\Gamma $ and $L_{\rm p}/L_{\rm Tot}$ as constrained through the comparison with observations of the photon fluxes from pion decay computed within the framework of our model. The allowed regions for each of the considered value of $\mu $ are the ones lying under the curve of corresponding type: $\mu =5 \rightarrow $ solid, $\mu =10 \rightarrow $ dashed, $\mu =20 \rightarrow $ dot-dashed. No constraints can be put on the wind parameters for $\mu =1$.

In Fig. 10 we use the information that can be extracted from high energy data (Fig. 8) to restrict the allowed parameter space for the values of $L_{\rm p}/L_{\rm Tot}$, $\Gamma $ and $\mu $. The allowed regions, for each choice of $\mu $, are those lying below the corresponding curve, as explained in the figure caption.

It is interesting to notice that, if $\mu =10$ or larger, the possibility that 60% of the pulsar wind energy might be in protons with a Lorentz factor $\Gamma=10^6$ would be excluded by the data.

We now use the constraints in Fig. 10 to estimate the maximum number of neutrino events that one could expect to detect from the Crab Nebula in 1 yr of observation with a detector of area $1\ {\rm km}^2$, such as IceCube will be. The result of this calculation is shown in Fig. 11 where we plot the maximum number of muons per year as a function of the wind Lorentz factor. $N_\mu$ is computed taking into account the probability for a neutrino to interact with the matter in the detector and create a muon. For each value of $\Gamma $ we consider the maximum allowed value of $f^\pi L_{\rm p}$ with the requirement that $\mu \le 20$ and $L_{\rm p} \le 0.9 L_{\rm Tot}$. This determines $J_\nu(\epsilon_\nu)$. $N_\mu$ is then found as:

 \begin{displaymath}N_\mu={1~ {\rm yr} \times 1~ {\rm km}^2 \over 4\ \pi\ {\rm d}...
...u(\epsilon_\nu) P_{\nu \mu}(\epsilon_\nu) {\rm d} \epsilon_\nu
\end{displaymath} (32)

where the maximum neutrino energy is $\epsilon_{\rm Max}=m_{\rm p} c^2 \Gamma/4$ and $P_{\nu \mu}=1.3 \times 10^{-6} \epsilon_{\nu,{\rm TeV}}^\beta$ with $\beta=1$ if $1\ {\rm TeV}<\epsilon_\nu<100\ {\rm TeV}$, and $\beta=0.5$ otherwise (Gaisser et al. 1995).

\par\includegraphics[width=8.8cm,clip]{2873.f11}\end{figure} Figure 11: The maximum number of muons per year that would be produced by neutrinos coming from the Crab Nebula in a ${\rm km}^2$ detector. Also shown, as a dot-dashed line, are the background counts, assuming an angular resolution of $(1^\circ )^2$, such as IceCube should have.

The background counts we show in Fig. 11 are those computed for IceCube by Alvarez-Muñiz & Halzen (2001) at the position of the Crab Nebula. The maximum number of neutrino events allowed by the $\gamma$-ray data is above the background for basically the entire range of values of $\Gamma $ we considered. The effect of neutrino oscillations would be to reduce these fluxes. For the case of maximal mixing, equal abundances of neutrinos of the three flavors are obtained, so that the flux of muon neutrinos is about half of that calculated above.

5 Conclusions

In this paper we have analysed the observational consequences that the presence in pulsar winds of a proton component, energetically dominant, might have.

We have found that when the plerion is very young (up to a few tens of years after the supernova explosion) the main loss mechanism for these protons is likely to be pion production due to nuclear collisions. A large flux of neutrinos, $\gamma$-rays and secondary pairs can be expected during these early stages. The synchrotron emission from secondaries is likely to be completely negligible, while the TeV photon and neutrino emission from a galactic object could be detectable even by present- day telescopes during the first few tens of years after the supernova explosion. Detection of neutrinos from older plerions should await the upcoming ${\rm km}^2$ detectors.

We have also estimated the present-day TeV photon and neutrino emission from the Crab Nebula. The existing TeV observations already exclude part of the parameter space, so that not all combinations of $\Gamma $, $L_{\rm p}$ and $\mu $ are allowed.

Finally we have computed the maximum number of neutrino events allowed by the $\gamma$-ray data according to our modeling and found it to be larger than the atmospheric neutrino background (as already suggested by Arons, 1998). Such a flux is likely to be detected by a ${\rm km}^2$ detector such as IceCube. Neutrino detection can be therefore considered as a promising means to help answer the question of how large a fraction of the pulsar spin-down energy goes into relativistic baryons.

We thank F. Pacini and M. Salvati for fruitful discussions. We are very grateful to the referee, Yves Gallant, for valuable suggestions. This work was partly supported by MIUR under grant Cofin 2001-02-10.


Copyright ESO 2003