Parameterization of γ, e± and Neutrino Spectra Produced by p −p Interaction in Astronomical Environment Tuneyoshi Kamae1, Niklas Karlsson2, Tsunefumi Mizuno3, Toshinori Abe4, Tatsumi Koi Stanford Linear Accelerator Center, Menlo Park, CA 94025 [email protected] 7 0 0 ABSTRACT 2 n a ± We present the yield and spectra of stable secondary particles (γ, e , ν , ν¯ , J e e 6 νµ, andν¯µ)ofp−p interactioninparameterized formulaetofacilitatecalculations 1 involving them in astronomical environments. The formulae are derived from the 3 up-to-datep−pinteractionmodelbyKamae et al.(2005),whichincorporatesthe v logarithmically rising inelastic cross section, the diffraction dissociation process, 1 8 and the Feynman scaling violation. To improve fidelity to experimental data 5 5 in lower energies, two baryon resonance contributions have been added: one 0 representing ∆(1232)andtheotherrepresenting multipleresonances around1600 6 0 MeV/c2. The parametrized formulae predict that all secondary particle spectra / h beharder by about 0.05in power-law indices thanthatof theincident protonand p their inclusive cross-sections be larger than those predicted by p−p interaction - o models based on the Feynman scaling. r t s a Subject headings: cosmic rays — galaxies: jets — gamma-rays: theory — ISM: : v general — neutrinos — supernovae: general i X r a 1. Introduction Gamma-ray emission due to neutral pions produced by proton-proton interaction has long been predicted from the Galactic ridge, supernova remnants (SNRs), active galactic nu- cleus(AGN)jets,andotherastronomicalsites(Hayakawa1969;Stecker1971;Murthy & Wolfendale 1Also with Kavli Institute for Particle Astrophysics and Cosmology, Stanford University, Menlo Park, CA 94025 2Visiting scientist from Royal Institute of Technology, SE-10044 Stockholm, Sweden 3Present address: Department of Physics, Hiroshima University, Higashi-Hiroshima, Japan 739-8511 4 Present address: Department of Physics, University of Tokyo, Tokyo, Japan 113-0033 – 2 – 1986; Sch¨onfelder 2001; Schlickeiser 2002; Aharonian 2004; Aharonian et al. 2004). High en- ergyneutrinosproducedbyp−pinteractioninAGNjetswillsoonbedetectedwithlarge-scale neutrino detectors that are under construction (Halzen 2005). Spectra of these gamma-rays, neutrinos, and other secondaries depend heavily on the incident proton spectrum, which is unknown and needs to be derived, in almost all cases, from the observed spectra them- selves. Such analyses often involve iterative calculations with many trial proton spectra. The parameterized model presented here is aimed to improve accuracy of such calculations. Among the secondaries of p −p interaction in astronomical environment, gamma-rays have been best studied. The p-p interaction is one of the two dominant gamma-ray emission mechanisms in the sub-GeV to multi-TeV range, the other being Compton up-scattering of low energy photons by high energy electrons. Gamma-rays in this energy range have been detectedfrompulsars, theGalacticRidge, SNRs, blazarsandothersourcecategories(Stecker 1971; Murthy & Wolfendale 1986; Ong 1998; Sch¨onfelder 2001; Weekes 2003; Aharonian 2004). High energy gamma-rays from AGN jets are interpreted as mostly due to the inverse Compton up-scattering of low-energy photons by multi-TeV electrons. Observed radio and X-ray spectra match those of synchrotron radiation by these electrons. Synchronicity in variability between the observed X-ray and gamma-ray fluxes has given strong support for the inverse Compton up-scattering scenario (Ong 1998; Scho¨nfelder 2001; Schlickeiser 2002; Aharonian 2004). For some AGN jets, the above scenario does not work well, and p − p interaction has been proposed as an alternative mechanism (Mu¨cke & Protheroe 2001; Mu¨cke et al. 2003; B¨ottcher & Reimer 2004). High energy gamma-rays detected by COS-B and EGRET from the Galactic Ridge, on the other hand, are interpreted as predominantly due to neutral pions produced by interac- tion of protons andnuclei with the interstellar matter (ISM; Stecker 1973; Strong et al. 1978, 1982,2000;Stephens & Badhwar1981;Mayer-Hasselwander et al.1982;Bloemen et al.1984; Bloemen 1985; Dermer 1986a; Stecker 1990; Hunter et al. 1997; Mori 1997; Stanev 2004). The measured gamma-ray flux and spectral shape (Hunter et al. 1997) have been viewed as the key attestation to this interpretation. It is also known that inverse Compton scattering contributes significantly to the Galactic Ridge gamma-ray emission (Murthy & Wolfendale 1986; Strong et al. 2000; Sch¨onfelder 2001). In the past two years, several SNRs, including RX J1713-3946 and RX J0852-4622, ◦ have been imaged in the TeV band with an angular resolution around 0.1 by H.E.S.S. (Aharonian et al. 2004b,c; Aharonian 2005). A smooth featureless spectrum, suggestive of synchrotron radiation by multi-TeV electrons, has been detected in the X-ray band from RX J1713-3946 (Koyama et al. 1997; Slane et al. 1999; Uchiyama et al. 2003) and RX J0852- – 3 – 4622 (Tsunemi et al. 2000; Iyudin et al. 2005). The measured TeV gamma-ray fluxes and spectra, however, do not agree well with those predicted by the inverse Compton scenario (see, e.g., the analysis in Uchiyama et al. (2003)). Several authors have proposed that the TeV gamma-rays are possibly due to interaction of accelerated protons with the ISM (Berezhko & Volk 2000; Enomoto et al. 2002; Aharonian 2004; Katagiri et al. 2005). Higher precision data are expected, in the GeV range, from GLAST Large Area Tele- scope (GLAST-LAT 2005)1 and, in the TeV range, from the upgraded Air Cherenkov Tele- scopes (Aharonian et al. 2004): they will soon test applicability of the inverse Compton up- scattering and the proton interaction with ISM for various astronomical gamma-ray sources. In many objects, secondary electrons and positrons may produce fluxes of hard X-rays and low-energy gamma-rays detectable with high-sensitivity instruments aboard Integral, Swift, Suzaku, and NuSTAR.2 The formulae given here will give fluxes and spectra of these sec- ondary particles for arbitrary incident proton spectrum. This work is an extension of that by Kamae et al. (2005), where up-to-date knowl- edge of the π0 yield in the p − p inelastic interaction has been used to predict the Galac- tic diffuse gamma-ray emission. The authors have found that past calculations (Stecker 1970, 1973, 1990; Strong et al. 1978; Stephens & Badhwar 1981; Dermer 1986a,b; Mori 1997; Strong et al. 2000) had left out the diffractive interaction and the Feynman scaling violation in the non-diffractive inelastic interaction. Another important finding by them is that most previous calculations have assumed an energy-independent p − p inelastic cross-section of about 24 mb for T ≫ 10 GeV, whereas recent experimental data have established a log- p arithmic increase with the incident proton energy. Updating these shortfalls has changed the prediction on the gamma-ray spectrum in the GeV band significantly: the gamma-ray power-law index is harder than that of the incident proton; and the GeV−TeV gamma-ray flux is significantly larger than that predicted on the constant cross-section and Feynman scaling (Kamae et al. 2005). The model by Kamae et al. (2005) will hereafter be referred to as model A. Model A does not model p−p interaction accurately near the pion production thresh- old. To improve prediction of gamma-rays, electrons, and positrons produced near the pion production threshold, two baryon resonance excitation components have bee added to model A: ∆(1232), representing the ∆ resonance, and res(1600), representing resonances around 1GLAST Large Area Telescope, http://www-glast.stanford.edu. 2See INTEGRAL Web site, http://sci.esa.int/esaMI/Integral, the NuSTAR Web site, http://www.nustar.caltech.edu, the Suzaku Web site http://www.isas.jaxa.jp/e/enterp/missions/astro-e2, and the Swift Web site, http://swift.gsfc.nasa.gov/docs/swift. – 4 – 1600 MeV/c2. We note here that ∆(1232) is the most prominent and lightest Baryon res- onance excited in p − p interaction. It has a mass of 1.232 GeV/c2 and a width of about 0.12 GeV/c2, and decays to a nucleon (proton or neutron) and a pion (π+,0,−). The other resonance, res(1600), is assumed to decay to a nucleon and two pions. Introduction of these contributions have necessitated adjustment of the model A at lower energies as described below. The readjusted model will be referred to as the “readjusted model A”. The parameterized model presented here exhibits all features of model A at higher energies (proton kinetic energy, T > 3 GeV) and reproduces experimental data down to the p pion production threshold. The inclusive gamma-ray and neutrino cross-section formulae can be used to predict their yields and spectra for a wide range of incident proton spectrum. Formulae for electrons and positrons predict sub-TeV to multi-TeV secondary electrons and positrons supplied by p − p interaction. We note that space-borne experiments such as PAMELA3 will soon measure the electron and positron spectra in the sub-TeV energy range, wherethesecondariesofp−pinteractionmaybecomecomparabletotheprimarycomponents (Mu¨ller 2001; Stephens 2001; DuVernois et al. 2001). Due to paucity of experimental data and widely accepted modeling, we have not param- eterized inclusive secondary cross-sections for α −p nor p−He nor α−He interactions. We note that α-particles are known to make up about 7% by number of cosmic-rays observed near the Earth (Schlickeiser 2002) and that He to make up ∼ 10% by number of interstellar gas. The total non p − p contribution is comparable to that of p − p contribution. The α-particle and He nucleus can be regarded, to a good approximation, as four independent protons beyond the resonance region (T > 3 GeV): the error introduced is expected to be p less than 10% for high energy light secondary particles (Kamae et al. 2005). Inclusion of α-particle as projectile and He nuclei as targets will change the positron- electron ratio significantly (about 10−15%) as discussed below. Fermi motion of nucleons and multiple nucleonic interactions in the nucleus are known to significantly affect pion production near the threshold and in the resonance region (T < 3 GeV; Crawford et al. p 1980; M˚artensson et al. 2000); we acknowledge need for separate treatment of p−He, α−p, and α−He interactions in the future. 3 See http://wizard.roma2.infn.it/pamela. – 5 – 2. Monte Carlo Event Generation ± The parameterization of the inclusive cross sections for γ, e , ν, and ν¯ has been carried out, separately, for non-diffractive, diffractive, and resonance-excitation processes, in three steps: First, the secondary particle spectra have been extracted out of events generated for mono-energetic protons (0.488 GeV < T < 512 TeV) based on the readjusted model A. We p then fit these spectra with a common parameterized function, separately for non-diffractive, diffractive and resonance-excitation processes. In the third step, the parameters determined for mono-energetic protons are fitted as functions of proton energy, again separately for the three processes. The above procedure has been repeated for all secondary particle types. The functional formulae often introduce tails extending beyond the energy-momentum conservation limits, which may produce artifacts when wide range spectral energy density (E2dflux(γ)/dE)isplotted. Toeliminatesuchartifacts, weintroduceanothersetoffunctions to impose the kinematic limits. Several simulation programs have been used in model A (Kamae et al. 2005): for the highenergy non-diffractive process (T > 52.6 GeV),Pythia 6.2 (Sj¨ostrand et al. 2001)with p the multi-parton-level scaling violation option (Sj¨ostrand & Skands 2004);4 for the lower en- ergy non-diffractive process, the parameterized model by (Blattnig et al. 2000); for the diffractive part, the program by T. Kamae (2004, personal communication)5. In the read- justed model A, two programs to simulate two resonance-excitation components have also been added. Modeling of the two resonance components will be explained below. 3. Non-Diffractive, Diffractive, and Resonance-Excitation Cross-Sections Experimental data on p − p cross-sections are archived for a broad range of the inci- dent proton energy and various final states. The total and elastic cross-sections have been compiled from those by Hagiwara et al. (2002), as shown in Figure 1. The two thin curves running through experimental data points in the figure are our eye-ball fits to the total and elastic cross sections. We then define the “empirical” inelastic cross section as the difference of the two curves to which the sum of non-diffractive, diffractive, ∆(1232)-excitation, and res(1600)-excitation components are constrained. Typical errors in the empirical inelastic cross-section are 20% for T < 3 GeV and 10% for T > 3 GeV. p p 4See http://cepa.fnal.gov/CPD/MCTuning1and http://www.phys.ufl.edu/∼rfield/cdf. 5 Diffractive process has been included in Pythia after the work began. – 6 – The four component cross sections of the readjusted model A and their sum are shown inFigures 1and 2. The empirical inelastic cross section isshown by aseries ofsmall circles in Figure2. Thecomponentcrosssectionstakeformulaegiveninequations(1)(non-diffractive), (2) (diffractive), (3) [∆(1232)] and (4) [res(1600)]. These are also shown in Figures 1 and 2. We note that there is no clear experimental method separating the four components, especially at lower energies (T < 20 GeV). This ambiguity does not significantly affect the p secondary particle fluxes, as long as the sum agrees with the total inelastic cross section and the total secondary inclusive cross sections agree with the corresponding experimental data. The secondary particle spectra for a mono-energetic proton are normalized to the com- ponent cross sections given in equations (1), (2), (3) and (4) at the corresponding proton energy. We note that the non-diffractive component for T < 52.6 GeV is based on the p formula by Blattnig et al. (2000), which is normalized to their π0 inclusive cross section formula not to the total inelastic cross section. In the readjusted model A, this compo- nent cross section is defined by equation (1) and the π0 inclusive cross section formula of Blattnig et al. (2000) has been redefined so that the sum of the four components reproduce the experimental π0 inclusive cross section. The positive and negative pion inclusive cross sections are also redefined as products of our π0 inclusive cross section and the ratio of the π+,− and π0 inclusive cross sections given in Blattnig et al. (2000). – 7 – 200 100 50 20 b℄ m [ 10 p p (cid:27) 5 2 1 0.5 1 1 2 3 4 5 10 10 10 10 10 Pp [GeV/ ℄ Fig. 1.— Experimental p − p cross sections, as a function of proton momentum, and that of readjusted model A: experimental total (squares), experimental elastic (triangles), total inelastic (thick solid line), non-diffractive (dashed line), diffractive process (dot-dashed line), ∆(1232) (dotted line), and res(1600) (thin solid line). The total inelastic is the sum of the four components. The thin solid and dot-dot-dashed lines running through the two experimental data sets are eye-ball fits to the total and elastic cross sections, respectively. – 8 – 40 30 b℄ m [ el 20 n p;i p (cid:27) 10 0 1 2 5 10 Pp [GeV/ ℄ Fig. 2.— Experimental p − p cross sections, as a function of proton momentum, and that of readjusted model A for T < 10 GeV. Small circles represent the empirical inelastic cross p section described in the text. Lines are the same as in Figure 1. – 9 – 0 P < 1 GeV/c, p  0.57(x/a )1.2(a +a x2 +a x3  0 2 3 4    +a5exp(−a6(x+a7)2)) 1 ≤ Pp ≤ 1.3 GeV/c,  σpp (x)[mb] =  (b |a −x|+b |a −x|)/(a −a ) 1.3 ≤ P ≤ 2.4GeV/c, (1) NonDiff  0 1 1 0 1 0 p a +a x2 +a x3 2 3 4   +a exp(−a (x+a )2) 2.4 ≤ P ≤ 10 GeV/c,  5 6 7 p   c +c x+c x2 P > 10 GeV/c,  0 1 2 p  0 P < 2.25 GeV/c, p  (x−d )/d  p 0 1  ×(d2 +d3log10(d4(x−0.25))  σpp (x)[mb] =  +d x2 −d x3) 2.25 ≤ P ≤ 3.2 GeV/c, (2) Diff  5 6 p d +d log (d (x−0.25)) 2 3 10 4   +d x2 −d x3 3.2 ≤ P ≤ 100 GeV/c,  5 6 p    e +e x P > 100 GeV/c,  0 1 p  0 E < 1.4 GeV, p   f0Ep10 1.4 ≤ Ep ≤ 1.6 GeV,  σpp (x)[mb] =  f exp(−f (E −f )2) 1.6 ≤ E ≤ 1.8 GeV, (3) ∆1232  1 2 p 3 p f E−10 1.8 ≤ E ≤ 10 GeV, 4 p p   0 E > 10 GeV,  p   0 E < 1.6 GeV, p   g0Ep14 1.6 ≤ Ep ≤ 1.9 GeV,  σpp (x)[mb] =  g exp(−g (E −g )2) 1.9 ≤ E ≤ 2.3 GeV, (4) Res(1600)  1 2 p 3 p g E−6 2.3 ≤ E ≤ 20 GeV, 4 p p   0 E > 20 GeV,  p   where x = log (P [GeV/c]) and E is the proton energy in GeV. 10 p p 3.1. Introduction of Resonance-Excitation Processes to Model A One or both of the projectile and target protons can be excited to baryon resonances in the p−p interaction. Here we use “baryon resonances” to represent both nucleon resonances – 10 – (iso-spin=1/2) and ∆ resonances (iso-spin=3/2). These excitations enhance the pion pro- duction (and hence secondary particle production) near the inelastic threshold. The most prominent resonance among them is ∆(1232), which has a mass of 1232 MeV/c2 and decays predominantly (> 99%) to a nucleon and a pion (Hagiwara et al. 2002). Stecker (1970) proposed a cosmic gamma-raymodel inwhich neutral pions areproduced only through the ∆(1232) excitation for T ≤ 2.2 GeV. The resonance is assumed to move p only in the direction of the incident proton. At higher energies, another process, the fireball process, sets in and produces pions with limited transverse momenta. Dermer (1986a) compared predictions of models on π0 kinetic energy distribution in the proton-proton center-of-mass (CM) system with experiments and noted that the model by Stecker (1970) reproduces data better than the scaling model by Stephens & Badhwar (1981) for T < 3 GeV. He proposed a cosmic gamma-ray production model that covers a p wider energy range by connecting the two models in the energy range T = 3−7 GeV. p Model A by Kamae et al. (2005) has been constructed primarily for the p−p inelastic interaction T ≫ 1 GeV and has left room for improvement for T < 3 GeV. The diffrac- p p tion dissociation component of model A has a resonance-excitation feature similar to that implemented in Stecker (1971) for T > 3 GeV where either or both protons can be excited p to nucleon resonances (iso-spin=1/2 and mass around 1600 MeV/c2) along the direction of the incident and/or target protons. What has not been implemented in model A is the en- hancement by baryon resonances in the inclusive pion production cross sections below T < 3 p GeV. We note here that the the models by Stecker (1970) and by Dermer (1986a, see also Dermer 1986b) used experimental data on the inclusive π0 yield (and that of charged pions) to guide their modeling, but not the total inelastic cross section. Model A by Kamae et al. (2005), on the other hand, has simulated all particles in each event (referred to as the “ex- clusive” particle distribution) for all component cross sections. One exception is simulation of the low-energy non-diffractive process (T < 52.6 GeV) by Blattnig et al. (2000). The p inclusive π0 (or gamma-ray) yield is obtained by collecting π0 (or gamma-rays) in simulated exclusive events. When readjusting model A by adding the resonance-excitation feature similar to that by Stecker (1970), overall coherence to model A has been kept. We adjusted the ∆(1232) excitation cross section to reproduce the total inelastic cross section given in Figure 2 and fixed the average pion multiplicities for + : 0 : − to those expected by the one-pion-exchange hypothesis, 0.73 : 0.27 : 0.0. As higher-mass resonances begin to con- tribute, the average pion multiplicity is expected to increase. To reproduce the experimental π0 inclusive cross section and total inelastic cross section for T < 3 GeV, we introduced a p second resonance, res(1600). This resonance does not correspond to any specific resonance