GALACTIC COSMIC RAY ORIGINS AND OB ASSOCIATIONS: EVIDENCE FROM SuperTIGER OBSERVATIONS OF ELEMENTS ₂₆Fe THROUGH ₄₀Zr

The SuperTIGER (Trans-Iron Galactic Element Recorder) instrument (Binns et al. 2014) was flown on a NASA long-duration balloon flight over Antarctica for 55 days in the 2012–2013 austral summer at altitudes from about 36.6 to 39.6 km and a mean atmospheric overburden of 4.4 g cm⁻². The instrument measured the elemental abundances of Galactic cosmic ray (GCR) nuclei with 10 ≤ Z ≤ 40 above ∼700 GeV nucleon⁻¹ (at the top of the atmosphere).

In this paper we present analysis of measurements of the elemental composition of "ultra-heavy" GCRs with atomic number 26 ≤ Z ≤ 40. These measurements are the first in which each element in the 30 ≤ Z ≤ 40 charge range has been measured with single-element resolution and good statistics.

In recent years, an explanation of GCR origins has emerged based on measurements from the Cosmic Ray Isotope Spectrometer (CRIS) (Stone et al. 1998a) onboard the NASA Advanced Composition Explorer (ACE) satellite (Stone et al. 1998b) and from the TIGER balloon-borne instrument (Link 2003; Rauch 2008; Rauch et al. 2009).

In this explanation, the GCR source material is thought to be a mixture of material from massive stars (supernova ejecta and stellar wind outflow), primarily within OB associations, and normal interstellar medium (ISM) material with solar system (SS) composition. These nuclei are then accelerated to cosmic ray energies by supernova shocks. The possibility of OB association origin of GCR was first discussed by Reeves (1973) and later developed by Hainebach et al. (1976), Cassé & Paul (1982), Cesarsky & Montmerle (1981), and others. The previously measured composition of cosmic ray isotopes and elements has been shown to be consistent with GCR origin in a source which is a mixture of ∼20% material from massive star outflow and supernova ejecta, and ∼80% material with SS abundances (Higdon & Lingenfelter 2003; Binns et al. 2005, 2013; Rauch et al. 2009 ). For the remainder of this paper, the combined massive star wind outflow and supernova ejecta will be referred to as massive star material (MSM).

Observed GCR abundances show that elements that are found in interstellar dust grains (refractory elements) are preferentially accelerated compared to those that exist primarily as interstellar gases (volatile elements) (Ellison et al. 1997; Meyer et al. 1997; Meyer & Ellison 1999). Elemental abundances measured by ACE at energies of hundreds of MeV nucleon⁻¹ (Binns et al. 2013), TIGER at GeV nucleon⁻¹ (Rauch et al. 2009), and CREAM at TeV nucleon⁻¹ energies (Ahn et al. 2010) show that this enhancement is mass dependent for both refractory and volatile elements, and that the ordering of these elements with atomic mass A is greatly improved by comparing GCR source abundances with a mixture of normal ISM and MSM, rather than normal ISM alone (we note that Ellison et al. 1997 explain this mass-dependent trend only for the volatile elements).

In this paper, we demonstrate that the elemental abundances are consistent with a GCR source that consists of a best-fit mixture of ${19}_{-6}^{+11}$% MSM mixed with ∼81% material with SS abundance, and an acceleration mechanism in which elements found in interstellar dust grains are preferentially accelerated over those found in interstellar gases.

Recent γ-ray observations of supernova remnants using the Fermi Large Area Telescope (Acero et al. 2016) and ground-based imaging atmospheric Cherenkov telescope arrays such as HESS (Abramowski et al. 2015), MAGIC (Aleksic et al. 2012), and VERITAS (Humensky et al. 2015) have provided evidence of particle acceleration to very high energies. Moreover, observation of the π⁰ turn-on feature indicates that at least supernova remnants (SNRs) W44, IC443, and W51C are accelerating protons (Ackermann et al. 2013; Jogler & Funk 2016), and all of these are believed to be the remnants of core-collapse supernovae. In addition, GeV and TeV emission has been observed from acceleration within the 30 Dor C superbubble in the Large Magellanic Cloud (Abramowski et al. 2015). Fermi has detected extended emission coinciding with a "cocoon"-like morphology in the Cygnus superbubble (Ackermann et al. 2011), extending ∼50 parsecs from the Cygnus OB2 association. TeV emission from the Cygnus superbubble was also observed by ARGO-YBJ (Bartoli et al. 2014). These γ-ray observations lend further support to a model in which OB associations are a significant source of GCRs.

The SuperTIGER instrument (Binns et al. 2014; Murphy 2015) consists of two nearly identical ∼1 m × ∼2 m modules, each consisting of a suite of seven detectors. Figure 1 shows a schematic side view of one module. Three scintillator detectors measure the differential energy loss $\tfrac{{dE}}{{dx}}$ within the instrument, which is a function of a particle's charge (Z) and velocity; two Cherenkov detectors, one with an aerogel radiator (C0), and one with an acrylic radiator (C1), which give signals that are different functions of Z and velocity for particles above the Cherenkov threshold, and a scintillating fiber hodoscope (consisting of two separate x, y planes) which measures particle trajectory. For one SuperTIGER module, the Cherenkov radiators in the aerogel (C0) detector have an index of refraction n = 1.04, while the other module has n = 1.04 aerogel in half of the module and n = 1.025 aerogel in the other half. The active area of each module measures approximately 1.16 m × 2.4 m, and the full geometry factor of both modules combined is ∼8.3 m² sr for particles whose trajectory zenith angle is less than 70 degrees. After accounting for losses due to nuclear interactions within the instrument, the "effective" geometry factor is ∼3.9 m² sr for ₃₄Se.

Zoom In Zoom Out Reset image size

The charge (Z) of each cosmic-ray nucleus detected in the instrument was determined by one of two complementary techniques. At low energies (above the C1 threshold of ∼350 MeV nucleon⁻¹ but below the C0 threshold of ∼2.5 GeV nucleon⁻¹ or ∼3.3 GeV nucleon⁻¹, depending on the half-module the event went through), the charge was determined using a combination of signals from the top two scintillator detectors (S1 and S2) and the acrylic (C1) Cherenkov detector. At energies above the aerogel (C0) threshold, the charge was determined with a combination of the C1 and C0 detector signals. This technique was also used to analyze data from the 1997, 2001, and 2003 TIGER flights (Sposato 1999; Link 2003; Rauch et al. 2009).

For each scintillator and Cherenkov detector, the signal was taken as the sum of the signals from all its photomultipliers. This sum was corrected for photomultiplier gain differences, temporal variations, and area nonuniformities using the ∼5 × 10⁶ ₂₆Fe nuclei detected during flight to map the detector response. The particle trajectory from the scintillating fiber hodoscope was used to determine the particle position for the mapping correction. The trajectory angle θ with respect to the normal to the detector surface was used to correct for the sec(θ) dependence of signal on path length within the detector. Particles that underwent a charge-changing nuclear interaction within the instrument were identified and rejected by requiring agreement to within approximately one charge unit in the S1 scintillator, C1, and C0 Cherenkov detectors for events above the C0 threshold, and in the S1 and S2 scintillator detectors, and the C1 Cherenkov detector for events below the C0 threshold. Figure 2 is a cross-plot of C1 versus C0 in which each point represents the signals from the two Cherenkov detectors for a single cosmic-ray nucleus (in this plot, only events and signals from the n = 1.04 aerogels are shown). For particles with energies above the aerogel Cherenkov threshold, the combination of the two Cherenkov signals gives a well resolved charge assignment for events to the right of the red line. The points in Figure 2 to the left of the line with low C0 signals represent events with energies near or below the aerogel Cherenkov threshold. Those events are analyzed using the signals from the scintillators and the acrylic Cherenkov detector.

Zoom In Zoom Out Reset image size

Figure 3 is a similar cross-plot of signals from the top (S1) scintillator detector versus the acrylic Cherenkov (C1) detector. The events to the right of the red (right) line have energies above the aerogel Cherenkov (C0) threshold, and were analyzed using signals from the two Cherenkov detectors as described previously. The events to the left of the left (blue) line are particles near or below the acrylic Cherenkov (C1) threshold of ∼350 MeV nucleon⁻¹, and are not included in the analysis. The events between the two lines were selected for analysis. For these lower-energy events, signals from the top two scintillator detectors (S1 and S2) and the acrylic Cherenkov (C1) detector were used to assign charge. The bottom (S3) scintillator detector was used to define a more restrictive data set for the development of these analysis techniques, but was not used in the final analysis. We fit curves of constant charge to each charge contour on the cross plot, and used those to fit an energy-independent form of the scintillator saturation model of Voltz et al. (1966) using the formalism given by Ahlen (1980). The resulting charge histogram was then renormalized so that charge peaks corresponded to their integer charges.

Zoom In Zoom Out Reset image size

A charge histogram is shown in Figure 4, which includes events analyzed with both techniques from the SuperTIGER flight. The 1σ charge resolution at ₂₆Fe is 0.18 charge units (c.u.). Figure 5 uses a coarser binning of the data in the charge range 30 ≤ Z ≤ 40, showing well defined, single-element peaks for every charge in that range. This histogram was fit with a multi-Gaussian function using a maximum likelihood method that was used to derive the measured abundances shown in column 3 of Table 1.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We have derived abundances at the top of the atmosphere by correcting for the charge-dependent probability of particles undergoing nuclear interactions within the instrument and atmosphere, and the charge-dependent energy losses in the atmosphere and instrument. The correction for nuclear interactions within the instrument accounts for those particles identified and discarded during analysis. The atmospheric correction includes both the fraction of particles interacting (∼36% for ₃₄Se) and secondary production. This was done using the same technique used for the TIGER data analysis, described in detail by Rauch (2008), using the total and partial charge-changing cross sections derived from accelerator data by Nilsen et al. (1995). The derived top-of-atmosphere abundances are listed in columns 7–9 of Table 1 and shown in Figure 6. That figure also shows GCR abundances measured in space by HEAO-3-C2 (Byrnak et al. 1983), ACE-CRIS (Binns et al. 2013), and top-of-atmosphere abundances measured by TIGER (Rauch et al. 2009). Also shown are SS elemental abundances (Lodders 2003). For the SuperTIGER points, the combined statistical and systematic error bars are shown in solid orange. The systematic uncertainties for the SuperTIGER abundances were obtained by adjusting the total and partial charge changing cross-sections up and down by the uncertainty in those cross-sections (estimated by Nilsen et al. 1995), and then comparing the results of the propagation calculation to the results obtained using unmodified cross sections. These systematic uncertainties are small compared to the statistical uncertainties. The ACE-CRIS error bars are purely statistical, while the TIGER error bars are calculated using top-of-instrument statistical uncertainties propagated through the atmosphere, and the HEAO-3-C2 error bars are statistical uncertainties, reported in Byrnak et al. (1983).

Zoom In Zoom Out Reset image size

The SuperTIGER data points generally agree with previous experiments, but have significantly smaller error bars. Of particular note are ₃₁Ga and ₃₂Ge. TIGER measurements indicate nearly equal abundances for these two elements, which was not expected in view of the high SS ₃₂Ge/₃₁Ga ratio (Rauch et al. 2009). While the SuperTIGER measurement of each of these elements is not in statistical disagreement with those of TIGER, SuperTIGER, with its much better statistics, shows that the abundances of these two elements are not equal; rather ₃₂Ge/₃₁Ga is approximately 1.5, with a statistical difference in ratios of nearly 5σ.

Cosmic ray source abundances were derived from the SuperTIGER top-of-atmosphere abundances using a leaky box propagation model (Wiedenbeck et al. 2007), which uses total destruction cross sections (a modified form of those from Webber et al. 1990) and partial cross sections from Silberberg et al. (1998). The interstellar propagation results were used as input to a spherically symmetric modulation model based on Fisk (1971), with a modulation level ϕ = 543 MV and a typical top-of-atmosphere energy of ∼3.1 GeV nucleon⁻¹, to obtain modulated values for comparison with the abundances observed at Earth. This ϕ was inferred from spectra observed by ACE/CRIS during the SuperTIGER flight, which were measured using the method described in Wiedenbeck et al. (2005). The assumed cosmic-ray source abundances were adjusted to yield agreement with the data. The derived source abundances are shown in columns 3–5 of Table 2. The uncertainties reported in Table 2 are the propagated top-of-atmosphere uncertainties. Figure 7 shows the top-of-atmosphere abundances reported in Table 1 compared to the calculated source abundances.

Zoom In Zoom Out Reset image size

Figure 8 is a plot of the ratio of Galactic cosmic ray source (GCRS) elemental abundances to SS abundances from Lodders (2003) as a function of atomic mass A. For elements with Z < 26, the source abundances are those derived by Engelmann et al. (1990) from HEAO-3-C2 data. For ₂₇Co and ₂₉Cu, the abundances from Rauch et al. (2009) were used. For all other elements shown with Z > 26, each point represents the source abundances calculated for SuperTIGER combined with TIGER abundances from Rauch et al. (2009), weighted by the statistics recorded with each experiment. The refractory elements have equilibrium condensation temperatures (Lodders 2003) greater than ∼1200 K and the volatile elements have condensation temperatures lower than ∼1200 K. As noted by Meyer et al. (1997) and Ellison et al. (1997), among others, the GCRS/SS ratio is generally higher for refractory elements than for volatile elements, especially at low A. However, at high A the two groups merge and there is significant scatter, as noted by Rauch et al. (2009).

Zoom In Zoom Out Reset image size

Figure 9 is a modification of Figure 8 where the GCRS abundances are compared to a mixture by mass of 81% material with SS composition (Lodders 2003) and 19% MSM averaged over an initial mass function (calculated by Woosley & Heger 2007), instead of pure SS material. This figure shows a significant improvement in the organization of data compared with Figure 8, with a clear separation of the refractory and volatile elements, each with a similar mass dependence. The value of 19% MSM was determined by comparing the GCRS abundances with source mixtures consisting of SS material with MSM mixed in 1% increments from 0 to 100%. For each source mixture, the refractory and volatile elements were each fit with a simple curve of the form $y={C}_{0}{A}^{{C}_{1}}$, and the combined χ² value for the mixture was calculated. The mixture with the minimum total χ² was selected as the best-fit mixture, with ±1σ uncertainty levels obtained by finding the mixtures with a χ² value of ${\chi }_{\min }^{2}+1$. The best-fit mixture was found to be ${19}_{-6}^{+11} \%$ MSM by mass, with the rest being normal ISM material with SS elemental abundances. Figure 10 shows the total χ² for refractory elements, volatile elements, and the combined total χ² as a function of the percentage of MSM in the source mixture. The vertical dashed lines show the ±1σ range in the percentage of MSM.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 11 uses the same data and curves as Figure 9, but with the addition of a ₂He datum from Cummings et al. (2016). This point was not included in the fit, but the point still falls very near to the best-fit line.

Zoom In Zoom Out Reset image size

We use SS abundances from Lodders (2003) instead of the more recent abundances from Lodders et al. (2009) because the massive star outflow and ejecta model of Woosley & Heger (2007) used the Lodders (2003) relative abundances as the input to their yield calculations. The more recent calculation of SN yields by Sukhbold et al. (2016) also uses Lodders (2003). We use the Woosley & Heger (2007) results here instead of Sukhbold et al. (2016) because, as noted in their paper, the Sukhbold et al. (2016) SN yields are "problematic" for elements in the ultra-heavy charge range. We note that yields have also been calculated by Chieffi & Limongi (2013). Their yields for these UH elements differ from Woosley & Heger (2007) and show strong dependence on the choice of mass cut, which determines how much material is ejected in the supernova explosion. A similar analysis should be done using their calculated yields. However, this is beyond the scope of the present paper.

This best fit had a reduced χ² value of 1.26. The C₁ parameter of the best fit is 0.583 ± 0.072 for refractory elements and 0.632 ± 0.119 for volatile elements. Figure 12 shows the fit C₁ values as a function of the percent of MSM in the source mixture.

Zoom In Zoom Out Reset image size

We interpret Figures 9 and 11 as strong evidence in support of the model of cosmic ray origin in which the source material is a mix of MSM with normal ISM, and refractory elements are preferentially accelerated over volatile elements. The contribution to the GCRS mixture from MSM indicates that OB associations are a significant source of GCR.

In addition, the recent detection of ⁶⁰Fe in cosmic rays (⁶⁰Fe is a radioactive primary cosmic ray with half-life 2.6 Myr that is primarily synthesized in core-collapse supernovae; Travaglio et al. 2004; Seitenzahl et al. 2013; Sukhbold et al. 2016) conclusively shows that recently synthesized material (within the last few million years) is accelerated to cosmic ray energies (Binns et al. 2016). The most natural place for this to occur is in OB associations.

These observations do not directly tell us where or how the normal ISM is injected and accelerated into cosmic rays. SNe Ia often explode into the normal ISM and it is estimated that 15% of core-collapse supernovae occur outside of superbubbles (Higdon & Lingenfelter 2003, 2005). γ-ray observations show that these supernovae accelerate high-energy cosmic rays (Wang et al. 2007; Dermer & Powale 2013). However, additional accelerators are required since an unreasonably large fraction of the supernova energy would be required to power cosmic rays if these were the only source of normal ISM acceleration. Since most supernovae are in OB associations, it appears that supernova shocks from stars in OB associations must also be accelerating cosmic rays from the normal ISM, perhaps from walls of superbubbles surrounding OB associations and residual ISM within the superbubble itself (Higdon & Lingenfelter 2003).

We have presented new SuperTIGER measurements of the elemental abundances of GCRs from ₂₆Fe to ₄₀Zr. Our results support a model of cosmic ray origin in a source mixture of ${19}_{-6}^{+11} \%$ MSM and ∼81% normal ISM material with SS abundances. This indicates that a significant fraction of GCR acceleration occurs in OB associations. We also find a preferential acceleration of refractory elements over volatile elements by a factor of between ∼4 and ∼4.5, ordered by atomic mass (A). Both the refractory and volatile elements show a mass-dependent enhancement with similar slopes.

We gratefully acknowledge the excellent and highly professional work of the NASA Columbia Scientific Balloon Facility, the NASA Balloon Program Office, and the NSF Office of Polar Programs, who together made possible the record long-duration balloon flight of SuperTIGER. NASA supported this research under the ROSES 2007 APRA program under grants NNX09AC17G to Washington University in St. Louis and NNX09AC18G to the California Institute of Technology and the Jet Propulsion Laboratory, and APRA07-0146 to NASA/GSFC. We thank Nasser Barghouty for his help calculating interaction cross sections for the Galactic propagation. We are also grateful for support from the Peggy and Steve Fossett Foundation and the McDonnell Center for the Space Sciences at Washington University. We thank the referee who identified himself as Don Ellison for his helpful comments that improved this paper.