Neutron monitor count rate increase as a proxy for dose rate assessment at aviation altitudes during GLEs

Radiation exposure due to cosmic rays, specifically at cruising aviation altitudes, is an important topic in the field of space weather. While the effect of galactic cosmic rays can be easily assessed on the basis of recent models, estimate of the dose rate during strong solar particle events is rather complicated and time consuming. Here we compute the maximum effective dose rates at a typical commercial flight altitude of 35 kft ( 11 000 m above sea level) during ground level enhancement events, where the necessary information, namely derived energy/rigidity spectra of solar energetic particles, is available. The computations are carried out using different reconstructions of the solar proton spectra, available in bibliographic sources, leading to multiple results for some events. The computations were performed employing a recent model for effective dose and/or ambient dose equivalent due to cosmic ray particles. A conservative approach for the computation was assumed. A highly significant correlation between the maximum effective dose rate and peak NM count rate increase during ground level enhancement events is derived. Hence, we propose to use the peak NM count rate increase as a proxy in order to assess the peak effective dose rate at flight altitude during strong solar particle events using the real time records of the worldwide global neutron monitor network.


Introduction
According to the generally accepted definition, space weather concerns dynamic conditions on the Sun and solar wind resulting in changes in the Earth's radiation environment, magnetosphere and ionosphere, which eventually can compromise the performance of spacecraft and ground-based systems and can endanger human health (e.g. Baker, 1998;Lilensten & Bornarel, 2009).
Active processes on the Sun, such as Coronal Mass Ejections (CMEs), solar flares and high-speed solar wind streams may lead to sequence(s) of disturbances in the Earth's magnetosphere and atmosphere, occasionally even impacting the ground level as geomagnetic storms and enhancements of relativistic electron populations in outer radiation belts. Trapped protons and ions can damage satellites. An important part of space weather effects is related to the variable and highly dynamic radiation environment in the near-Earth space and Earth's atmosphere (e.g. Vainio et al., 2009, and references therein). Sporadically, following such solar eruptive processes as solar flares and CMEs Solar Energetic Particles (SEPs) can be produced (e.g. Reames, 1999;Cliver et al., 2004;Reames, 2013;Desai & Giacalone, 2016;Klein & Dalla, 2017, and references therein). SEP events are distinct enhancements of particle fluxes originating from Sun.
In this work we focus on solar energetic protons. Their energy is usually of the order of a few tens of MeV/nucleon, but can occasionally reach about hundred MeV/nucleon or even above a GeV/nucleon. While less energetic SEPs are fully absorbed in the atmosphere, more energetic ones can initiate an atmospheric cascade, similarly to Galactic Cosmic Rays (GCRs), whose secondaries eventually reach the ground, leading to an enhancement of count rates of ground based detectors, in particular Neutron Monitors (NMs). This special class of SEP events is known as Ground Level Enhancements (GLEs). Their occurrence rate is about 10-12 per solar cycle with a higher probability during maximum and decline phase of the solar activity cycle (Shea & Smart, 1990;Stoker, 1995). The high energy SEPs can dramatically change the Earth's radiation environment (Matthiä et al., 2009a;Vainio et al., 2009;Sato et al., 2014). Therefore, strong SEP events form a potential space weather hazard, specifically at aviation altitudes.
The increased intensity of secondary cosmic rays at flight altitudes, specifically during SEP events, is an important space weather issue (e.g. Mewaldt, 2006;Pulkkinen, 2007;Shea & Smart, 2012, and references therein). Aircrews are exposed to an additional complex radiation field, particularly during intercontinental flights over the sub-polar and polar regions, where the magnetospheric shielding is marginal. Recently it was advised to consider the exposure to cosmic radiation of aircrew as occupational (ICRP, 2007). Accordingly, a health monitoring and assessment of the individual accumulated doses of the flight personnel was suggested in the EU (EURATOM, 2013).
The radiation environment, accordingly aircrew exposure depends on geographic position, altitude and solar activity (Spurny et al., 1996(Spurny et al., , 2003Shea & Smart, 2000) and is mainly defined by GCRs, which originate from the Galaxy and consist mostly of protons and a-particles with small abundance of heavier nuclei (e.g. Grieder, 2001;Gaisser & Stanev, 2010, and references therein). The GCR flux, specifically the lowenergy part, is modulated by the solar magnetic field in the Heliosphere inversely following the 11-year solar cycle (e.g. Dorman, 2006;Potgieter, 2013, and references therein). When cosmic ray particles penetrate into the atmosphere, they induce a complicated nuclear-electromagnetic-muon cascade that ionizes the ambient air (Bazilevskaya et al., 2008;Usoskin et al., 2009) and results in enhancement of radiation exposure compared to the sea level (Shea & Smart, 2000;Vainio et al., 2009).
The complex radiation environment at flight altitudes is governed by GCRs. It is determined by different types and energy ranges of the produced secondary particles. On the other hand, high energy SEP events may enhance the radiation exposure in the atmosphere at commercial flight altitudes, specifically in the polar regions (e.g. Shea & Smart, 2000;Spurny et al., 2003). During strong GLE events the radiation exposure is a superposition of the GCRs and SEPs contributions. In addition, the characteristics of SEPs, which determine the radiation environment during GLEs: spectra, angular distribution and duration, differ from event to event (Gopalswamy et al., 2012;Moraal & McCracken, 2012). Therefore, the radiation exposure during GLE events is usually assessed retrospectively using reconstructed spectra from ground based and/or space-borne instruments, individually for each event.
Several space weather issues can be addressed using spaceborne measurements and thorough modelling using improved models of solar activity, interaction of the solar wind with the magnetosphere and models of the radiation belts (e.g. Eastwood et al., 2017;Ganushkina et al., 2017;Koskinen et al., 2017, and references therein). However, assessments of the radiation hazard, specifically during strong SEPs events, are really challenging since it requires detailed information of SEP spectral and anisotropy characteristics. The low energy part of SEP characteristics can be estimated using space-borne particle measurements (e.g. Aguilar et al., 2010;Adriani et al., 2016), accordingly the higher energy part, namely GLE particles, can be studied by NMs (e.g. Dorman, 2004, and references therein). Unfortunately, estimates of the GLE spectral and angular characteristics using measurements from the global NM network require accurate modelling of their propagation through the Earth's atmosphere and magnetosphere, which is time consuming (e.g. Shea & Smart, 1982;Humble et al., 1991;. Besides, possible systematic errors of the retrieved spectral and angular characteristics of SEPs can lead to significant uncertainties of the estimated radiation exposure (Bütikofer & Flückiger, 2013. In this work, we perform a statistical analysis of 34 GLEs out of 72 registered, where the information about the spectral characteristics is available (e.g. Tuohino et al., 2018). For several events there exist different reconstructions of the solar proton spectra. Therefore, for those events we employ various spectra leading to several results for one event. We compute conservatively the maximal effective dose rate during those evens employing a recently proposed model and the derived GLE spectra. As a result we propose a convenient proxy for the maximum effective dose rate due to SEPs at commercial aviation altitudes during GLEs, which is suitable for operational purposes.
The sensitivity of NMs to primary CR is governed by both geomagnetic and atmospheric shielding, leading to the effective cut-off, which is determined by the geomagnetic location and the altitude above sea level. The latter determines the thickness of atmospheric layer above the monitor, since the primary CR must possess a minimum energy necessary to induce an atmospheric shower (%430 MeV/n at the sea level), which can reach the ground (e.g. Grieder, 2001;Dorman, 2004). The geomagnetic cut-off rigidity is marginal in the polar regions, where the atmospheric cut-off dominates the shielding. Therefore, polar NMs, particularly high-altitude ones such as SOPO/SOPB and DOMC/DOMB ( Fig. 1 and Table 1) are more sensitive to primary CR, specifically SEPs, than mid and highcut-off rigidity NMs. In addition, polar NMs possess a better angular resolution, which is important for the GLE analysis (Bieber & Evenson, 1995). An illustration of the asymptotic directions at quiet magnetospheric conditions for several low cut-off rigidity NMs is shown in Figure 2. The asymptotic directions are plotted in the range of maximal response of NMs to GLE particles. A. Mishev et al.: J. Space Weather Space Clim. 2018, 8, A46 For the computations of the NM cut-off rigidity and asymptotic directions we employed the MAGNETOCOSMICS code . We used a combination of the internal geomagnetic model IGRF (epoch 2015) (Macmillan et al., 2003) and the external Tsyganenko-89 model (Tsyganenko, 1989), which offers a good balance between simplicity and precision (Kudela & Usoskin, 2004;Nevalainen et al., 2013).

Model for dose rate computation at flight altitude
In order to estimate the radiation exposure at flight altitudes one needs precise information of the spectral and angular parameters of energetic particles and to posses a precise model for their propagation in the atmosphere. A convenient way to  compute the radiation exposure at a typical flight altitude is based on a Monte Carlo simulation of CR particles propagation and interaction in the Earth's atmosphere (e.g. Ferrari et al., 2001;Roesler et al., 2002). Over the years several models of this kind and/or based on other methods were proposed (Schraube et al., 2000;Ferrari et al., 2001;Roesler et al., 2002;Lewis et al., 2005;Takada et al., 2007;Matthiä et al., 2008;Sato et al., 2008;Latocha et al., 2009;Mertens et al., 2013;Matthiä et al., 2014;Mishev et al., 2014a;Wilson et al., 2014;Copeland, 2017), and a reasonable agreement between several models was achieved (Bottollier-Depois et al., 2009). While the models agree well between each other, it was recently shown that differences in input spectral and angular characteristics of SEPs can lead to significant, up to an order of magnitude, differences of the computed radiation exposure(s) (for details see Bütikofer & Flückiger, 2013. Herein, we employed a recent numerical model for computation of the effective and/or ambient dose equivalent at flight altitudes. The model is based on pre-computed yield functions (see details in . The yield function represents the effective dose produced by a monoenergetic unit flux of primary CR particle, which enters in the atmosphere. It is obtained on the basis of high-statistics Monte Carlo simulations. The effective dose rate at a given atmospheric depth h induced by a primary CR particle is computed by convolution of the yield function with a corresponding primary CR particle spectrum: where J i (T) is the differential energy spectrum of the primary CR arriving at the top of the atmosphere for i-th component (proton or a-particle) and Y i is the effective dose yield function for this type of particles. The integration is over the kinetic energy T above Tcut (Pc), which is defined by the local cut-off rigidity P c for a nuclei of type i the proton's rest mass. Accordingly, the effective dose yield function Yi is defined as: where C j (T*) is the coefficient converting the fluence of secondary particles of type j (neutron, proton, c, e À , e + , l À , l + p À , p + )with energy T* to the effective dose, F i,j (h, T, T*, h, u) is the fluence of secondary particles of type j, produced by a primary particle of type i (proton or a-particle) with a given primary energy T arriving at the top of the atmosphere from zenith angle h and azimuth angle u. The conversion coefficients C j (T*) are considered according to Pelliccioni (2000); Petoussi-Henss et al. (2010). Similar expressions are used for the ambient dose equivalent, the details and look-up tables are given in .
Herein, for computations of the ambient dose equivalent (Fig. 3) or effective dose we employ the force field model of GCR spectrum (Gleeson & Axford, 1968;Caballero-Lopez & Moraal, 2004;Usoskin et al., 2005), where the solar modulation parameter was considered according to Usoskin et al. (2011, Fig. 2. Asymptotic directions of polar NMs computed for the epoch 2015 at quiet magnetospheriC conditions. The abbreviations are given in Table 1. The color lines depict asymptotic directions plotted in the rigidity range 1-5 GV, for DOMC and SOPO from 0.7 to 5 GV respectively. A. Mishev et al.: J. Space Weather Space Clim. 2018, 8, A46 2017). We consider a realistic mass composition of GCRs with the nucleonic ratio of heavier particles including a-particles to protons in the interstellar medium as 0.3 similarly to (Mishev & Velinov, 2011;Kovaltsov et al., 2012). The local interstellar spectrum was taken according to Burger et al. (2000).
Applied for computation of the radiation exposure due to GCR, the model demonstrated a very good agreement with the reference data  and other model (Mertens et al., 2013), the details are given elsewhere by . In this study, the model is compared with recent measurements and two other widely used models at two different locations and three altitudes (Fig. 3) (Schraube et al., 2000;Matthiä et al., 2014;Meier et al., 2016). One can see the good agreement, specifically between Oulu and PANDOCA models, with recent experimental data. Note, that the EPCARD model provides effective dose rate, while the other models are applied here for ambient dose equivalent computation for the given altitude and location. However, a good agreement is observed, which is consistent with (Mertens et al., 2013; results.

Assessment of effective dose at flight altitude
Assessment of the effective dose at a flight altitude during GLE events is challenging, since events possess different features, viz. energy spectrum, duration and time evolution (Gopalswamy et al., 2012;Moraal & McCracken, 2012) and also because of the large diversity of secondaries at the flight altitude (e.g. Spurny et al., 1996). Therefore it is necessary to study each GLE event individually. For correct computations, it is necessary to derive precise spectral and angular characteristics of GLE particles in order to compute the effective dose rate (Eq. (1)). In general, this is possible using the NM data, but it requires such time consuming operations as computations of asymptotic cones, modelling of the NM response and fulfilling an optimization (e.g. Shea & Smart, 1982;Humble et al., 1991;Bombardieri et al., 2006;Vashenyuk et al., 2006;Mishev & Usoskin, 2016;. Moreover, different sets of the derived spectra result in considerably different assessments of the dose rate (Bütikofer & Flückiger, 2013. Here we computed the effective dose rate at the altitude of 35 kft (%11 km a.s.l.) during GLE events, using rigidity spectra derived from the NM data Lockwood et al., 1990a;Humble et al., 1991;Smart et al., 1993;Cramp et al., 1997b;Lovell et al., 1998;Deeley et al., 2002;Bombardieri et al., 2006Bombardieri et al., , 2007Bombardieri et al., , 2008Bütikofer et al., 2009;Matthiä et al., 2009b;Vashenyuk et al., 2011;Bieber et al., 2013;Mishev et al., 2014bPlainaki et al., 2014;Kravtsova & Sdobnov, 2016;Mishev & Usoskin, 2016;Kocharov et al., 2017;.
We conservatively assumed an isotropic distribution of the GLE particles similarly to Copeland et al. (2008) in order to assess the maximum radiation exposure. An important anisotropy is observed in most of GLE events, particularly during The measurements with HAWK environmental radiation monitoring system are adopted from (Meier et al., 2016), the estimations with PANDOCA model are interpolated from (Matthiä et al., 2014). Note that the EPCARD results are in effective dose.
A. Mishev et al.: J. Space Weather Space Clim. 2018, 8, A46 the event onset and initial phase (e.g. Bombardieri et al., 2008;Bütikofer et al., 2009). The anisotropy reveals nonsymmetric solar proton flux over the globe. Therefore, the anisotropy effects are important specifically during the event onset and lead to non-uniform distribution of the exposure (e.g. Bütikofer et al., 2008;Matthiä et al., 2009a, b;Velinov et al., 2013;Mishev & Velinov, 2015. The explicit consideration of the anisotropy would result on underestimation of the effective dose rate compared to an isotropic distribution. In addition, the effective dose rate rapidly decreases at regions with higher cut-off rigidity because of considerably softer spectrum of SEPs compared to GCRs. Therefore, in order to provide a conservative approach, we performed all the computations in a region with the geomagnetic cut-off rigidity P c < 1 GV, where the expected exposure is maximal. Hence, we computed the maximum effective dose rate dueto GLE particles, as shown in Table 2. In cases when several spectral reconstructions of GLE particles are available for a given event, resulting in considerably different effective dose estimates (more than 30-40%), we presented an effective dose range (e.g., for GLEs 30,42,45,69). In this study we considered only the maximum radiation exposure in order to provide a conservative approach. Therefore, mainly the prompt hard component of the GLE particles was taken into account, since it results in the maximum radiation exposure, while the soft component will be studied in a separate work.
One can see that the proposed fit (Eq. (3)) encompasses, within the 95% confidence level, most of the data points.  NM count rate increase can be used as a convenient proxy for an assessment of the maximum effective dose at a flight altitude during a GLE event due to high energy SEPs.
A comparison of computed effective dose rates during several GLEs using Oulu  and PANDOCA (Matthiä et al., 2014) models is performed (see for details Table 4). Note, that here we assume the same GLE spectra (Matthiä et al., 2009a, b) in order to minimize the differences and avoid discrepancy due to different spectral reconstruction (e.g. Bütikofer & Flückiger, 2013. One can see that an agreement within about 10% is achieved, the small difference is most likely due to anisotropy effects, GCR parametrization and/or model features as discussed in (Bottollier-Depois et al., 2009).
About 65% of the maximum NM increases are registered by low rigidity high-altitude stations (Table 2). However, several of those stations are not active (MTWS, SLPM and VSTK). Thus, only SOPO/B and DOMC/B, both stations located in South hemisphere, are low rigidity high-altitude NMs (Fig. 1). In order to provide convenient and effective proxy suitable for operational purposes and to avoid usage of a single NM, which may lead to a bias in a case of highly anisotropic events (Fig. 2), we performed a similar study but for sea level NMs (Table 2). We derived similar correlation as in the previous case, with the same quality. However, considering only sea level NMs, the fit of the distribution leads to more conservative assessment of the effective dose during GLEs (see below). Considering only the sea level NMs, the assessed effective dose is about 50% greater than the previous case i.e. a more conservative approach is proposed. In addition, it is more convenient for operational purposes, due to the large number of sea level uniformly distributed stations.
The distribution shown in Fig. 4 was fitted by the equation: where E(x) is the maximum effective dose rate at the altitude of 35 kft a.s.l. due to SEPs, x corresponds to peak NM count increase in % during the event, a and b are the fitted parameters, the details are given in Table 5. Note, that the GCR contribution to the exposure is routinely computed with model(s), explicitly considering the solar activity, and shall be superposed with Equation (3) calculations.

Conclusions
In this work a number of GLE events have been studied, where the necessary information was available. On the basis of the SEP rigidity spectra derived from NM records, the maximum effective dose rates at the aviation altitude of 35 kft during these events was calculated using a recently proposed model. A highly significant correlation between the maximum effective dose rate due to solar protons and the peak NM count rate increase was found. We propose to use the NM count rate increase as a proxy to assess the effective dose at a flight altitude. Thismakes it possible to obtain a quick estimate of the effective dose rate due to SEPs during GLE events on the basis of records from the global NM network, specifically the low cut-off rigidity stations (see Table 1 and Fig. 2).
Note, that in most cases, the GLE effective dose rate profile, consists of rapid rising at the event onset at initial phase, followed by a decay. Therefore, the peak effective dose rate is due mostly to the hard -prompt component. In addition, the most energetic solar protons from the prompt component arrive in the vicinity of Earth before the bulk of SEPs. Therefore, approximations of the dose rates for operational purposes can be achieved instantly (cf. Latocha et al., 2009). Hence, the global NM network can be used to assess an important space weather effect, namely the radiation exposure of aircrew due to high energy particles of solar origin.