The most common proxy of solar activity is formed by the data on cosmogenic radionuclides (e.g., 10Be and 14C), which are produced by cosmic rays in the Earth’s atmosphere (e.g, Stuiver and Quay, 1980; Beer et al., 1990; Bard et al., 1997; Beer, 2000). Other cosmogenic nuclides, which are used in geological and paleomagnetic dating, are less suitable for studies of solar activity (see e.g., Beer, 2000; Beer et al., 2012). Cosmic rays are the main source of cosmogenic nuclides in the atmosphere (excluding anthropogenic factors during the last decades) with the maximum production being in the upper troposphere/stratosphere. After a complicated transport in the atmosphere, the cosmogenic isotopes are stored in natural archives such as polar ice, trees, marine sediments, etc. This process is also affected by changes in the geomagnetic field and climate. Cosmic rays experience heliospheric modulation due to solar wind and the frozen-in solar magnetic field. The intensity of modulation depends on solar activity and, therefore, cosmic-ray flux and the ensuing cosmogenic isotope intensity depends inversely on solar activity. An important advantage of the cosmogenic data is that primary archiving is done naturally in a similar manner throughout the ages, and these archives are measured nowadays in laboratories using modern techniques. If necessary, all measurements can be repeated and improved, as has been done for some radiocarbon samples. In contrast to fixed historical archival data (such as sunspot or auroral observations) this approach makes it possible to obtain homogeneous data sets of stable quality and to improve the quality of data with the invention of new methods (such as accelerator mass spectrometry). Cosmogenic isotope data is the main regular indicator of solar activity on the very long-term scale but it cannot resolve the details of individual solar cycles. The redistribution of nuclides in terrestrial reservoirs and archiving may be affected by local and global climate/circulation processes, which are, to a large extent, unknown for the past. However, a combined study of different nuclides data, whose responses to terrestrial effects are very different, may allow for disentangling external and terrestrial signals.
The flux of cosmic rays (highly energetic fully ionized nuclei) is considered roughly constant (at least at the time scales relevant for the present study) in the vicinity of the Solar system. However, before reaching the vicinity of Earth, galactic cosmic rays experience complicated transport in the heliosphere that leads to modulation of their flux. Heliospheric transport of GCR is described by Parker’s theory (Parker, 1965; Toptygin, 1985) and includes four basic processes: the diffusion of particles due to their scattering on magnetic inhomogeneities, the convection of particles by out-blowing solar wind, adiabatic energy losses in expanding solar wind, drifts of particles in the magnetic field, including the gradient-curvature drift in the regular heliospheric magnetic field, and the drift along the heliospheric current sheet, which is a thin magnetic interface between the two heliomagnetic hemispheres. Because of variable solar-magnetic activity, CR flux in the vicinity of Earth is strongly modulated (see Figure 3). The most prominent feature in CR modulation is the 11-year cycle, which is in inverse relation to solar activity. The 11-year cycle in CR is delayed (from a month up to two years) with respect to the sunspots (Usoskin et al., 1998). The time profile of cosmic-ray flux as measured by a neutron monitor (NM) is shown in Figure 3 (panel b) together with the sunspot numbers (panel a). Besides the inverse relation between them, some other features can also be noted. A 22-year cyclicity manifests itself in cosmic-ray modulation through the alteration of sharp and flat maxima in cosmic-ray data, originated from the charge-dependent drift mechanism. One may also note short-term fluctuations, which are not directly related to sunspot numbers but are driven by interplanetary transients caused by solar eruptive events, e.g., flares or CMEs. An interesting feature is the increase of CR flux in 2009, when it was the highest ever recorded by NMs (Moraal and Stoker, 2010), as caused by the favorable heliospheric conditions (unusually weak heliospheric magnetic field and the flat heliospheric current sheet) (McDonald et al., 2010). For the previous 50 years of high and roughly-stable solar activity, no trends have been observed in CR data; however, as will be discussed later, the overall level of CR has changed significantly on the centurial-millennial timescales.
Full solution of the CR transport problems is a complicated task and requires sophisticated 3D time-dependent self-consistent modelling. However, the problem can be essentially simplified for applications at a long-timescale. An assumption on the azimuthal symmetry (requires times longer that the solar-rotation period) and quasi-steady changes reduces it to a 2D quasi-steady problem. Further assumption of the spherical symmetry of the heliosphere reduces the problem to a 1D case. This approximation can be used only for rough estimates, since it neglects the drift effect, but it is useful for long-term studies, when the heliospheric parameters cannot be evaluated independently. Further, but still reasonable, assumptions (constant solar-wind speed, roughly power-law CR energy spectrum, slow spatial changes of the CR density) lead to the force-field approximation (Gleeson and Axford, 1968), which can be solved analytically. The differential intensity of the cosmic-ray nuclei of type with kinetic energy at 1 AU is given in this case asGarcia-Munoz et al., 1975; Burger et al., 2000; Webber and Higbie, 2003, 2009) agree with each other for energies above 20 GeV but may contain uncertainties of up to a factor of 1.5 around 1 GeV. These uncertainties in the boundary conditions make the results of the modulation theory slightly model-dependent (see discussion in Usoskin et al., 2005; Herbst et al., 2010) and require the LIS model to be explicitly cited. This approach gives results, which are at least dimensionally consistent with the full theory and can be used for long-term studies1 (Usoskin et al., 2002b; Caballero-Lopez and Moraal, 2004). Differential CR intensity is described by the only time-variable parameter, called the modulation potential , which is mathematically interpreted as the averaged rigidity (i.e., the particle’s momentum per unit of charge) loss of a CR particle in the heliosphere. However, it is only a formal spectral index whose physical interpretation is not straightforward, especially on short timescales and during active periods of the sun (Caballero-Lopez and Moraal, 2004). Despite its cloudy physical meaning, this force-field approach provides a very useful and simple single-parametric approximation for the differential spectrum of GCR, since the spectrum of different GCR species directly measured near the Earth can be perfectly fitted by Equation 3 using only the parameter in a wide range of solar activity levels (Usoskin et al., 2011). Therefore, changes in the whole energy spectrum (in the energy range from 100 MeV/nucleon to 100 GeV/nucleon) of cosmic rays due to the solar modulation can be described by this single number within the framework of the adopted LIS. The concept of modulation potential is a key concept for the method of solar-activity reconstruction by cosmogenic isotope proxy as it makes it possible to parameterize the GCR with one single parameter.
Cosmic rays are charged particles and therefore are affected by the Earth’s magnetic field. Thus the geomagnetic field puts an additional shielding on the incoming flux of cosmic rays. The shielding effect of the geomagnetic field is usually expressed in terms of the cutoff rigidity , which is the minimum rigidity a vertically incident CR particle must posses (on average) in order to reach the ground at a given location and time (Cooke et al., 1991). Neglecting such effects as the East-West asymmetry, which is roughly averaged out for the isotropic particle flux, or nondipole magnetic momenta, which decay rapidly with distance, one can come to a simple approximation, called the Störmer’s equation, that describes the vertical geomagnetic cutoff rigidity :25 G cm3), is the Earth’s mean radius, is the distance from the given location to the dipole center, and is the geomagnetic latitude. The cutoff concept works like a Heaviside step-function so that all cosmic rays whose rigidity is below the cutoff are not allowed to enter the atmosphere while all particles with higher rigidity can penetrate. This approximation provides a good compromise between simplicity and reality, especially when using the eccentric dipole description of the geomagnetic field (Fraser-Smith, 1987). The eccentric dipole has the same dipole moment and orientation as the centered dipole, but the dipole’s center and consequently the poles, defined as crossings of the axis with the surface, are shifted with respect to geographical ones.
The shielding effect is the strongest at the geomagnetic equator, where the present-day value of may reach up to 17 GV in the region of India. There is almost no cutoff in the geomagnetic polar regions (). However, even in the latter case the atmospheric cutoff becomes important, i.e., particles must have rigidity above 0.5 GV in order to initiate the atmospheric cascade which can reach ground (see Section 3.1.3).
The geomagnetic field is seemingly stable on the short-term scale, but it changes essentially on centurial-to-millennial timescales (e.g., Korte and Constable, 2006). Such past changes can be evaluated based on measurements of the residual magnetization of independently-dated samples. These can be paleo- (i.e., natural stratified archives such as lake or marine sediments or volcanic lava) or archaeological (e.g., clay bricks that preserve magnetization upon baking) samples. Most paleo-magnetic data preserve not only the magnetic field intensity but also the direction of the local field, while archeo-magnetic samples provide information on the intensity only. Using a large database of such samples, it is possible to reconstruct (under reasonable assumptions) the large-scale magnetic field of the Earth. Data available provides good global coverage for the last 3 millennia, allowing for a reliable paleomagnetic reconstruction of the true dipole moment (DM) or virtual dipole moment2 (VDM) and its orientation (the ArcheoInt collaboration – Genevey et al., 2008). Less precise, but still reliable reconstructions of the DM and its orientation are possible for the last seven millennia (the CALS7K.2 model by Korte and Constable, 2005), however they may slightly underestimate the dipole moment, especially in the earlier part of the period (Korte and Constable, 2008). Directional paleomagnetic reconstruction are less reliable on a longer timescale, because of the spatial sparseness of the paleo/archeo-magnetic samples in the earlier part of the Holocene (Korte et al., 2011). Some paleomagnetic reconstructions are shown in Figure 4. All paleomagnetic models depict a similar long-term trend – an enhanced intensity during the period between 1500 BC and 500 AD and a significantly lower field before that.
Changes in the dipole moment inversely modulate the flux of CR at Earth, with strong effects in tropical regions and globally. The migration of the geomagnetic axis, which changes the geomagnetic latitude of a given geographical location is also important; while not affecting the global flux of CR, it can dramatically change the CR effect regionally, especially at middle and high latitude. These changes affect the flux of CR impinging on the Earth’s atmosphere both locally and globally and must be taken into account when reconstructing solar activity from terrestrial proxy data (Usoskin et al., 2008, 2010). Accounting for these effects is quite straightforward provided the geomagnetic changes in the past are known independently, e.g., from archeo and paleo-magnetic studies (Donadini et al., 2010). However, because of progressively increasing uncertainties of paleomagnetic reconstructions back in time, it presently forms the main difficulty for the proxy method on the long-term scale (Snowball and Muscheler, 2007), especially in the early part of the Holocene. On the other hand, the geomagnetic field variations are relatively well known for the last few millennia (Genevey et al., 2008; Korte and Constable, 2008).
When an energetic CR particle enters the atmosphere, it first moves straight in the upper layers, suffering mostly from ionization energy losses that lead to the ionization of the ambient rarefied air and gradual deceleration of the particles. However, after traversing some amount of matter (the nuclear interaction mean-free path is on the order of 100 g/cm2 for a proton in the air) the CR particle may collide with a nucleus in the atmosphere, producing a number of secondaries. These secondaries have their own fate in the atmosphere, in particular they may suffer further collisions and interactions forming an atmospheric cascade (e.g., Dorman, 2004). Because of the thickness of the Earth’s atmosphere (1033 g/cm2 at sea level) the number of subsequent interactions can be large, leading to a fully-developed cascade (also called an air shower) consisting of secondary rather than primary particles. A schematic view of the atmospheric cascade is shown in Figure 5.Three main components can be separated in the cascade:
- The “hadronic” nucleonic component is formed by the products of nuclear collisions of primary cosmic rays and their secondaries with the atmospheric nuclei, and consists mostly of superthermal protons and neutrons.
- The “soft” or electromagnetic component consists of electrons, positrons and photons.
- The “hard” or muon component consists mostly of muons; pions are short lived and decay almost immediately upon production, feeding muons and the “soft” component.
The development of the cascade depends mostly on the amount of matter traversed and is usually linked to residual atmospheric depth, which is very close to the static barometric pressure, rather than to the actual altitude, that may vary depending on the exact atmospheric density profile.
Cosmogenic isotopes are a by-product of the hadronic branch of the cascade (details are given below). Accordingly, in order to evaluate cosmic-ray flux from the cosmogenic isotope data, one needs to know the physics of cascade development. Several models have been developed for this cascade, in particular its hadronic branch with emphasis on the generation of cosmogenic isotope production. The first models were simplified quasi-analytical (e.g., Lingenfelter, 1963; O’Brien and Burke, 1973) or semi-empirical models (e.g., Castagnoli and Lal, 1980). With the fast advance of computing facilities it became possible to exploit the best numerical method suitable for such problems – Monte-Carlo (e.g., Masarik and Beer, 1999, 2009; Webber and Higbie, 2003; Webber et al., 2007; Usoskin and Kovaltsov, 2008a; Kovaltsov and Usoskin, 2010; Kovaltsov et al., 2012; Argento et al., 2013). The fact that models, based on different independent Monte-Carlo packages, namely, a general GEANT tool and a specific CORSIKA code, yield similar results provides additional verification of the approach.
A scheme for the transport and redistribution of the two most useful cosmogenic isotopes, 14C and 10Be, is shown in Figure 6. After a more-or-less similar production, the two isotopes follow quite different fates, as discussed in detail in Sections 3.2.3 and 3.3.3. Therefore, expected terrestrial effects are quite different for the isotopes and comparing them with each other can help in disentangling solar and climatic effects (see Section 3.7.3). A reader can find great detail also in a book by Beer et al. (2012).
The most commonly used cosmogenic isotope is radiocarbon 14C. This radionuclide is an unstable isotope of carbon with a half-life of about 5730 years. Since the radiocarbon method is extensively used in other science disciplines where accurate dating is a key issue (e.g., archeology, paleoclimatology, quaternary geology), it was developed primarily for this task. The solar-activity–reconstruction method, based on radiocarbon data, was initially developed as a by-product of the dating techniques used in archeology and Quaternary geology, in an effort to improve the quality of the dating by means of better information on the 14C variable source function. The present-day radiocarbon calibration curve, based on a dendrochronological scale, uninterruptedly covers the whole Holocene (and extending to 50,000 BP – Reimer et al., 2009) and provides a solid quantitative basis for studying solar activity variations on the multi-millennial time scale.
Radiocarbon is usually measured in tree rings, which allows an absolute dating of the samples by means of dendrochronology. Using a complicated technique, the 14C activity3 is measured in an independently dated sample, which is then corrected for age as13C isotope) of the sample, the radiocarbon value of 14C is calculated (see details in Stuiver and Pollach, 1977). 13C is the per mille deviation of the 13C content in the sample from that in the standard belemnite sample calculated similarly to Equation 6. The value of 14C (measured in per mille ) is further used as the index of radiocarbon relative activity. The series of 14C for the Holocene is presented in Figure 7A as published by the INTCAL04 collaboration of 21 dating laboratories as a result of systematic precise measurements of dated samples from around the world (Reimer et al., 2004). The most recent INTCAL09 dataset is available at http://www.radiocarbon.org/IntCal09.htm.
A potentially interesting approach has been made by Lal et al. (2005), who measured the amount of 14C directly produced by CR in polar ice. Although this method is free of the carbon-cycle influence, the first results, while being in general agreement with other methods, are not precise.
The main source of radioisotope 14C (except anthropogenic sources during the last decades) is cosmic rays in the atmosphere. It is produced as a result of the capture of a thermal neutron by atmospheric nitrogen3.1.3) but their flux varies in time along with the modulation of cosmic-ray flux. This provides continuous source of the isotope in the atmosphere, while the sinks are isotope decay and transport into other reservoirs as described below (the carbon cycle).
The connection between the cosmogenic-isotope–production rate, , at a given location (quantified via the geomagnetic latitude ) and the cosmic-ray flux is given by3.1.2), is the differential energy spectrum of CR (see Section 3.1.1) and is the differential yield function of cosmogenic isotope production, calculated using a Monte-Carlo simulation of the cosmic-ray–induced atmospheric cascade (Kovaltsov et al., 2012). Because of the global nature of the carbon cycle and its long attenuation time, the radiocarbon is globally mixed before the final deposition, and Equation 9 should be integrated over the globe. The yield function of the 14C production is shown in Figure 8A together with those for 10Be (see Section 3.3.2) and for a ground-based neutron monitor (NM), which is the main instrument for studying cosmic-ray variability during the modern epoch. One can see that the yield function increases with the energy of CR. On the other hand, the energy spectrum of CR decreases with energy. Accordingly, the differential production rate (i.e., the product of the yield function and the spectrum, – the integrand of Equation 9), shown in Figure 8B, is more informative. The differential production rate reflects the sensitivity to cosmic rays, and the total production rate is simply an integral of over energy above the geomagnetic threshold.
Thanks to the development of atmospheric cascade models (Section 3.1.3), there are numerical models that allow one to compute the radiocarbon production rate as a function of the modulation potential and the geomagnetic dipole moment . The overall production of 14C is shown in Figure 9.
The production rate of radiocarbon, , can vary as affected by different factors (see, e.g., Damon and Sonett, 1991):
- Variations of the cosmic-ray flux on a geological timescale due to the changing galactic background (e.g., a nearby supernova explosion or crossing the dense galactic arm).
- Secular-to-millennial variations are caused by the slowly-changing geomagnetic field. This is an important component of the variability, which needs to be independently evaluated from paleo and archeo-magnetic studies.
- Modulation of cosmic rays in the heliosphere by solar magnetic activity. This variation is the primary aim of the present method.
- Short-term variability of CR on a daily scale (suppression due to interplanetary transients or enhancement due to solar energetic-particle events) can be hardly resolved in radiocarbon data.
Therefore, the production rate of 14C in the atmosphere can be modelled for a given time (namely, the modulation potential and geomagnetic dipole moment) and location. The global production rate is then obtained as a result of global averaging.
There is still a small unresolved discrepancy in the absolute value of the modeled 14C production rate. Different models yield the global-average production rate of 1.7 – 2.3 atoms cm–2 s–1 (see, e.g., O’Brien, 1979; Masarik and Beer, 1999; Goslar, 2001; Usoskin et al., 2006b; Kovaltsov et al., 2012, and references therein), which is consistent with a direct estimate of the radiocarbon reservoir, based on analyses of the specific 14C activity on the ground, 1.76 – 2.0 (Suess, 1965; Damon et al., 1978; O’Brien, 1979). On the other hand, the steady-state production calculated from the 14C inventory in the carbon-cycle model (see Section 3.2.3) typically yields 1.6 – 1.7 atoms cm–2 s–1 for the pre-industrial period (e.g., Goslar, 2001, and references therein). Thus, results obtained from the carbon cycle models and production models agree only marginally in the absolute values, and this needs further detailed studies. The situation has been improved recently, when the newest numerical model (Kovaltsov et al., 2012) yields the global average production of 1.64 and 1.88 atoms cm–2 s–1 for the modern and pre-industrial periods, respectively. In 14C-based reconstructions, the pre-industrial steady-state production is commonly used.
Upon production cosmogenic radiocarbon gets quickly oxidized to carbon dioxide CO2 and takes part in the regular carbon cycle of interrelated systems: atmosphere-biosphere-ocean (Figure 6). Because of the long residence time, radiocarbon becomes globally mixed in the atmosphere and involved in an exchange with the ocean. It is common to distinguish between an upper layer of the ocean, which can directly exchange CO2 with the air and deeper layers. The measured 14C comes from the biosphere (trees), which receives radiocarbon from the atmosphere. Therefore, the processes involved in the carbon cycle are quite complicated. The carbon cycle is usually described using a box model (Oeschger et al., 1974; Siegenthaler et al., 1980), where it is represented by fluxes between different carbon reservoirs and mixing within the ocean reservoir(s), as shown in Figure 10. Production and radioactive decay are also included in box models. Free parameters in a typical box model are the 14C production rate , the air-sea exchange rate (expressed as turnover rate ), and the vertical–eddy-diffusion coefficient , which quantifies ocean ventilation. Starting from the original representation (Oeschger et al., 1974), a variety of box models have been developed, which take into account subdivisions of the ocean reservoir and direct exchange between the deep ocean and the atmosphere at high latitudes. More complex models, including a diffusive approach, are able to simulate more realistic scenarios, but they require knowledge of a large number of model parameters. These parameters can be evaluated for the present time using the bomb test – studying the transport and distribution of the radiocarbon produced during the atmospheric nuclear tests. However, for long-term studies, only the production rate is considered variable, while the gas-exchange rate and ocean mixing are kept constant. Under such assumptions, there is no sense in subdividing reservoirs or processes, and a simple carbon box model is sufficient.
Using the carbon cycle model and assuming that all its parameters are constant in time, one can evaluate the production rate from the measured 14C data. This assumption is well validated for the the Holocene (Damon et al., 1978; Stuiver et al., 1991) as there is no evidence of considerable oceanic change or other natural variability of the carbon cycle (Gerber et al., 2002), and accordingly all variations of 14C predominantly reflect the production rate. This is supported by the strong similarity of the fluctuations of the 10Be data in polar ice cores (Section 3.3) compared to 14C, despite their completely different geochemical fate (Bard et al., 1997; Steinhilber et al., 2012). However, the changes in the carbon cycle during the last glaciation and deglaciation were dramatic, especially regarding ocean ventilation; this and the lack of independent information about the carbon cycle parameters, make it hardly possible to qualitatively estimate solar activity from 14C before the Holocene.
First attempts to extract information on production-rate variations from measured 14C were based on simple frequency separations of the signals. All slow changes were ascribed to climatic and geomagnetic variations, while short-term fluctuations were believed to be of solar origin. This was done by removing the long-term trend from the 14C series and claiming the residual as being a series of solar variability (e.g., Peristykh and Damon, 2003). This oversimplified approach was natural at earlier times, before the development of carbon cycle models, but later it was replaced by the inversion of the carbon cycle (i.e., the reconstruction of the production rate from the measured 14C concentration). Although mathematically this problem can be solved correctly as a system of linear differential equations, the presence of fluctuating noise with large magnitude makes it not straightforward, since the time derivative cannot be reliably identified leading thus to possible amplification of the high-frequency noise in 14C data. One traditional approach (e.g., Stuiver and Quay, 1980) is based on an iterative procedure, first assuming a constant production rate, and then fitting the calculated 14C variations to the actual measurements using a feedback scheme. A concurrent approach based on the presentation of the carbon cycle as a Fourier filter (Usoskin and Kromer, 2005) produces similar results. Roughly speaking, the carbon cycle acts as an attenuating and delaying filter for the 14C signal (see Figure 11). The higher the frequency is, the greater the signal is attenuated. In particular, the large 11-year solar cycle expected in the 14C is attenuated by a factor of hundred in the measured 14C data, making it hardly detectable. Because of the slow oceanic response, the 14C data is also delayed with respect to the production signal. The production rate for the Holocene is shown in Figure 7 and depicts both short-term fluctuations as well as slower variations, mostly due to geomagnetic field changes (see Section 3.2.5).
Unfortunately, cosmogenic 14C data cannot be easily used for the last century, primarily because of the extensive burning of fossil fuels. Since fossil fuels do not contain 14C, the produced CO2 dilutes the atmospheric 14CO2 concentration with respect to the pre-industrial epoch. Therefore, the measured 14C cannot be straightforwardly translated into the production rate after the late 19th century, and a special correction for fossil fuel burning is needed. This effect, known as the Suess effect (e.g., Suess, 1955), can be up to in 14C in 1950 (Tans et al., 1979), which is an order of magnitude larger than the amplitude of the 11-year cycle of a few . Moreover, while the cosmogenic production of 14C is roughly homogeneous over the globe and time, the use of fossil fuels is highly nonuniform (e.g., de Jong and Mook, 1982) both spatially (developed countries, in the northern hemisphere) and temporarily (World Wars, Great Depression, industrialization, etc.). This makes it very difficult to perform an absolute normalization of the radiocarbon production to the direct measurements. Sophisticated numerical models (e.g., Sabine et al., 2004; Mikaloff Fletcher et al., 2006) aim to account for the Suess effect and make good progress. However, the results obtained indicate that the determination of the Suess effect does not yet reach the accuracy required for the precise modelling and reconstruction of the 14C production for the industrial epoch. As noted by Matsumoto et al. (2004), “…not all is well with the current generation of ocean carbon cycle models. At the same time, this highlights the danger in simply using the available models to represent state-of-the-art modeling without considering the credibility of each model.” Note that the atmospheric concentration of another carbon isotope 13C is partly affected by land use, which has also been modified during the last century.
Another anthropogenic activity greatly disturbing the natural variability of 14C is related to the atmospheric nuclear bomb tests actively performed in the 1960s. For example, the radiocarbon concentration nearly doubled in the early 1960s in the northern hemisphere after nuclear tests performed by the USSR and the USA in 1961 (Damon et al., 1978). On one hand, such sources of momentary spot injections of radioactive tracers (including 14C) provide a good opportunity to verify and calibrate the exchange parameters for different carbon-cycle reservoirs and circulation models (e.g., Bard et al., 1987; Sweeney et al., 2007). Thus, the present-day carbon cycle is more-or-less known. On the other hand, the extensive additional production of isotopes during nuclear tests makes it hardly possible to use the 14C as a proxy for solar activity after the 1950s (Joos, 1994).
These anthropogenic effects do not allow one to make a straightforward link between pre-industrial data and direct experiments performed during more recent decades. Therefore, the question of the absolute normalization of 14C model is still open (see, e.g., the discussion in Solanki et al., 2004, 2005; Muscheler et al., 2005).
As discussed in Section 3.1.2, knowledge of geomagnetic shielding is an important aspect of the cosmogenic isotope method. Since radiocarbon is globally mixed in the atmosphere before deposition, its production is affected by changes in the geomagnetic dipole moment , while magnetic-axis migration plays hardly any role in 14C data.
The crucial role of paleomagnetic reconstructions has long been known (e.g., Elsasser et al., 1956; Kigoshi and Hasegawa, 1966). Many earlier corrections for possible geomagnetic-field changes were performed by detrending the measured 14C abundance or production rate (Stuiver and Quay, 1980; Voss et al., 1996; Peristykh and Damon, 2003), under the assumption that geomagnetic and solar signals can be disentangled from the production in the frequency domain. Accordingly, the temporal series of either measured 14C or its production rate is decomposed into the slow changing trend and faster oscillations. The trend is supposed to be entirely due to geomagnetic changes, while the oscillations are ascribed to solar variability. Such a method, however, obliterates all information on possible long-term variations of solar activity. Simplified empirical correction factors were also often used (e.g., Stuiver and Quay, 1980; Stuiver et al., 1991). The modern approach is based on a physics-based model (e.g., Solanki et al., 2004; Vonmoos et al., 2006) and allows the quantitative reconstruction of solar activity, explicitly using independent reconstructions of the geomagnetic field. In this case the major source of errors in solar activity reconstructions is related to uncertainties in the paleomagnetic data (Snowball and Muscheler, 2007). These errors are insignificant for the last millennium (Usoskin et al., 2006a), but become increasingly important for earlier times.
The cosmogenic isotope 10Be is useful for long-term studies of solar activity because of its long half-life of around 1.5 × 106 years. Its concentration is usually measured in stratified ice cores allowing for independent dating. Because of its long life, the beryllium concentration is difficult to measure by the decay rate (Beer, 2000). Accordingly, the 10Be/9Be ratio needs to be precisely measured at an accuracy better than 10–13. This can be done using AMS (Accelerator Mass Spectrometry) technique, which makes the measurements complicated and expensive. Correction for the decay is straightforward and does not include isotope fractionating. From the measured samples, first the 10Be concentration is defined, usually in units of 104 atoms/g. Sometimes, a correction for the snow precipitation amount is considered leading to the observable 10Be flux, which is the number of atoms, precipitating to the surface per cm2 per second.
There exist different 10Be series suitable for studies of long-term solar activity, coming from ice cores in Greenland and Antarctica. They have been obtained from different cores with different resolutions, and include data from Milcent, Greenland (Beer et al., 1983), Camp Century, Greenland (Beer et al., 1988), Dye 3, Greenland (Beer et al., 1990), Dome Concordia and South Pole, Antarctica (Raisbeck et al., 1990), GRIP, Greenland (Yiou et al., 1997), GISP2, Greenland (Finkel and Nishiizumi, 1997), Dome Fuji, Antarctica (Horiuchi et al., 2007, 2008), Dronning Maud Land, Antarctica (Ruth et al., 2007), etc. We note that data on 10Be in other archives, e.g., lake sediments, is usually more complicated to interpret because of the potential influence of the climate (Horiuchi et al., 1999; Belmaker et al., 2008).
The isotope 10Be is produced as a result of spallation of atmospheric nitrogen and oxygen (carbon is less abundant by far in the atmosphere and makes a negligible contribution) by the nucleonic component of the cosmic-ray–induced atmospheric cascade (Section 3.1.3).
A small contribution may also exist from photo-nuclear reactions (Bezuglov et al., 2012). The cross section (a few mb) of the spallation reactions is almost independent of the energy of impacting particles and has a threshold of about 15 MeV. Thus, the production of 10Be is defined mostly by the multiplicity of the nucleonic component, which increases with the energy of primary cosmic rays (see Figure 8). Maximum production occurs at an altitude of 10 – 15 km due to a balance between the total energy of the cascade (which increases with altitude) and the number of secondaries (decreasing with altitude). Most of the global 10Be is produced in the stratosphere (55 – 70%) and the rest in the troposphere (Lal and Peters, 1967; Masarik and Beer, 1999, 2009; Usoskin and Kovaltsov, 2008a; Kovaltsov and Usoskin, 2010).
Computation of 10Be isotope production is straightforward, provided a model of the atmospheric cascade is available. The first consistent model was developed by D. Lal et al. (Bhandari et al., 1966; Lal and Peters, 1967; Lal and Suess, 1968), using an empirical approach based on fitting simplified model calculations to measurements of the isotope concentrations and “star” (inelastic nuclear collisions) formations in the atmosphere. Next was an analytical model by O’Brien (1979), who solved the problem of the GCR-induced cascade in the atmosphere using an analytical stationary approximation in the form of the Boltzmann equation. Those models were based on calculating the rate of inelastic collisions or “stars” and then applying the mean spallation yield per “star”. A new step in the modelling of isotope production was made by Masarik and Beer (1999), who performed a full Monte-Carlo simulation of a GCR-initiated cascade in the atmosphere and used cross sections of spallation reactions directly instead of the average “star” efficiency. Modern models (Webber and Higbie, 2003; Webber et al., 2007; Usoskin and Kovaltsov, 2008a; Kovaltsov and Usoskin, 2010) are based on a full Monte-Carlo simulation of the atmospheric cascade, using improved cross sections. The global production rate of 10Be is about 0.02 – 0.03 atoms cm–2 s–1 (Masarik and Beer, 1999; Webber et al., 2007; Kovaltsov and Usoskin, 2010), which is lower than that for 14C by two orders of magnitude (about 2 atoms cm–2 s–1; see Section 3.2.2). The yield function of 10Be production is shown in Figure 8A and the differential production rate in Figure 8B. One can see that the peak of 10Be sensitivity, especially in polar regions, is shifted towards lower energies (below 1 GeV) compared with a neutron monitor. This implies that the 10Be isotope is relatively more sensitive to less energetic CR and is, therefore, more affected by solar energetic particles (Usoskin et al., 2006b). Comparison of model computations with direct beryllium production experiments (Usoskin and Kovaltsov, 2008a; Kovaltsov and Usoskin, 2010), and also the results of modelling of the short-living 7Be isotope (Usoskin et al., 2009a) suggest that some numerical models (Masarik and Beer, 1999; Webber and Higbie, 2003; Webber et al., 2007) tend to underestimate the production.
Although the production of 10Be can be more or less precisely modelled, a simple normalization “surface”, similar to that shown in Figure 9 for 14C, is not easy to produce because of partial mixing in the atmosphere (see Section 3.3.3). Simplified models, assuming either only global (e.g., Beer, 2000) or polar production (Bard et al., 1997; Usoskin et al., 2004), have been used until recently. However, it has been recognized that a more realistic model of the limited atmospheric mixing should be used. Without detailed knowledge of 10Be transport in the atmosphere, it is impossible to relate the quantitatively-measured concentration to the production (as done for 14C using the carbon cycle), and one has to assume that the measured abundance is proportional (with an unknown coefficient) to the production rate in a specific geographical region (see Section 3.3.3).
After production, the 10Be isotope has a seemingly simple (Figure 6) but difficult-to-account-for fate in the atmosphere. Its atmospheric residence time depends on scavenging, stratosphere-troposphere exchange and inter-tropospheric mixing (e.g., McHargue and Damon, 1991). Soon after production, the isotope becomes attached to atmospheric aerosols and follows their fate. In addition, it may be removed from the lower troposphere by wet deposition (rain and snow). The mean residence time of the aerosol-bound radionuclide in the atmosphere is quite different for the troposphere, being a few weeks, and stratosphere, where it is one to two years (Raisbeck et al., 1981). Accordingly, 10Be produced in the troposphere is deposited mostly locally, i.e., in the polar regions, while stratospheric 10Be can be partly or totally mixed. In addition, because of the seasonal (usually Spring) intrusion of stratospheric air into the troposphere at mid-latitudes, there is an additional contribution of stratospheric 10Be. Therefore, the measured 10Be concentration (or flux) in polar ice is modulated not only by production but also by climate/precipitation effects (e.g., Steig et al., 1996; Bard et al., 1997). This led Lal (1987) to the extreme conclusion that variations of polar 10Be reflect a meteorological, rather than solar, signal. However, comparison between Greenland and Antarctic 10Be series and between 10Be and 14C data (e.g., Bard et al., 1997; Horiuchi et al., 2008; Beer et al., 2012) suggests that the beryllium data mostly depicts production variations (i.e., solar signal) on top of which some meteorological effects can be superposed (see also Section 3.7.3).
Since both assumptions of the global and purely-local polar production of 10Be archived in polar ice are over-simplified, several attempts have been made to overcome this problem. For instance, McCracken (2004) proposed several simple mathematical models of partial atmospheric mixing (without division in the troposphere and stratosphere) and compared them with observed data. From this semi-empirical approach McCracken concluded that M2 (full mixing above 60° latitude and a limited mixing between 40° and 60° latitude) is a reasonable model for Antarctica. Vonmoos et al. (2006) assumed that the production of 10Be recorded in Greenland is related to the entire hemisphere in the stratosphere (i.e, global stratospheric mixing) but is limited to latitudes above 40° latitude in the troposphere (partial tropospheric mixing). This approach uses either semi-empirical or indirect arguments in choosing the unknown degree of mixing.
Recent efforts in employing modern atmospheric 3D circulation models for simulations of 10Be transport and deposition, including realistic air-mass transport and dry-vs-wet deposition (Field et al., 2006; Heikkilä et al., 2008, 2009), look more promising. An example of 10Be deposition computed on the world grid using the NASA GISS model (Field et al., 2006) is shown in Figure 12. Precision of the models allows one to distinguish local effects, e.g., for Greenland (Heikkilä et al., 2008). A simulation performed by combining a detailed 10Be-production model with an air-dynamics model can result in an absolute model relating production and deposition of the radionuclide. We may expect this breakthrough to occur in the near future. The validity and usefulness of this approach has been recently demonstrated by Usoskin et al. (2009a), who directly modeled production (using the CRAC model – Usoskin and Kovaltsov, 2008a) and transport (using the GISS ModelE – Koch et al., 2006) of a short-living beryllium isotope 7Be and showed that such a combined model is able to correctly reproduce both the absolute level and temporal variations of the 7Be concentration measured in near ground air around the globe. Keeping in mind the similarity between production and transport of the two beryllium isotopes, 7Be and 10Be, this serves as support for the advanced modelling of 10Be transport. A similar general agreement between measured and modelled seasonal variability has been recently found for 10Be in an Antarctic ice core (Pedro et al., 2011).
In order to properly account for geomagnetic changes (Section 3.1.2), one needs to know the effective region in which the radionuclide is produced before being stored in the archive analyzed. For instance, if the concentration of 10Be measured in polar ice reflects mainly the isotope’s production in the polar atmosphere (as, e.g., assumed by Usoskin et al., 2003c), no strong geomagnetic signal is expected to be observed, since the geographical poles are mostly related to high geomagnetic latitudes. On the other hand, assuming global mixing of atmospheric 10Be before deposition in polar ice (e.g., Masarik and Beer, 1999), one expects that only changes in the geomagnetic dipole moment affect will the signal. However, because of partial mixing, which can be different in the stratosphere and troposphere, taking into account migration and displacement of the geomagnetic dipole axis may be essential for a reliable reconstruction of solar variability from 10Be data (McCracken, 2004). Therefore, only a full combination of the transport and production models, the latter explicitly including geomagnetic effects estimated from paleomagnetic reconstructions, can adequately account for geomagnetic changes and separate the solar signal. These will form the next generation of physics-based models for the cosmogenic-isotope proxy method. We note that paleomagnetic data should ideally not only provide the dipole moment (VADM or VDM) but should also provide estimates of the geomagnetic axis attitude and displacement of the dipole center (Korte et al., 2011).
An interesting new potential proxy for solar activity (or cosmic ray) variability on the long-term centennial-to-millennial time scale has been proposed recently by Traversi et al. (2012). This is the nitrate content in a polar ice core Talos Dome in Antarctica, which has a favorite location in the sense of snow accumulation and conservation of such volatile specie as nitrate. Nitrate-related species are partly produced in the stratosphere/troposphere as a result of the ionization of the atmospheric air by cosmic rays and, partly, via terrestrial sources (e.g., lightnings) and are subject to air transport (Rozanov et al., 2012). As shown by Traversi et al. (2012), the nitrate concentration/flux measured in the Talos Dome ice core for the Holocene period agrees well with the cosmogenic data of 14C in tree rings and 10Be in both Antarctic and Greenland ice cores, on the time scales from centennia to millennia. Due to the large errors of the ice core dating, 200 – 300 years (Schüpbach et al., 2011), shorter time scales cannot be considered. The level of the nitrate variability is generally consistent with that predicted by theoretical models assuming its production by GCR in the atmosphere (Semeniuk et al., 2011; Rozanov et al., 2012). Thus, the nitrate in an ice core provides a potential to become a new proxy of long-term solar activity, with independent atmospheric fate, which would strengthen the robustness of the reconstructions. However, an independent confirmation of the result and a more detailed model are needed before it can serve as a new quantitative proxy. Note that the mechanism of the nitrate production and transport is not related to the possible nitrate peaks claimed to be caused by strong solar energetic-particle events (see Section 5.3).
Several methods have been developed historically to convert measured cosmogenic-isotope data into a solar activity index, ranging from very simple regressions to physics-based models. A new step in long-term solar-activity reconstruction has been made recently, which is the development of the proxy method in which physics-based models are used, instead of a phenomenological regression, to link SN with cosmogenic-isotope production (Usoskin et al., 2003c, 2007; Solanki et al., 2004; Vonmoos et al., 2006; Muscheler et al., 2007; Steinhilber et al., 2008). Due to recent theoretical developments, it is now possible to construct a chain of physical models to model the entire relationship between solar activity and cosmogenic data.
The physics-based reconstruction of solar activity (in terms of sunspot numbers) from cosmogenic proxy data includes several steps:
- Computation of the isotope’s production rate in the atmosphere from the measured concentration in the archive (Sections 3.2.2 and 3.3.2);
- Computation, considering independently-known secular geomagnetic changes (see Section 3.2.5) and a model of the CR-induced atmospheric cascade, of the GCR spectrum parameter quantified via the modulation potential (Section 3.5.2), some reconstructions being terminated at this point;
- Computation of a heliospheric index, whether of the open solar magnetic flux or of the average HMF intensity at the Earth’s orbit (Section 3.5.2)
- Computation of a solar index (sunspot number series), corresponding to the above-derived heliospheric parameter (Section 3.5.3).
Presently, all these steps can be completed using appropriate models. Some models stop after computations of the modulation potential as its translation into the solar index may include additional uncertainties. Although the uncertainties of the models may be considerable, the models allow a full basic quantitative reconstruction of solar activity in the past. However, much needs to be done, both theoretically and experimentally, to obtain an improved reconstruction.
Mathematical regression is the most apparent and often used (even recently) method of solar-activity reconstruction from proxy data (see, e.g., Stuiver and Quay, 1980; Ogurtsov, 2004). The reconstruction of solar activity is performed in two consecutive steps. First, a phenomenological regression (either linear or nonlinear) is built between a proxy data set and a direct solar-activity index for the available “training” period (e.g., since 1750 for WSN or since 1610 for GSN). Then this regression is extrapolated backwards to evaluate SN from the proxy data. The main shortcoming of the regression method is that it depends on the time resolution and choice of the “training” period. The former is illustrated by Figure 13, which shows the scatter plot of the 10Be concentration vs. GSN for the annual and 11-year smoothed data.One can see that the slope of the 10Be-vs-GSN relation (about -500 g/atom) within individual cycles is significantly different from the slope of the long-term relation (about -100 g/atom), i.e., individual cycles do not lie on the line of the 11-year averaged cycles. Moreover, the slope of the regression for individual 11-year cycles varies essentially depending on the solar activity level. Therefore, a formal regression built using the annual data for 1610 – 1985 yields a much stronger GSN-vs-10Be dependence than for the cycle-averaged data (see Figure 13b), leading to a potentially-erroneous evaluation of the sunspot number from the 10Be proxy data.
It is equally dangerous to evaluate other solar/heliospheric/terrestrial indices from sunspot numbers, by extrapolating an empirical relation obtained for the last few decades back in time. This is because the last decades (after the 1950s), which are well covered by direct observations of solar, terrestrial and heliospheric parameters, correspond to a very high level of solar activity. After a steep rise in activity level between the late 19th and mid 20th centuries, the activity remained at a roughly constant high level, being totally dominated by the 11-year cycle without a long-term trends. Accordingly, all empirical relations built based on data for this period are focused on the 11-year variability and can overlook possible long-term trends (Mursula et al., 2003). This may affect all regression-based reconstructions, whose results cannot be independently (directly or indirectly) tested. In particular, this may be related to solar irradiance reconstructions, which are often based on regression-like models, built and verified using data for the last three solar cycles, when there was no strong trend in solar activity.
As an example let us consider an attempt (Belov et al., 2006) to reconstruct cosmic-ray intensity since 1610 from sunspot numbers using a (nonlinear) regression. The regression between the count rate of a neutron monitor and sunspot numbers, established for the last 30 years, yields an agreement at a 95% confidence level for the period 1976 – 2003 (see Figure 3). Based on that, Belov et al. (2006) extrapolated the regression back in time to produce a reconstruction of cosmic-ray intensity (quantified in NM count rate) to 1560 (see Figure 14). One can see that there is no notable long-term trend in the reconstruction, and the fact that all CR maxima essentially lie at the same level, from the Maunder minimum to modern times, is noteworthy. It would be difficult to dispute such a result if there was no direct test for CR levels in the past. Independent reconstructions based on cosmogenic isotopes or theoretical considerations (e.g., Usoskin et al., 2002b; Scherer et al., 2004; Scherer and Fichtner, 2004) provide clear evidence that cosmic-ray intensity was essentially higher during the Maunder minimum than nowadays. This example shows how easy it is to overlook an essential feature in a reconstruction based on a regression extrapolated far beyond the period it is based on. Fortunately, for this particular case we do have independent information that can prevent us from making big errors. In many other cases, however, such information does not exist (e.g., for total or spectral solar irradiance), and those who make such unverifiable reconstructions should be careful about the validity of their models beyond the range of the established relations.
The modulation potential (see Section 3.1.1) is directly related to cosmogenic isotope production in the atmosphere. It is a parameter describing the spectrum of galactic cosmic rays (see the definition and full description of this index in Usoskin et al., 2005) and is sometimes used as a stand-alone index of solar (or, actually, heliospheric) activity. We note that, provided the isotope production rate is estimated and geomagnetic changes can be properly accounted for, it is straightforward to obtain a time series of the modulation potential, using, e.g., the relation shown in Figure 9. Several reconstructions of modulation potential for the last few centuries are shown in Figure 15. While being quite consistent in the relative changes, they differ in the absolute level and fine details. Reconstructions of solar activity often end at this point, representing solar activity by the modulation potential, as some authors (e.g., Beer et al., 2003; Vonmoos et al., 2006; Muscheler et al., 2007) believe that further steps (see Section 3.5.3) may introduce additional uncertainties. However, since is a heliospheric, rather than solar, index, the same uncertainties remain when using it as an index of solar activity. Moreover, the modulation potential is a model-dependent quantity (see discussion in Section 3.1.1) and therefore does not provide an unambiguous measure of heliospheric activity. In addition, the modulation potential is not a physical index but rather a formal fitting parameter to describe the GCR spectrum near Earth (Usoskin et al., 2005) and, thus, is not a universal solar-activity index.
Modulation of GCR in the heliosphere (see Section 3.1.1) is mostly defined by the turbulent heliospheric magnetic field (HMF), which ultimately originates from the sun and is thus related to solar activity. It has been shown, using a theoretical model of the heliospheric transport of cosmic rays (e.g., Usoskin et al., 2002b), that on the long-term scale (beyond the 11-year solar cycle) the modulation potential is closely related to the open solar magnetic flux Fo, which is a physical quantity describing the solar magnetic variability (e.g., Solanki et al., 2000; Krivova et al., 2007).
Sometimes, instead of the open magnetic flux, the mean HMF intensity at Earth orbit, B, is used as a heliospheric index (Caballero-Lopez and Moraal, 2004; McCracken, 2007; Steinhilber et al., 2010). Note that B is linearly related to Fo assuming constant solar-wind speed, which is valid on long-term scales. An example of HMF reconstruction for the last 600 years is shown in Figure 16. In addition, the count rate of a “pseudo” neutron monitor (i.e., a count rate of a neutron monitor if it was operated in the past) is considered as a solar/heliospheric index (e.g., Beer, 2000; McCracken and Beer, 2007).
The open solar magnetic flux Fo described above is related to the solar surface magnetic phenomena such as sunspots or faculae. Modern physics-based models allow one to calculate the open solar magnetic flux from data of solar observation, in particular sunspots (Solanki et al., 2000, 2002; Krivova et al., 2007). Besides the solar active regions, the model includes ephemeral regions. Although this model is based on physical principals, it contains some unknowns like the decay time of the open flux, which cannot be measured or theoretically calculated and has to be found by means of fitting the model to data. This free parameter has been determined by requiring the model output to reproduce the best available data sets for the last 30 years with the help of a genetic algorithm. Inversion of the model, i.e., the computation of sunspot numbers for given Fo values is formally a straightforward solution of a system of linear differential equations, however, the presence of noise in the real data makes it only possible in a numerical-statistical way (see, e.g., Usoskin et al., 2004, 2007). By inverting this model one can compute the sunspot-number series corresponding to the reconstructed open flux, thus forging the final link in a chain quantitatively connecting solar activity to the measured cosmogenic isotope abundance. A sunspot-number series reconstructed for the Holocene using 14C isotope data is shown in Figure 17. While the definition of the grand minima (Section 4.2) is virtually insensitive to the uncertainties of paleomagnetic data, the definition of grand maxima depends on the paleomagnetic model used (Usoskin et al., 2007). Since the Y00 paleomagnetic model forms an upper bound for the true geomagnetic strength (Section 3.1.2), the corresponding solar-activity reconstructions may underestimate the solar-activity level. Accordingly, the grand maxima defined using the Y00 model are robust and can be regarded as “maximum maximorum” (see Section 4.3).
As very important for the climate research, the variations of the total solar irradiance (TSI) are sometimes reconstructed from the solar proxy data (Steinhilber et al., 2009; Vieira et al., 2011). However, the absolute range of the TSI variability on the centennial-millennial time scales still remains unknown (Schmidt et al., 2012).
Detailed computational models of cosmogenic isotope production in the atmosphere (e.g., Masarik and Beer, 1999) have opened up a new possibility for long-term solar-activity reconstruction (e.g., Beer, 2000). The first quantitative reconstructions of solar activity from cosmogenic proxy appeared in the early 2000s based on 10Be deposited in polar ice (Beer et al., 2003; Usoskin et al., 2003c).
Beer et al. (2003) reconstructed the modulation potential on a multi-millennial timescale using the model computations by Masarik and Beer (1999) and the 10Be data from the GISP2 core in Greenland. This result has been extended, even including the 14C data set, and presently covers the whole Holocene (Vonmoos et al., 2006; Steinhilber et al., 2008, 2010). Usoskin et al. (2003c) presented the reconstruction of sunspot activity over the last millennium, based on 10Be data from both Greenland and Antarctica, using a physics-based model described in detail in Usoskin et al. (2004). This result reproduces the four known grand minima of solar activity – Maunder, Spörer, Wolf and Oort minima (see Section 4.2). Later Solanki et al. (2004) reconstructed 10-year–averaged sunspot numbers from the 14C content in tree rings throughout the Holocene and estimated its uncertainties (see Figure 18). This result was disputed by Muscheler et al. (2005), whose concurrent model, however, rested on an erroneous normalization, as argued in Solanki et al. (2005). The reconstruction of Solanki et al. (2004) has been recently updated by Usoskin et al. (2006a), using a newer paleomagnetic reconstruction by Korte and Constable (2005), and was later slightly revised (Usoskin et al., 2007), considering an updated model of the solar open magnetic flux by Krivova et al. (2007). Reconstruction of the HMF from 10Be data has been performed by Caballero-Lopez and Moraal (2004), using a model of CR modulation in the heliosphere and a 10Be production model by Webber and Higbie (2003). Recently it was revised (McCracken, 2007) to present a detailed reconstruction of HMF intensity since 1428. The most recent reconstruction of the heliospheric modulation potential was done by Steinhilber et al. (2012) using the combined principal component analysis of several data sets.
The obtained results are discussed in Section 4.
Because of the diversity of the methods and results of solar-activity reconstruction, it is vitally important to verify them. Even though a full verification is not possible, there are different means of indirect or partial verification, as discussed below. Several solar-activity reconstructions on the millennium timescale, which differ from each other to some degree and are based on terrestrial cosmogenic isotope data, have been published recently by various groups. Also, they may suffer from systematic effects. Therefore, there is a need for an independent method to verify/calibrate these results in order to provide a reliable quantitative estimate of the level of solar activity in the past, prior to the era of direct observations.
The most direct verification of solar-activity reconstruction is a comparison with the actual GSN sunspot data for the last few centuries. However, regression-based models (see Section 3.5.1) cannot be tested in this way, since it would require a long set of independent direct data outside the “training” interval. It is usual to include all available data into the “training” period to increase the statistics of the regression, which rules out the possibility of testing the model. On the other hand, such a comparison to the actual GSN since 1610 can be regarded as a direct test for a physics-based model since it does not include phenomenological links over the same time interval. The period of the last four centuries is pretty good for testing purposes since it includes the whole range of solar activity levels from the nearly spotless Maunder minimum to the modern period of a very active sun. As an example, a comparison between the observed GSN and the 14C-based (Solanki et al., 2004) and 10Be-based (Usoskin et al., 2003c) reconstructions is shown in Figure 18. The agreement between the actual and reconstructed sunspot numbers is quite good, the correlation coefficient for the 14C-based series is r = 0.93 with the RMS deviation between the two series being six for the period of 1610 – 1900 (Solanki et al., 2004). We want to stress that this reconstruction is fully physics based and does not include any fitting to the whole GSN data series; thus this comparison verifies the approach in both absolute level and relative variations. The agreement between GSN and 10Be-based reconstructions (Usoskin et al., 2003c) is also good (r = 0.78, RMS = 10 for 1700 – 1985). In this case, however, the comparison can only test the relative variation because of the unknown proportionality coefficient between the measured concentration of 10Be and the production rate (Section 3.3.3), which is fitted to match the overall level of the reconstructed solar activity. One can see that the reconstructed sunspot series generally follows the real GSN series, depicting the same main features, namely, the Maunder minimum, the tiny Dalton minimum, a slight decrease of activity around 1900 (sometimes called the modern minimum) as well as a steep rise in the first half of the 20th century. This validates the reliability of the physics-based reconstruction of sunspot numbers. Note, however, that individual 11-year cycles are poorly resolved in these reconstructions.
Models focused on the reconstruction of heliospheric parameters (HMF or the modulation potential ) cannot be verified in this manner since no heliospheric data exists before the middle of the 20th century. Comparison to direct cosmic-ray data after the 1950s (or, with caveats, after the 1930s – McCracken and Beer, 2007) is less conclusive, since the latter are of shorter length and correspond to a period of high solar activity, leading to larger uncertainties during grand minima.
It is important that some (semi)empirical relations forming the basis for the proxy method are established for the recent decades of high solar activity. The end of the Modern grand maximum of activity and the current low level of activity, characterized by the highest ever observed cosmic ray flux as recorded by ground-based neutron monitors, the very low level of the HMF and geomagnetic activity, should help to verify the connections between solar activity, cosmic ray fluxes, geomagnetic activity, the heliospheric magnetic field, and open field. Since some of these connections are somewhat controversial, these extreme conditions should help to quantify them better.
Another more-or-less direct test of solar/heliospheric activity in the past comes from cosmogenic isotopes measured in lunar rock or meteorites. Cosmogenic isotopes, produced in meteoritic or lunar rocks during their exposure to CR in interplanetary space, provide a direct measure of cosmic-ray flux. Uncertainties due to imprecisely known terrestrial processes, including the geomagnetic shielding and redistribution process, are naturally avoided in this case, since the nuclides are directly produced by cosmic rays in the body of the rock, where they remain until they are measured, without any transport or redistribution. The activity of a cosmogenic isotope in meteorite/lunar rock represents an integral of the balance between the isotope’s production and decay, thus representing the time-integrated CR flux over a period determined by the mean life of the radioisotope. The results of different analyses of measurements of cosmogenic isotopes in meteoritic and lunar rocks show that the average GCR flux remained roughly constant – within 10% over the last million years and within a factor of 1.5 for longer periods of up to 109 years (e.g., Vogt et al., 1990; Grieder, 2001).
By means of measuring the abundance of relatively short-lived cosmogenic isotopes in meteorites, which fell through the ages, one can evaluate the variability of the CR flux, since the production of cosmogenic isotopes ceases after the fall of the meteorite. A nearly ideal isotope for studying centurial-scale variability is 44Ti with a half-life of 59.2 ± 0.6 yr (a lifetime of about 85 years). The isotope is produced in nuclear interactions of energetic CR with nuclei of iron and nickel in the body of a meteorite (Bonino et al., 1995; Taricco et al., 2006). Because of its mean life, 44Ti is relatively insensitive to variations in cosmic-ray flux on decade (11-year Schwabe cycle) or shorter timescales, but is very sensitive to the level of CR flux and its variations on a centurial scale. Using a full model of 44Ti production in a stony meteorite (Michel and Neumann, 1998) and data on the measured activity of cosmogenic isotope 44Ti in meteorites, which fell during the past 235 years (Taricco et al., 2006), Usoskin et al. (2006c) tested, in a straightforward manner, several recent reconstructions of heliospheric activity after the Maunder minimum. First, the expected 44Ti activity has been calculated from the reconstructed series of the modulation potential, and then compared with the results of actual measurements (see Figure 19). It has been shown that 44Ti data can distinguish between various reconstructions of past solar activity, allowing unrealistic models to be ruled out. Since the life-time of the 44Ti is much longer than the 11-year cycle, this method does not allow for the reconstruction of solar/heliospheric activity, but it serves as a direct way to test existing reconstructions independently. Most of the reconstructions appear consistent with the measured 44Ti activity in meteorites, including the last decades, thus validating their veracity. The only apparently-inconsistent model is the one by Muscheler et al. (2005), which is based on erroneous normalization (as discussed in Solanki et al., 2005). In particular, the 44Ti data confirms significant secular variations of the solar magnetic flux during the last century (cf. Lockwood et al., 1999; Solanki et al., 2000; Wang et al., 2005).
As an indirect test of the solar-activity reconstruction, one can compare different isotopes. The idea behind this test is that two isotopes, 14C and 10Be, have essentially different terrestrial fates, so that only the production signal, namely, solar modulation of cosmic rays, can be regarded as common in the two series. Processes of transport/deposition are completely different (moreover, the 14C series is obtained as an average of the world-wide–distributed samples). The effect of changing geomagnetic fields is also different (although not completely) for the two isotopes, since radiocarbon is globally mixed, while 10Be is only partly mixed before being stored in an archive. Even comparison between data of the same 10Be isotope, but measured in far-spaced ice cores (e.g., Greenland and Antarctica), may help in separating climatic and extraterrestrial factors, since meteorology in the two opposite polar areas is quite different.
The first thorough consistent comparison between 10Be and 14C records for the last millennium was performed by Bard et al. (1997). They assumed that the measured 10Be concentration in Antarctica is directly related to CR variations. Accordingly, 14C production was considered as proportional to 10Be data. Then, applying a 12-box carbon-cycle model, Bard et al. (1997) computed the expected 14C synthetic record. Finally, these 10Be-based 14C variations were compared with the actual measurements of 14C in tree rings, which depicted a close agreement in the profile of temporal variation (coefficient of linear correlation r = 0.81 with exact phasing). Despite some fine discrepancies, which can indicate periods of climatic influence in either (or both) of the series, that result has clearly proven the dominance of solar modulation of cosmogenic nuclide production variations during the last millennium. This conclusion has been confirmed (e.g., Usoskin et al., 2003c; Muscheler et al., 2007) in the sense that quantitative solar-activity reconstructions, based on 10Be and 14C data series for the last millennium, yield very similar results, which differ only in small details. However, a longer comparison over the entire Holocene timescale suggests that, while centennial variations of solar activity reconstructed from the two isotopes are very close to each other, there might be a discrepancy in the very long-term trend (Vonmoos et al., 2006; Muscheler et al., 2007), whose nature is not clear (climate changes, geomagnetic effects or model uncertainties).
Recently, Usoskin et al. (2009b) studied the dominance of the solar signal in different cosmogenic isotope data on different time scales. They compared the expected 10Be variations computed from 14C-based reconstruction of cosmic ray intensity with the actually measured 10Be abundance at the sites and found that: (1) There is good agreement between the 14C and 10Be data sets, on different timescales and at different locations, confirming the existence of a common solar signal in both isotope data; (2) The 10Be data are driven by the solar signal on timescales from about centennial to millennial time scales; (3) The synchronization is lost on short (< 100 years) timescales, either due to local climate or chronological uncertainties (Delaygue and Bard, 2011) but the solar signal becomes important even at short scales during periods of Grand minima of solar activity, (4) There is an indication of a possible systematic uncertainty in the early Holocene (cf., Vonmoos et al., 2006), likely due to a not-perfectly-stable thermohaline circulation. Overall, both 14C- and 10Be-based records are consistent with each other over a wide range of timescales and time intervals.
Thus, comparison of the results obtained from different sources implies that the variations of cosmogenic nuclides on the long-term scale (centuries to millennia) during the Holocene are primarily defined by the solar modulation of CR.
Most of the earlier solar activity reconstructions are based on single proxy records, either 14C or 10Be. Although they are dominated by the same production signal, viz. solar activity, (see Section 3.7.3), they still contain essential fractions of noise.
A promising first step in the direction of extracting the common solar signal from different proxy records was made recently by Steinhilber et al. (2012). They combined, in a composite reconstruction, different 10Be ice core records from Greenland and Antarctica with the global 14C tree ring record. The composite was made in a mathematical way, using the principal component analysis as a numerical tool. This analysis formally finds the common variability in different series, that is assumed to be the solar signal. The resultant reconstruction is shown in Figure 20 and is the most consistent reconstruction of the cosmic ray (and thus solar magnetic activity) variability over the Holocene, available up-to-date. However, since the used mathematical tool can only work with the relative variability, the reconstruction also yields the relative values rather than absolute values, and it is not available in the terms of sunspot numbers.
A full physics-based multi-proxy composite reconstruction of the solar activity on the millennial time scale is still pending.
In this section, a proxy method of past–solar-activity reconstruction is described in detail.
This method is based on the use of indirect proxies of solar activity, i.e., quantitative parameters, which can be measured now, but represent signatures, stored in natural archives, of the different effects of solar magnetic activity in the past. Such traceable signatures can be related to nuclear or chemical effects caused by cosmic rays in the Earth’s atmosphere, lunar rocks or meteorites. This approach allows one to obtain homogeneous data sets with stable quality and to improve the quality of data when new measurement techniques become available. It provides the only possible regular indicator of solar activity on a very long-term scale.
The most common proxy of solar activity is formed by data of the cosmogenic radionuclides, 10Be and 14C, produced by cosmic rays in the Earth’s atmosphere. After a complicated transport in the atmosphere, these cosmogenic isotopes are stored in natural archives such as polar ice, trees, marine sediments, from where they can now be measured. This process is also affected by changes in the geomagnetic field and the climate.
Radioisotope 14C, measured in independently-dated tree rings, forms a very useful proxy for long-term solar-activity variability. It participates in the complicated carbon cycle, which smoothes out spatial and short-term variability of isotope production. For the Holocene period, with its stable climate, it provides a useful tool for studying solar activity in the past. Existing models allow the quantitative conversion between the measured relative abundance of 14C and the production rate in the atmosphere. The use of radiocarbon for earlier periods, the glacial and deglaciation epochs, is limited by severe climate and ocean ventilation changes. Radiocarbon data cannot be used after the end of the 19th century because of the Suess effect and atmospheric nuclear tests.
Another solar activity proxy is the cosmogenic 10Be isotope measured in stratified polar ice cores. Atmospheric transport of 10Be is relatively straightforward, but its details are as of yet unresolved, leading to the lack of a reliable quantitative model relating the measured isotope concentration in ice to the atmospheric production. Presently, it is common to assume that the production rate is proportional, with an unknown coefficient, to the measured concentration. However, a newly-developed generation of models, which include 3D atmospheric-circulation models, will hopefully solve this problem soon.
Recently a new proxy, nitrate concentration measured in an Antarctic ice core, has been proposed for long-term solar activity reconstructions, but it still needs verification and model support.
Modern physics-based models make it possible to build a chain, which quantitatively connects isotope production rate and sunspot activity, including subsequently the GCR flux quantified via the modulation potential, the heliospheric index, quantified via the open solar magnetic flux or the average HMF intensity at the Earth’s orbit, and finally the sunspot-number series. Presently, all these steps can be made using appropriate models allowing for a full basic quantitative reconstruction of solar activity in the past. The main uncertainties in the solar-activity reconstruction arise from paleo-magnetic models and the overall normalization.
An independent verification of the reconstructions, including direct comparison with sunspot numbers, cosmogenic isotopes in meteorites and the comparison of different models with each other, confirms their veracity in both relative variations and absolute level. It also implies that the variations in cosmogenic nuclides on the long-term scale (centuries to millennia) during the Holocene are primarily defined by the solar modulation of CR.