Issue 
J. Space Weather Space Clim.
Volume 5, 2015



Article Number  A8  
Number of page(s)  18  
DOI  https://doi.org/10.1051/swsc/2015009  
Published online  02 April 2015 
Research Article
Assimilation of realtime riometer measurements into models of 30 MHz polar cap absorption
Space Plasma Environment and Radio Science Group (SPEARS), Department of Physics, Lancaster University, Lancaster
LA1 4YB, UK
^{*} Corresponding author: ncrogers@physics.org
Received:
13
November
2014
Accepted:
9
March
2015
Space weather events may adversely affect high frequency (HF) radio propagation, hence the ability to provide nowcasting and forecasting of HF radio absorption is key for industries that rely on HF communications. This paper presents methods of assimilating 30 MHz radio absorption measurements into two types of ionospheric polar cap absorption (PCA) model to improve their performance as nowcasting tools. Type 1 models calculate absorption as m times the square root of the flux of solar protons above an energy threshold, E_{t}. Measurements from 14 riometers during 94 solar proton events (1995–2010) are assimilated by optimising the day and night values of m by linear regression. Further nonlinear optimisations are demonstrated in which parameters such as E_{t} are also optimised and additional terms characterise local time and seasonal variations. These optimisations reduce RMS errors by up to 36%. Type 2 models incorporate altitude profiles of electron and neutral densities and electron temperatures. Here the scale height of the effective recombination coefficient profile in the Dregion is optimised by regression. Furthermore, two published models of the rigidity cutoff latitude (CL) are assessed by comparison with riometer measurements. A small improvement in performance is observed by introducing a 3h lag in the geomagnetic index K_{p} in the CL models. Assimilating data from a single riometer in the polar cap reduces RMS errors below 1 dB with less than 0.2 dB bias. However, many highlatitude riometers now provide absorption measurements in nearreal time and we demonstrate how these data may be assimilated by fitting a loworder spherical harmonic function to both the measurements and a PCA model with optimised parameters.
Key words: Ionosphere (polar) / Proton precipitation / SEP / Aeronomy / Radio Sciences
© N.C. Rogers and F. Honary, Published by EDP Sciences 2015
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Introduction
The HF (3–30 MHz) radio band is important for longdistance communications and overthehorizon surveillance, particularly in polar regions where ground infrastructure and satellite coverage may be limited. However, during space weather events such as solar proton events (SPEs), the highlatitude ionosphere becomes an efficient absorber of HF radio waves. There is therefore a need to develop accurate HF propagation prediction services that provide maps of highlatitude absorption on a realtime or forecast basis. Such online services would be of particular benefit to airlines operating on polar routes, for example.
Skywave HF propagation is facilitated by refraction in the highaltitude Fregion of the ionosphere where resonant free electrons contribute significantly to the onward propagation of the radio wave. However, the radio wave attenuates in the lower Dregion of the ionosphere (50–90 km) since in this region free electrons collide frequently with neutral atoms and molecules and are either captured to form relatively immobile negative ions or subsequently reradiate with a random incoherent phase. Radio waves are most strongly attenuated during SPEs in which energetic protons precipitate into the Dregion ionosphere causing widespread HF absorption across the polar caps (Bailey 1964; Potemra 1972). The latitudinal extent of the absorption region depends on the particles’ rigidity (momentum per unit charge) and the shape of the deflecting geomagnetic field.
In this study, ionospheric absorption has been recorded by 14 riometers which measure ionospheric absorption of the 30 MHz cosmic noise background (see Sect. 1). Historically, riometer measurements (e.g. Potemra 1972; Sellers et al. 1977) have been used to define parameters of climatological polar cap absorption (PCA) models. However, riometer measurements may now be made available on a realtime basis, allowing measurements to be assimilated into a nowcast model of absorption which could be made available online. The data would enhance existing models such as the DRegion Absorption Prediction (DRAP) model (Sauer & Wilkinson 2008), which currently include realtime proton flux measurements from a nearEarth orbiting satellite and nowcasts of geomagnetic indices.
This paper analyses some approaches to assimilating data from riometers into global absorption models, and assesses the performance for the most practicable scenario in which a single polar cap riometer provides a realtime stream of cosmic noise absorption (CNA) measurements.
Two generic types of PCA model are assessed and riometer measurements are assimilated by optimising parameters of these models by regression. Applying higher weights to more recent measurements improves the optimisation in a nowcast tool for radio absorption.
In Section 1 the riometer measurements and data selections are described. Section 2 describes the two types of PCA model and compares the performance of two published models of the rigidity cutoff latitude. Section 3 describes methods of optimising the parameters in the PCA model based on regression to riometer measurements and quantifies the potential improvements in the associated error and correlation statistics. Section 4 describes how the optimisation may be conducted in real time by weighting the regression according to the age of the measurements. Errors and correlation statistics are quantified for an example in which the model parameters are optimised at a single riometer in the polar cap and applied over a range of latitudes. Finally in Section 5 we describe methods by which data from multiple distributed riometers may be assimilated into an optimised model of CNA over the highlatitude region.
1. Measurements
1.1. Riometer measurements
Absorption measurements were obtained from the Canadian Space Agency’s NORSTAR riometer array (Rostoker et al. 1995; Danskin et al. 2008) which, prior to 2005, formed part of a programme called CANOPUS. Figure 1 provides the locations of these riometers, giving their abbreviated place names. The “Churchill line” of seven riometers labelled in Figure 1 lies close to the 94° W geographic meridian (±2°) from Taloyoak (talo) in the north to Pinawa (pina) in the south.
Fig. 1. Riometers of the NORSTAR array. Contours (red) indicate corrected geomagnetic latitudes for epoch 2000. 
These data were complemented by measurements from the zenithal beam of the Imaging Riometer for Ionospheric Studies (IRIS) (Browne et al. 1995; Honary et al. 2011) operated by Lancaster University at Kilpisjärvi, Finland in conjunction with the Sodankylä Geophysical Observatory (SGO). The data were processed using the MultiInstrument Analysis (MIA) software package (Marple & Honary 2004).
Riometer locations are given in Table 1, together with the invariant latitude at 50 km altitude (mean values determined from the International Geomagnetic Reference Field (IGRF) (Finlay et al. 2010) for the data period covered (1995–2010)). Details for five further widebeam riometers from the SGO array in Finland are also included in Table 1 as these were used in a separate assessment of the SPE of July 2000 presented in Section 5.
The Finland riometer measurements were standardised to 30 MHz (as used in the Canadian array) assuming the frequency dependence (Sauer & Wilkinson 2008)(1)
The exponent of 1.5 is based on J. K. Hargreaves’s analysis of multifrequency absorption measurements by Parthasarathy et al. (1963) and agrees with exponents from 1.4 to 1.6 presented in Figure 10 of Patterson et al. (2001) which were determined from three SPEs in 1989 and 1990. Whilst these measurements are based on daytime events, we follow (Sauer & Wilkinson 2008) in also adopting the 1.5 exponent during nighttime measurements.
With the exception of IRIS, the riometers employed widebeam antennas (60° or more between halfpower points) which measured CNA over a broad range of slant angles through the ionosphere. These measurements were corrected to give the absorption values expected for a narrow zenithal beam following the method of Hargreaves & Detrick (2002) and using archived measurements of the antenna directivity gain.
Median absorption values in 5min periods were recorded for all 94 SPEs in the period 1995–2010. An SPE is defined as a period for which the flux of >10 MeV protons exceeds 10 particle flux units (pfu) (1 pfu = 1 cm^{−2} s^{−1} sr^{−1}).
1.2. Satellite measurements
Proton and Xray flux measurements were taken from the energetic particles sensor (EPS) and Xray sensor (XRS) of the Space Environment Monitor (SEM) subsystem (SSL 1996) on the NASA GOES 8 satellite (up to 17 June 2003) and the GOES 11 satellite (from 19 June 2003). (The GOES8 and GOES11 proton sensors are intercalibrated to better than 10%, based on an analysis of SPEs in the period 1998–2013 by Rodriguez et al. (2014).) SPE commencement times were determined from a catalogue published by the NOAA Space Weather Prediction Center (NOAA 2014). The stated SPE commencement time for 24 April 1999 was corrected from 18:04 to 18:40 UT. The end time of each SPE (not published by NOAA) was determined as the beginning of the first period in which J(>10 MeV) remained below 10 pfu continuously for 2 h, based on 5min flux averages.
1.3. Excluded data
Sudden Commencement Absorption (SCA) is a distinct absorption phenomenon which coincides with a geomagnetic storm sudden commencement (SSC) or sudden impulse (SI) and has been studied extensively by Ritchie & Honary (2009a, 2009b), although individual SCA events have been reported since the 1960s. The absorption events begin simultaneously with the SSC, last for 5–10 min and are confined to the auroral zones on both the dayside and nightside. PCA models do not incorporate the SCA mechanism, so following Neal et al. (2013), periods up to 15 min before and up to 6 h after an SI or SSC were excluded from the data set. The times of SI and SSC were taken from the IAGA databases of SSC 1994–2012 (Ebre 2014), retaining only those events for which the majority of the reporting magnetometer stations qualified the SI/SSC event with confidence factor of two or more (on a scale of 0–3 as defined in the reference above).
Whilst riometers are designed to measure the cosmic background noise, extraneous noise sources such as solar radio bursts and manmade noise may contaminate the measurements and there may be errors in the determination of the quietday curve (QDC) – an average diurnal profile of noise determined by manual inspection of the data in the absence of absorption events. Thus to exclude periods of extraneous noise, any absorption values that were negative (noise exceeding the QDC) or less than 0.1 dB were excluded from the data set, as were periods within 15 min of these times. However, it is possible that this procedure fails to exclude periods of weaker extraneous noise that occur during strong absorption events (i.e. where the measured absorption remains greater than 0.1 dB). Details of the generation of the QDC are described for the NORSTAR riometers in a University of Calgary report (NORSTAR 2014) and for the IRIS riometer at Kilpisjärvi using a variant of the percentile technique of Browne et al. (1995) described by Rodger et al. (2013, p. 7815).
Calibration noise signals (of approx. 2min duration every 4 h) were excised from the NORSTAR riometer measurements, as were 26 short periods for which data was clearly corrupted by artefacts. Furthermore, NORSTAR riometer measurements derived from raw signal strengths below 0.2 V were excluded from the data set due to an apparent insensitivity of the calibration procedure at these very weak signal levels. This condition excluded 3.8% of the complete data set.
2. Modelling HF radiowave absorption in the ionosphere
2.1. Photoionisation effects
One source of ionisation in the Dregion ionosphere is solar irradiation, particularly at Xray and extreme ultraviolet (EUV) wavelengths which penetrate to greater depths in the atmosphere. Thus HF absorption is a function of the solarzenith angle, χ, which varies with the timeofday, latitude and season. Xray and EUV flux also varies over the solar activity cycle (with a period of approximately 11 years) which may be parameterised directly by the 12monthsmoothed 10.7cm solar radio flux (F_{10.7}) or by correlation to the International Sunspot Number, R_{12}. It should be noted, however, that in processing riometer absorption measurements the regular diurnal and longerterm variations are effectively removed since absorption is determined relative to the QDC which is recalculated approximately every 2 weeks to account for changing conditions local to the riometer (snow cover, etc.) and the annual cyclic change in the sidereal diurnal profile of cosmic noise. However, individual intense solar flares may cause spikes in the Xray flux and simultaneous peaks in HF absorption called shortwave fadeouts (SWF) which can affect the whole sunlit hemisphere.
Following the procedure used in the DRAP model, absorption in the sunlit ionosphere due to solar Xray flux is determined using the empirical model of Schumer (2009) who determined that on a vertical incidence ionogram the “Highest Affected Frequency” (HAF) affected by at least 1 dB of absorption was(2)where F (W m^{−2}) is the solar Xray flux in the 0.1–0.8 nm band. This is converted to oneway riometer absorption at 30 MHz using Eq. (1) (3)where the factor of 0.5 assumes that absorption on ionosonde measurements is equally divided between the up and down propagation paths.
For the SPE periods in our data set this is a small correction (a 0.02 dB mean and 0.38 dB maximum was observed), and throughout this paper A_{x} has been added to the predicted absorption from proton precipitation models prior to comparison with riometer measurements. Since QDCs are obtained during quiet ionospheric conditions (in the absence of solar flares), it is assumed that A_{x} has a negligible effect on the QDC level.
2.2. Auroral absorption events
In the region of aurorae (approximately 60–80° magnetic latitude), precipitation of energetic electrons during substorms gives rise to “auroral absorption” (AA) enhancements (Foppiano & Bradley 1985; Kavanagh et al. 2004a). This is typically less than 1 dB at 30 MHz, although peaks of up to 6 dB have been observed (Nielsen & Honary 2000). AA events are less likely to occur during the first day of an SPE and are characterised by more irregular and impulsive variations particularly in the hours around local midnight where electron precipitation is more localised and transient in nature.
The probability and median levels of AA have been mapped and modelled on a climatological basis (e.g. Foppiano & Bradley 1985). However, the timing and geographical extent of individual absorption events cannot be determined reliably without localised and timely measurements from, for example, magnetometer arrays which detect substorm onsets (Beharrell et al. in press), or energetic electron sensors distributed in the Earth’s magnetosphere.
In this paper we have neglected AA due to electron precipitation on the assumption that highenergy solar protons will dominate ionisation of the Dregion during SPEs. This is a reasonable assumption for riometers located inside the polar cap and for all riometers in the first two days of the SPE (before the possible arrival of solarcoronal mass ejecta). Short periods of spikey AA could have been identified subjectively by visual inspection but we have not taken this approach. It is envisaged that a combined PCA and AA model would require altitude profiles of both the proton and electronproduced ionisation rates with a common profile for the effective recombination rate.
2.3. Type 2 (fullprofile) PCA models
A radio wave of frequency ω (rad s^{−1}) propagating along the vertical axis, z, has an electric field component , where n is the refractive index and c is the freespace velocity. The wave attenuates at a rate of(4)where ψ is the imaginary part of n, which may be determined using the Appleton formula (Davies 1990), (5)
Here , where is the plasma frequency, with N_{e }= electron number density (m^{−3}), e = electron charge, ε_{0} = freespace permittivity and m = electron mass. Z = ν/ω, where ν is the electronneutral collision rate (s^{−1}).
In Eq. (5), and throughout this paper geomagnetic field effects (birefringence) are neglected since the ordinary and extraordinary waves are not resolved in our riometer measurements and there is negligible difference in the absorption of these components at 30 MHz.
In the Dregion (<85 km) and lower Eregion, μ^{2} ≫ ψ^{2} and so Eq. (4) is given approximately by (Davies 1990)(6) (7)and the quantity μ ≅ 1. The electronneutral collision rate, ν, depends on the temperature of the electrons and the densities of the various neutral constituents in the plasma. Since ω^{2} ≫ ν^{2} it is reasonable to substitute in Eq. (7) an effective electronneutral collision frequency given by the sum(8)where 〈ν_{e,n}〉 is the mean rate of electron collisions with neutral species n when the electrons have a Maxwellian velocity distribution. These have been determined as(9a) (9b) (9c) (9d) (9e)where N_{s} denotes the density of species s (m^{−3}) and T_{e} is the electron temperature (K). Following Beharrell & Honary (2008) we adopt the mean rate given in Morrison et al. (1997) for molecular nitrogen (9a) and that of Thomas & Nesbet (1975) for atomic oxygen (9b), whilst for all other prevalent species (O_{2}, H and He) (9c–9e) we use the rates determined in Banks (1966). The contributions from atomic species (O, H, He) are negligible at Dregion altitudes but we have included these terms for completeness.
Altitude profiles of N_{s} and T_{e} were determined at each riometer location from the climatological neutral atmosphere model NRLMSISE00 (Picone et al. 2002), and it was assumed that the free electrons are in thermal equilibrium with the neutral gas.
Equation (7) shows that variations in radio wave absorption are linearly dependent on changes in N_{e}, which are determined by the rate of ionisation, q, and the combined rate of electron attachment to neutrals and ions and detachment from negative ions. The latter is expressed as where α_{eff} is the effective recombination coefficient. Examples of the altitude profile α_{eff}(z) during solar proton events are presented in Figure 2 based on measurements from various sources (Vickrey et al. 1982; Gledhill 1986; Hargreaves & Birch 2005). The profile is broadly exponential with distinct scale heights of h_{E }~ 51 km in the Eregion (90–120 km) (Vickrey et al. 1982), and h_{Dd }~ 6.1 km and h_{Dn }~ 4.3 km in the day and night Dregion, respectively (Gledhill 1986). These profiles, which were applied in the PCA model of Patterson et al. (2001), are bestfit values from observations that can range over several orders of magnitude at 50 km altitude. It will be shown in Section 3 how the Dregion scale heights h_{Dd} and h_{Dn} may be used as variable parameters of the model which can be optimised by regression to riometer measurements.
Fig. 2. Altitude profiles of α_{eff} based on regression fits to measurements from three published sources. 
Riometers measure the heightintegrated absorption over a vertical path given by(10)where κ is the specific absorption rate defined in Eq. (7), which is proportional to N_{e}. Multifrequency riometers may be used to determine N_{e}(z) directly (Parthasarathy et al. 1963; Lavergnat & Berthelier 1973; Kero et al. 2014) although these systems are not widely implemented. Instead N_{e}(z) is determined by considering the rate of freeelectron production (the ionisation rate) and the effective recombination rate. Under steadystate conditions (i.e. slowly varying energetic particle flux and negligible plasma transport) these rates are equal,(11) and hence (12)
The ionisation rate profile is given by(13)
where q′(E, z) is the rate of ionisation at altitude, z, as a function of the energy incident at the top of the ionosphere, E, and is determined (following Patterson et al. 2001) from measurements of the proton energy loss rate (eV m^{−1}) at sea level (Janni 1966; Eq. (4) of Patterson et al. 2001), modified by height profiles of atmospheric density relative to that at sea level (using the NRLMSISE00 model) and using an assumption that each generated electronion pair results in the energetic proton losing 1/35 eV (Porter et al. 1976). ∂J(E)/∂E is the differential isotropic proton flux (pfu MeV^{−1}) measured by a nearEarth satellite such as GOES. Following Patterson et al. (2001), the factor of (sr) in Eq. (12) is the effective solid angle of proton flux propagating in the downward direction assuming an isotropic angular distribution at the satellite. The energetic proton pitch angle distribution typically isotropises within a few hours of an SPE commencement (Kouznetsov et al. 2014).
2.4. Type 1 (energy threshold) PCA models
During an SPE, the integral proton energy spectrum has been demonstrated to follow a powerlaw distribution (Potemra 1972 and references therein; Smart & Shea 1985),(14)where E_{t} is a threshold energy, such that(15)
Under this approximation, the absorptionflux relation may be expressed as(16)where(17)using the expression for κ in Eq. (7) with μ ≅ 1 and substituting Eq. (12) for N_{e}(z). Choosing appropriate values of m and E_{t} forms the basis of the Type 1 PCA model.
Potemra (1972) calculated m as a function of γ and E_{t} for a number of daytime PCA events and found the values to be least sensitive to changes in γ when E_{t} = 7 MeV. The lower Dregion (40–80 km altitude) contributes little to absorption at night whilst the higher regions (70–80 km) continue to contribute, with most significant absorption in response to 3–5 MeV protons, rather than 5–20 MeV for daytime (Sellers et al. 1977; Kavanagh et al. 2004b; Hargreaves 2005). Sellers et al. (1977) extended the Potemra (1972) analysis, finding appropriate constants separately for night and daytime PCA events which were measured on a riometer at Thule, Greenland. Specifically these are m_{n} = 0.020, m_{d} = 0.115, E_{tn} = 2.2 MeV and E_{td} = 5.2 MeV where the subscripts _{d} and _{n} refer to day and night ionospheres, respectively. These parameters form the basis of the DRAP (Sauer & Wilkinson 2008) PCA model which takes the form(18)where(19)and(20)are the absorption (dB) due to proton collisional ionization in night and daytime ionospheres, respectively (n = 0.5 provides the squareroot relation in Eq. (16)). The terms E_{tn} and E_{td} in Eqs. (19) and (20) define night and daytime energy thresholds, respectively, and E_{c} represents the rigidity cutoff energy at the riometer location (as discussed in Sect. 2.5).(21)where χ is the solarzenith angle. In (Sauer & Wilkinson 2008), χ_{c} = 90° and the halfwidth of the daynight transition, ∆χ = 10°. Thus the absorption prediction in the region of the solar terminator (χ = 80 to 100°) is a linear interpolation between absorption predictions in the fullydeveloped day and night ionospheres.
Sauer & Wilkinson (2008) presented examples of model performance at Thule during 11 large SPEs from 1992 to 2002, and a further validation at Thule and five other riometer locations was presented by Akmaev et al. (2010), extending the data set to cover two further SPE events in 2005. Akmaev et al. (2010) showed that the model substantially overpredicted absorption for some SPEs, particular examples being the July 2000 “Bastille Day” storm and the April 21, 2002 storm. (This is assuming that the riometer measurements were accurately calibrated at the low signal levels associated with these periods of very strong absorption.)
2.5. Rigidity cutoff latitude models
Geomagnetic shielding prevents solar wind particles with insufficient rigidity from directly accessing the ionosphere at magnetic latitudes below a diffuse cutoff latitude (CL) (Störmer 1955; Smart & Shea 2005). Thus at any given location, protons with energy less than cutoff energy, E_{c}, do not contribute to ionospheric absorption (Smart et al. 1999; Leske et al. 2001; Sauer & Wilkinson 2008; Dmitriev et al. 2010; Nesse Tyssøy et al. 2013).
The CL model of Smart et al. (1999) used in DRAP (Sauer & Wilkinson 2008) was based on particle ray tracing through a Tsyganenko (1989) model of the magnetosphere parameterised by the K_{p} index for K_{p} = 0 to 5, and extended by Boberg et al. (1995) for D_{st} index values from 0 to −500 nT.
Dmitriev et al. (2010) determined the CL from proton flux measurements on the NASA Polar Orbiting Operational Environment Satellites (POES) fitting ellipses to the observed CL, parameterised by magnetic local time (MLT), the dipole tilt angle, and geomagnetic indices D_{st}, K_{p} and (in the case of energetic electrons) the AE index. They found that whilst the MLT variation in CL was minimal for highenergy protons (>36 MeV), for medium energies (2.5–6.9 MeV) it varied by up to 8° depending on MLT and displayed a pronounced duskward shift of the CL oval depending on the D_{st} index.
An illustration of the invariant cutoff latitudes predicted by the Dmitriev CL model for 10 MeV protons is presented in Figure 3 for three values of the K_{p} index over a 24h period of MLT. The MLTindependent values used in the DRAP model (from Smart et al. 1999) are also plotted in Figure 3 for comparison. The poleward excursion of CL in the hours around noon MLT has been independently observed and modelled in POES and METOP2 satellite measurements (Birch et al. 2005; Nesse Tyssøy et al. 2013), occurring particularly after a sudden dayside magnetopause compression which results from increases in solar wind dynamic pressure, density and speed. Leske et al. (2001) observed a <1° shift of the CL towards approximately 22 MLT in 8–16 MeV/nucleon alpha particles detected on the SAMPEX satellite and a peak CL around 11 am magnetic equatorial time was observed in OGO4 satellite measurements for lower energy protons (Fanselow & Stone 1972). The poleward excursion around local noon has been cited as a cause of the “midday recovery” in absorption observed on mid to highlatitude riometers (λ = 60–65°) near the CL boundary (Hargreaves et al. 1993; Ranta et al. 1995; Uljev et al. 1995) although Leinbach (1967) presented evidence to indicate that increased anisotropy of the energetic particle pitch angle distribution could be an alternative cause.
Fig. 3. An example comparison of Smart et al. (1999) and Dmitriev et al. (2010) models of the 10 MeV proton cutoff latitudes for three levels of K_{p}. (Date: 22/4/2000, D_{st} = 0.) 
Figure 4 illustrates how the Dmitriev CL model was applied in predicting HF absorption at the Fort McMurray riometer for the SPE of April 1998. It presents the timevarying proton fluxenergy spectrum interpolated to 0.1 MeV resolution from the corrected differential proton flux measurements of the EPS on GOES 8 assuming a powerlaw spectrum. The Dmitriev model CL and their equivalent Lshells, L_{cutoff}(E), were also calculated in 0.1 MeV steps and in the calculation of absorption the proton flux J(E) was zeroed where L_{cutoff} (E) was less than the L of the riometer. The pronounced diurnal variation of the cutoff portion of the proton flux is observed in Figure 4 at energies below approximately 40 MeV, with the highest flux cut off near noon MLT (20 UT).
Fig. 4. The proton flux spectrum, J(E, t) (interpolated between GOES energy channels at 0.1 MeV resolution) during the SPE of April 1998. The rigidity cutoff region at low energies (the black region) is determined from Dmitriev et al. (2010) at the Fort McMurray riometer. The dashed line represents the Smart et al. (1999) cutoff energies. 
Figure 5 presents an example of the 30 MHz absorption during this SPE, comparing the riometer measurements with the DRAP model (solid line) and DRAP with a substituted Dmitriev CL model (crosses). The DRAP/Dmitriev CL model gives a closer agreement with riometer measurements than the original DRAP/Smart CL model, although the midday recoveries observed near local noon (20 UT) on 21 and 22 April 1998 are not clearly observed in the measurements for this event.
Fig. 5. 30 MHz absorption at Fort McMurray during the April 1998 SPE, plotted with superposed predictions of Sauer & Wilkinson (2008) (S&W), which uses the Smart et al. (1999) CL model, S&W with a 2° equatorward shift of the CL, and S&W with CL defined by Dmitriev et al. (2010). 
A comparison of the error and correlation statistics using DRAP with the two CL models is given in Figure 6 for the seven riometers in the Churchill line and for all 94 SPEs. The rootmeansquared error (RMSE), bias, mean absolute error and the Pearson correlation coefficient are shown for the Smart CL model (solid lines) and the Dmitriev CL model (dashed lines). Whilst riometers at or poleward of Fort Churchill are unaffected by the change in CL model, the substitution of Dmitriev’s CL model improves all the error and correlation statistics for riometers at the lower latitude stations.
Fig. 6. Errors and correlation coefficients of the DRAP (Sauer & Wilkinson 2008) model for riometers in the meridional “Churchill line”. Values are compared for both the DRAP CL model (Smart et al. 1999) (solid lines), and an alternative CL model (Dmitriev et al. 2010) (dashed lines). 
Smart et al. 1999 noted that their 29–64 MeV proton CLs were systematically poleward of those reported by Leske et al. (1997) by around 1.5°. Comparison of DRAP model performance in the Canadian and European sectors by Akmaev et al. (2010) suggested possible errors in the location of the rigidity cutoff at high geomagnetic latitudes using the Smart CL model. Further analysis by Neal et al. (2013) using proton flux measurements from the NASA POES satellites suggested that a 1–2° equatorward shift in the Smart CL would improve the model’s performance.
An example of such a correction is presented as the dashed line in Figure 5 representing predictions of the DRAP model after the CL boundary of Smart et al. (1999) is shifted equatorward by 2° of invariant latitude. This provides a reasonable fit to the measurements in this example. Figure 7 presents DRAP model error and correlation statistics (as for Fig. 6) that result from applying a correction to the Smart CL over the range ±4° for three lowlatitude riometer stations in the NORSTAR array. These tests indicate that an equatorward shift, Δλ, of approximately 1–3° improves model errors and correlation. The optimum CL shift depends on the choice of error statistic and to some extent on the riometer chosen.
Fig. 7. Errors and correlation coefficients of the DRAP (Sauer & Wilkinson 2008) model with a range of corrections to the CL based on all 94 SPEs. (a) Fort McMurray riometer (L = 5.5), (b) Island Lake riometer (L = 5.5), (c) Pinawa riometer (L = 4.3). 
Neal et al. (2013) further recommended the use of a 3h time lag in the 3h k_{p} index used in the Smart CL model. Applying this technique results in marginal (<0.2 dB) improvements in RMSE (not shown) for the riometer stations presented in Figure 7.
3. General optimisation of PCA model parameters
3.1. Type 1 PCA models
Parameters of the two types of PCA model described in Sections 2.3 and 2.4 have been optimised so as to minimise the rootmeansquared error (RMSE) of the absorption estimate relative to riometer measurements for all 94 SPEs in the period 1995–2005. A trustregionreflective algorithm (Coleman & Li 1996) was used to search for the optimum parameters of the Type 1 model (18), with initial parameter estimates set as the fixed parameters of the DRAP model (Sellers et al. 1977; Sauer & Wilkinson 2008) with upper and lower parameter search limits set as in Table 2. (These constrain the model to physically realistic values during the optimisation.)
Initial parameter estimates and trustregion bounds used in the optimisation process.
As an extension of this analysis, the PCA model in Eq. (18) was generalised by (i) introducing a variable shift, Δλ, in the Smart CL model, and (ii) introducing additional terms that characterise Magnetic Local Time (MLT) and seasonal variations, where Eq. (18) is modified to:(22)where S, C and D are scalar coefficients, MLT is in hours, and δ is the solar declination, which varies between 23.44° at the June solstice and −23.44° at the December solstice. C quantifies the noonmidnight MLT variation in absorption predictions whilst S quantifies the dawndusk component and D represents the variation between solstices. However, these three additional terms represent only those MLT and seasonal variations that are not encapsulated in the regularly updated QDC measurements at the riometers and in this study they quantify the seasonal or MLT trends in the differences between the quiet ionospheric conditions before a SPE and the conditions during the SPE.
Table 3 presents the results of optimising parameters by regression to measurements from IRIS and the 13 NORSTAR riometers. Column 2 lists the fixed parameters of the Sauer & Wilkinson (2008) model (18) with the associated RMS error, bias (model – measurement), mean percentage overestimate, mean absolute error and the Pearson correlation coefficient. The third column of Table 3 lists the error and correlation statistics when the two linear scaling parameters m_{n} and m_{d} are optimised. The optimal parameter values are smaller than those in the DRAP model and they reduce the RMS error by 22% from 0.96 dB to 0.76 dB.
Error and correlation statistics for the Sauer & Wilkinson, 2008 PCA model (DRAP) based on 14 riometers and all 94 SPEs (1995–2005) for both original and optimised parameter sets. The Smart et al. (1999) CL model is used (with an optional equatorward shift of Δλ).
A marginal improvement on these figures is obtained by also optimising the night and daytime threshold energies E_{tn} and E_{td} (column 4 of Table 3). However, including the correction, Δλ, to the cutoff latitude in the optimisation (column 5) reduces RMS error more significantly to 0.68 dB (30% less than the DRAP RMSE), yielding a much smaller negative bias and an increase in the correlation coefficient.
The righthand column of Table 3 presents the result of further optimising a very wide selection of parameters. These include the linear scaling parameters S, C and D of the generalised model (22), the exponent n in Eqs. (19) and (20) and the limits of the twilight transition region defined by χ_{c} and Δχ in Eq. (21). Optimising all of these parameters further reduces the RMS error to 0.62 dB (36% less than DRAP) with a negligible positive bias (0.007 dB) and also provides the highest correlation coefficient (0.87). The optimal value of the exponent n is reduced from 0.5 to 0.38, and whilst the energy thresholds, E_{tn} and E_{td} increase, which tends to reduce J(>E_{t}), the scaling parameters m_{n} and m_{d} increase, in part to compensate this change. Figure 8 presents the model predictions using the original DRAP parameters (circles) and using the optimised parameters of the generalised model (crosses) given in the righthand column of Table 3.
Fig. 8. Predictions vs. measurements of 30 MHz zenithal CNA at all 14 riometers using the original DRAP model (circles) and a modified generalised bestfit model (crosses) (righthand column of Table 3). The dashed line indicates equality. 
The optimal values of the zenith angles χ_{c} and Δχ are coincidentally very close to the initial estimates used in DRAP (90° and 10°, respectively). This is unlikely to be artefact of the optimisation procedure since further trials over a wide range of initial estimates and bounds exhibited convergence to the same optimal values. Introducing MLT and seasonal dependences using the parameters S, C and D resulted in only small corrections of up to 0.12 dB in magnitude.
The process of assimilating riometer measurements into an optimised model should ideally give consideration to the likely measurement errors so as to weight the measurements appropriately in the regression. Riometer measurements below approximately 1 dB may be particularly subject to errors associated with inaccuracies in the QDC due to changing environmental conditions and to the effects of extraneous noise and interference in the receiver. Similarly, under conditions of high absorption (>10 dB) there may be greater inaccuracy in measuring the weak cosmic noise signal. To assess this, in Table 4 the effect of zeroweighting (i.e. excluding) measurements outside the range 1–10 dB is presented in the same form as for Table 3. The mean errors for the unoptimised S&W model are higher using this extra filter – as may be expected when a biased distribution is truncated – but the errors after optimisation follow a similar pattern of improvement (reading left to right for each error metric) to that observed for the full range of measurements.
In Table 5 the statistics for columns 1–4 of Table 3 are presented for the case in which the Dmitriev CL model replaced the Smart CL model. The principal effect of this is to increase the correlation coefficient to 0.85, and reduce RMS errors to 0.7 dB (a 29% reduction) and bias to −0.12 dB, values which are similar to the case in which the Smart CL was optimised by introducing an optimal 3.2° equatorward shift (column 5 of Table 3).
Optimal parameters and error statistics as for Table 3 (columns 1–4), substituting the Dmitriev et al. (2010) CL model.
In Table 6, error and correlation statistics are presented for the case in which measurements are recorded at Taloyoak only. At this polar cap location the cutoff rigidity is effectively zero, so optimal parameters did not depend on the rigidity CL model. A similar pattern of improvement in the statistics is observed as for Table 3, with RMS errors reduced from 1.2 dB to less than 0.6 dB (a 54% reduction), and bias reduced from 0.5 dB to −0.1 dB or better, even with the optimisation of just two parameters, m_{n} and m_{d}. The distribution of model predictions vs. Taloyoak riometer measurements is presented in Figure 9 for the original DRAP model parameters (circles) and for the generalised, optimised model (crosses) (using parameters in the righthand column of Table 6). The plot illustrates the improvement in bias and error variance of the model during solar proton events.
Fig. 9. Predictions vs. measurements of 30 MHz zenithal CNA at the Taloyoak riometer using the original DRAP model (circles) and a modified generalised bestfit model (crosses) (righthand column of Table 6). The dashed line indicates equality. 
3.2. Optimisation of Type 2 PCA models
For Type 2 models, a LevenbergMarquardt nonlinear regression algorithm (Seber & Wild 2003) – a “damped leastsquares” method – was implemented, and this successfully converged to a solution in the vast majority of cases. The parameters to be optimised were the scale heights h_{Dd} and h_{Dn} of α_{eff} in the day and nighttime Dregion, respectively. The initial estimates of h_{Dd} and h_{Dn} were those determined by Gledhill (1986). An example of initial and optimised parameters and the associated model error and correlation statistics is presented in Table 7 for optimisation using data for all 94 SPEs recorded at the Taloyoak riometer. Comparison with Table 6 indicates that this fullprofile (Type 2) model underperforms relative to the Type 1 model, but is considerably improved by optimisation of just two parameters, reducing the RMS error by 66% from 2.2 dB to 0.75 dB and the bias to −0.2 dB.
Initial and optimised parameters of a Type 2 PCA model for the Taloyoak riometer (all 94 SPEs).
4. Realtime optimisation of model parameters
The chemical constitution and temperature profiles of the mesosphere vary significantly over time, such that fixed (day or night) recombination coefficient and ionisation rate profiles also vary with respect to the climatological average. α_{eff}(z) has been observed to vary from day to day during SPEs, and gradually throughout the day (e.g. Reagan & Watt 1976) or during the night (e.g. Hargreaves et al. 1993). It is therefore beneficial to optimise PCA model parameters over similar timescales. To accomplish this in the parameter optimisation process, for a given time, t, the set of N riometer measurements recorded at times t_{i} (i = 1 to N) was assigned weights, w_{i}, that diminish exponentially with the age of the measurement, t − t_{i}, thus:(23)where τ is a characteristic decay time. This weighting ensures that optimisation of parameters at times for which there are no recent measurements (where t – t_{i} ≫ τ for all i) will revert to the optimum parameters for uniformly weighted measurements from the complete data set. The examples below include both past (t_{i} < t) and future (t_{i} > t) measurements in the optimisation so the “global” optimum parameters do not change with the date of each SPE.
Examples of optimum values of h_{Dd} and h_{Dn} are presented in Figure 10 during the SPE of April 1998, based on regression to riometer measurements at the polar cap riometer in Taloyoak, Canada, updated every hour. Values are presented for two values of the characteristic decay time, τ, and illustrate the greater effective smoothing of the optimal parameters at τ = 24 h compared with τ = 6 h.
Fig. 10. h_{Dn} and h_{Dd} optimised by regression to Taloyoak measurements and updated at 1h intervals during the SPE of April 1998. Values are shown for two values of the characteristic decay time, τ = 6 h and 24 h. 
Figure 11 presents the associated model predictions using the h_{Dn} and h_{Dd} parameters calculated for τ = 6 h plotted against the Taloyoak riometer measurements. The optimum parameters are recalculated every hour so the minimum age of highly weighted samples used in the optimisation varies between 0 and 1 h. Using these optimised parameters the model provides a reasonably close fit to the measurements.
Fig. 11. Riometer absorption at Taloyoak (“+”) with Type 2 model predictions (solid line) using h_{Dd} and h_{Dn} optimised with τ = 6 h updated at 1h intervals. 
Figure 12 presents the absorption measured for the same SPE at Rankin Inlet, 750 km south of Taloyoak (crosses). The dashed line represents the Type 2 PCA model prediction using the α_{eff}(z) profiles of Gledhill (1986), which overestimates absorption by several dB on the second and third day of the SPE. The positive bias during that period is reduced by substituting the h_{Dd} and h_{Dn} scale heights that were optimised for Taloyoak (Fig. 10) (solid line in Fig. 12), although errors are increased at other times.
Fig. 12. Absorption measured at Rankin Inlet riometer (“+”) with Type 2 model predictions (solid line) using the h_{Dd} and h_{Dn} (from Fig. 10) optimised for the Taloyoak riometer with τ = 6 h, updated at 1h intervals. The dashed line represents predictions using fixed h_{Dd} and h_{Dn} from Gledhill (1986). 
The performance of this approach has been determined statistically by repeating the exercise above for all riometers in the Churchill line of riometers and for 13 of the largest, multiday SPEs. The selection of SPEs listed in Table 8 matches that in the DRAP model evaluation by Akmaev et al. (2010) and represents 35% of the absorption measurements in the full data set.
Thirteen large, multiday SPE events in the period 1995–2005 (following the selection of Akmaev 2010).
The four panels of Figure 13 present (a) RMS errors, (b) error biases (model – measurements), (c) mean ratios (model dB/measurement dB) and (d) Pearson correlation coefficients. Each group of four bars represents one of the seven riometers in the Churchill line from Taloyoak to Pinawa in order of descending latitude. The first (lefthand) bar in each group represents the statistic averaged over all 13 selected SPEs for the Type 2 PCA model using the Gledhill (1986) profiles of α_{eff}, whilst the second bars represent the errors obtained by substituting the h_{Dd} and h_{Dn} scale heights optimised for Taloyoak.
Fig. 13. (a) RMS errors, (b) biases, (c) mean ratios and (d) correlation coefficients for Type 2 and Type 1 PCA model predictions, with original and realtime optimised parameters, averaged over 13 SPEs of Table 8. Optimisations (second and fourth bars in each group) use the Taloyoak riometer measurements only. All models apply the Dmitriev cutoff latitude model. 
The process was repeated for the Sauer & Wilkinson (2008) Type 1 PCA model (third bar in each group), and the fourth (righthand) bars represent the errors of this model after substituting values of coefficients m_{n} and m_{d} optimised for Taloyoak. In all cases a Dmitriev CL model was applied, optimisations were updated every 1h and a 6h characteristic decay time, τ, was used in Eq. (23) to weight the regression.
In all cases except the results for the most distant, lowest latitude riometer (Pinawa) each statistical error metric is improved by applying realtime optimisation based on the Taloyoak riometer observations. RMS errors are reduced to 1 dB or less (except at Pinawa), the magnitude of the bias is reduced to less than 0.2 dB (except at Pinawa and Eskimo Lake) and mean ratios are brought closer to unity (except at Pinawa). The correlation coefficient is little changed, except at Pinawa where correlation is reduced and at Taloyoak where there is an improved correlation (since the models are optimised for this site). These results demonstrate that a simple, singleriometer optimisation can improve PCA model outputs across a wide region of the polar cap. The poor correlation at Pinawa may arise from this station’s location well within the auroral zone which is subject to significant auroral absorption (AA).
5. Spatial interpolation and assimilation of multiple riometer measurements
Optimising PCA models using a single riometer risks extending local measurement errors over the entire area of predictions. For example, an unidentified period of extraneous noise or a sudden change to the antenna gain pattern (for example due to change in snow cover since the last QDC generation) may impair the generalised predictions. Whilst every effort should be made to monitor and detect riometer inaccuracies soon after they occur, the impact of any such failures can be mitigated through the incorporation of multiple distributed riometers. In recent years, an increasing number of riometers have been fitted with the capability of supplying online nearreal time (<15 min latency) measurements, with 23 operational stations in the Canadian region alone (Danskin et al. 2008). Consequently, algorithms have been developed which will assimilate data from multiple stations in realtime into PCA models.
Figure 14a presents a prediction of 30 MHz CNA using the DRAP model at 08:00 UT on 15 July 2000. Absorption measurements from the 13 NORSTAR riometers in Canada and 5 SGO riometers in Finland (see Table 1) are represented as coloured discs. This example is selected for a time near the peak of a major PCA event at which the model overestimated absorption at several riometer sites by several decibels and it is assumed (without validation) that the riometer measurements are accurate. Figure 14b presents the situation after the m_{n} and m_{d} coefficients of the model had been optimised by linear regression to all riometer measurements in a prior 30min period. This technique provides an improved fit to the measurements whilst retaining the spatial structure of the original model. The RMSE (relative to measurements from all riometers) at the time shown is reduced by 51% from 2.87 dB to 1.41 dB and the bias reduces from 1.41 dB to 0.13 dB (see Table 9). Using this technique, any riometer recording an outlier measurement at variance with measurements from neighbouring stations will not heavily influence the model prediction.
Fig. 14. Disc colours indicate absorption measured on riometers from the NORSTAR and SGO arrays on 15/7/2000 08:00 UT. Background colour represents (a) DRAP predictions, (b) DRAP modified with parameters m_{n} = 0.065, m_{d} = 0.022 found by regression to all riometers in the prior 30min period. 
An alternative method is to fit an interpolating function to both measurements and model predictions using a weighted regression. One such function is a truncated series of spherical harmonic basis functions,(24)where are associated Legendre polynomials, θ = colatitude and φ = longitude. The coefficients A_{lm} are found by a weighted linear regression, which allows higher weights to be assigned to the riometer measurements than are assigned to the regularly spaced model outputs.
An example of a fitted function (24) is presented in Figure 15 for L_{max} = 6 using riometer measurements from the NORSTAR and SGO arrays (see Table 1) and a regular grid of DRAP predictions using (a) original DRAP parameters and (b) modified values of m_{d} and m_{n} as determined by the previous method of regression to the riometer measurements. Since there are no riometers at low or southern latitudes, in this example θ is upscaled to a maximum colatitude of 60°. This is because the absorption at this latitude is relatively uniform (and negligible). The DRAP model estimates were incorporated at fixed intervals of 2° latitude and 4° longitude and assigned weights of sin(θ) to correct for the closer gridpoint spacing near the pole. This assimilation method allows spatial variations in absorption to be defined by both the riometer measurements and the model outputs, the contributions from each depending on the weights assigned to the riometers relative to the PCA model outputs. (A weight of 100 was applied for each riometer for the results shown in Fig. 15.) This technique can yield significant improvements to the RMSE and bias of the interpolating function (see Table 9). However, the technique introduces artificial structure into the absorption map arising from the interpolating spherical harmonic function.
Fig. 15. As for Figure 14 but with background colours representing (a) a 6th order sum of spherical harmonics fit by regression to the riometer data and DRAP data points; (b) as for (a) but using optimised DRAP parameters m_{n} and m_{d} as in Figure 14b. 
Conclusions
In recent years there has been a considerable increase in the number of riometers across the arctic region that provide realtime online measurements of 30 MHz radio absorption. We have demonstrated techniques by which these measurements could be assimilated into models of polar cap absorption – models which are currently based on realtime measurements of energetic particle precipitation and estimates of the rigidity cutoff latitude based on geomagnetic indices (or forecasts of their projected values). The assimilation is based on the assumption of slow changes in the altitude profile of the Dregion ionospheric constitution and temperature, which define profiles of ionisation rate and the effective recombination coefficient, which determine the altitude profile of specific absorption. In the case of Type 1 models such as DRAP, riometer measurements can be used to optimise the linear scaling coefficients, m_{n}, m_{d}, for night and day ionospheric profiles, respectively, and the proton energy thresholds E_{tn} and E_{td} may also be optimised by regression. Further parameters may easily be introduced and optimised, such as an offset to the rigidity cutoff latitude and factors based on MLT and season. These optimisations reduce RMS errors by 22–36% based on all riometers, or 54–63% for the polar cap Taloyoak riometer alone. In the case of the Type 2 (full altitude profile) model, optimising the scale height of the effective recombination coefficient in the Dregion ionosphere for night and day time ionospheres provides an effective means of optimisation, reducing RMS errors by 66% at Taloyoak, compared with a fixedparameter climatological model. In this analysis it is assumed that the riometers are perfectly calibrated for all signal levels above 0.2 V. During periods of the most intense absorption the riometer measurements may need to be deweighted in the regression due to inaccuracies associated with measuring very weak cosmic noise. Similarly, the relative accuracy of riometer measurements may be degraded under low absorption conditions due to inaccuracies in the QDC and the presence of extraneous noise. Type 1 model optimisations based exclusively on measurements between 1 and 10 dB from all riometers yielded 22–44% improvements to the RMSE.
By assigning weights to the measurements that decay exponentially with age over many hours, even a single riometer may be used to provide updates over a full range of MLTs. It has been demonstrated that a riometer located in the polar cap could have improved both Type 1 and Type 2 PCA model predictions for 13 major solar proton events between 1995 and 2005, reducing RMS errors to less than 1 dB in most cases, reducing biases to less than 0.2 dB and marginally improving correlation. However, the optimisation of the model using polar cap measurements reduced error performance in the lowest latitude station (Pinawa), suggesting poor correlation in the ionospheric profiles, potentially due to increased auroral absorption due to electron precipitation at this latitude. Assimilation of data from multiple riometers follows a similar approach: PCA model parameters are first optimised using recent measurements, and then, optionally, a loworder spherical harmonic interpolating function is fitted through both the riometer locations and regularly spaced model points which are required to constrain the fit over datasparse areas. The weights assigned to riometer measurements may be increased relative to the model data points to improve the fit to riometer measurements.
The location of the rigidity cutoff is the most important factor in determining absorption at lower latitude stations. Comparison with riometers principally in the Canadian region suggests that the Smart et al. (1999) model used in the DRAP model benefits from a simple correction that shifted the CL boundary 2–3° equatorward. Performance of the Sauer & Wilkinson (2008) model with the CL model substituted by that of Dmitriev et al. (2010) also reduces errors and improves correlation.
Realtime dataassimilative models of HF radio absorption should provide an essential component of online space weather services which are of particular interest to airtraffic controllers and airlines operating services on polar routes. Future research will extend the nowcast PCA models described in this paper to include other ionospheric radio absorption mechanisms such as auroral absorption, and also aim to provide a forecasting capability.
Acknowledgments
This work was supported by the UK Engineering and Physical Sciences Research Council, Grant No. EP/K007971/1. Operational support for the CANOPUS/NORSTAR instrument array was provided by the Canadian Space Agency. We acknowledge the efforts of Don Wallis who was principally responsible for the scientific operation of the CANOPUS riometer array and the NORSTAR team for providing the riometer data used in this study. We also acknowledge the Sodankylä Geophysical Observatory which provided measurements from five riometers in Finland for July 2000, and the US National Geophysical Data Centre for providing GOES satellite data. The editor thanks M. Friedrich and P. Stauning for their assistance in evaluating this paper.
References
 Akmaev, R.A., A. Newman, M. Codrescu, C. Schulz, and E. Nerney. DRAP model validation: I. Scientific report, 2010. Available online: www.ngdc.noaa.gov/stp/drap/DRAPVReport1.pdf. [Google Scholar]
 Bailey, D.K. Polar cap absorption. Planet. Space Sci., 12, 495–541, 1964, DOI: 10.1016/00320633(64)900406. [CrossRef] [Google Scholar]
 Banks, P. Collision frequencies and energy transfer electrons. Planet. Space Sci., 14, 1085–1103, 1966, DOI: 10.1016/00320633(66)900249. [CrossRef] [Google Scholar]
 Beharrell, M., and F. Honary. A new method for deducing the effective collision frequency profile in the Dregion. J. Geophys. Res., 113, A05303, 2008, DOI: 10.1029/2007JA012650. [Google Scholar]
 Beharrell, M.J., F. Honary, C.J. Rodger, and M.A. Clilverd. Substorm induced energetic electron precipitation: morphology and prediction. J. Geophys. Res., in press, DOI: 10.1002/2014JA020632. [Google Scholar]
 Birch, M.J., J.K. Hargreaves, A. Senior, and B.J.I. Bromage. Variations in cutoff latitude during selected solar energetic proton events. J. Geophys. Res., 110, A07221, 2005, DOI: 10.1029/2004JA010833. [Google Scholar]
 Boberg, P.R., A.J. Tylka, J.H. AdamsJr, E.O. Flückiger, and E. Kobel. Geomagnetic transmission of solar energetic protons during the geomagnetic disturbances of October 1989. Geophys. Res. Lett., 22, 1133–1136, 1995, DOI: 10.1029/95GL00948. [CrossRef] [Google Scholar]
 Browne, S., J.K. Hargreaves, and B. Honary. An imaging riometer for ionospheric studies. Electronics & Communication Engineering Journal, 7 (5), 209–217, 1995, DOI: 10.1049/ecej:19950505. [CrossRef] [Google Scholar]
 Coleman, T.F., and Y. Li. An interior, trust region approach for nonlinear minimization subject to bounds. SIAM Journal on Optimization, 6, 418–445, 1996, DOI: 10.1137/0806023. [Google Scholar]
 Danskin, D.W., D. Boteler, E. Donovan, and E. Spanswick. The Canadian riometer array, in Proceedings of the 12th International Ionospheric Effects Symposium, Alexandria,VA, USA, 13–15 May 2008, pp. 80–86, 2008 (Available from www.ntis.gov, PB2008112709). [Google Scholar]
 Davies, K. Ionospheric Radio, IEE Electromagnetic Series, 31, Peter Peregrinus ltd, on behalf of the Institution of Electrical Engineers, London, UK, ISBN: 086341186X, 1990. [Google Scholar]
 Dmitriev, A.V., P.T. Jayachandran, and L.C. Tasi. Elliptical model of cutoff boundaries for the solar energetic particles measured by POES satellites in December 2006. J. Geophys. Res., 115, A12244, 2010, DOI: 10.1029/2010JA015380. [CrossRef] [Google Scholar]
 Ebre Observatory. International service on rapid magnetic variations, 2014. Available online: http://www.obsebre.es/en/rapid. [Google Scholar]
 Fanselow, J.L., and E.C. Stone. Geomagnetic cutoffs for cosmicray protons for seven energy intervals between 1.2 and 39 MeV. J. Geophys. Res., 77 (22), 3999–4009, 1972, DOI: 10.1029/JA077i022p03999. [CrossRef] [Google Scholar]
 Finlay, C.C., S. Maus, C.D. Beggan, T.N. Bondar, A. Chambodut, et al. International Geomagnetic Reference Field: the eleventh generation. Geophys. J. Int., 183 (3), 1216–1230, 2010, DOI: 10.1111/j.1365246X.2010.04804.x. [Google Scholar]
 Foppiano, A.J., and P.A. Bradley. Morphology of background auroral absorption. J. Atmos. Terr. Phys., 47, 663–674, 1985, DOI: 10.1016/00219169(85)901023. [CrossRef] [Google Scholar]
 Gledhill, J.A. The effective recombination coefficient of electrons in the ionosphere between 50 and 150 km. Radio Sci., 21 (3), 399–408, 1986, DOI: 10.1029/RS021i003p00399. [CrossRef] [Google Scholar]
 Hargreaves, J.K. A new method of studying the relation between ionization rates and radiowave absorption in polarcap absorption events. Ann. Geophys., 23, 359–369, 2005, DOI: 10.5194/angeo233592005. [CrossRef] [Google Scholar]
 Hargreaves, J.K., and M.J. Birch. On the relations between proton influx and Dregion electron densities during the polarcap absorption event of 2829 October 2003. Ann. Geophys., 23, 3267–3276, 2005, DOI: 10.5194/angeo2332672005. [CrossRef] [Google Scholar]
 Hargreaves, J.K., and D.L. Detrick. Application of PCA events to the calibration of riometer systems. Radio Sci., 37 (3), 1035, 71 to 711, 2002, DOI: 10.1029/2001RS002465. [CrossRef] [Google Scholar]
 Hargreaves, J.K., A.V. Shirochkov, and A.D. Farmer. The polar cap absorption event of 19–21 March 1990: recombination coefficients, the twilight transition and the midday recovery. J. Atmos. Terr. Phys., 55 (6), 857–862, 1993, DOI: 10.1016/00219169(93)90026U. [CrossRef] [Google Scholar]
 Honary, F., S.R. Marple, K. Barratt, P. Chapman, M. Grill, and E. Nielsen. Digital beamforming imaging riometer systems. Rev. Sci. Instrum., 82 (3), 1–15, 2011, DOI: 10.1063/1.3567309. [CrossRef] [Google Scholar]
 Janni, J.F. Calculations of Energy Loss, Range, Path length, Straggling, Multiple, Scattering, and the Probability of Inelastic Nuclear Collisions for 0.1 to 1000MeV Protons, Air Force Weapons Laboratory, Kirtland Air Force Base, New Mexico, 1966. Available online: http://www.dtic.mil/dtic/tr/fulltext/u2/643837.pdf. [Google Scholar]
 Kavanagh, A.J., M.J. Kosch, F. Honary, A. Senior, S.R. Marple, E.E. Woodfield, and I.W. McCrea. The statistical dependence of auroral absorption on geomagnetic and solar wind parameters. Ann. Geophys., 22 (3), 877–887, 2004a, DOI: 10.5194/angeo228772004. [CrossRef] [Google Scholar]
 Kavanagh, A.J., S.R. Marple, F. Honary, I.W. McCrea, and A. Senior. On solar protons and polar cap absorption: constraints on an empirical relationship. Ann. Geophys., 22, 1133–1147, 2004b, DOI: 10.5194/angeo2211332004. [CrossRef] [Google Scholar]
 Kero, A., J. Vierinen, D. McKayBukowski, C.F. Enell, M. Sinor, L. Roininen, and Y. Ogawa. Ionospheric electron density profiles inverted form a spectral riometer measurement. Geophys. Res. Lett., 41, 5370–5375, 2014, DOI: 10.1002/2014GL060986. [CrossRef] [Google Scholar]
 Kouznetsov, A., D.J. Knudsen, E.F. Donovan, and E. Spanswick. Dynamics of the correlation between polar cap radio absorption and solar energetic proton fluxes in the interplanetary medium. J. Geophys. Res.: [Space Phys.], 119 (3), 1627–1642, 2014, DOI: 10.1002/2013JA019024. [CrossRef] [Google Scholar]
 Lavergnat, J., and J.J. Berthelier. An iterative mathematical technique for deriving electrondensity profiles from multifrequency riometer data. Radio Sci., 8 (7), 641–649, 1973, DOI: 10.1029/RS008i007p00641. [CrossRef] [Google Scholar]
 Leinbach, H.. Midday recoveries of polar cap absorption. J. Geophys. Res., 72 (21), 5473–5483, 1967, DOI: 10.1029/JZ072i021p05473. [CrossRef] [Google Scholar]
 Leske, R.A., R.A. Mewaldt, E.C. Stone, and T.T. Rosenvinge. Geomagnetic cutoff variations during solar energetic particle events – implications for the space station, in Proceedings of the 25th International Cosmic Ray Conference, Durban, South Africa, 30 July – 6 August, 2, 381, 1997. Available online: http://adsabs.harvard.edu/abs/1997ICRC....2..381L. [Google Scholar]
 Leske, R.A., R.A. Mewaldt, E.C. Stone, and T.T. Rosenvinge. Observations of geomagnetic cutoff variations during solar energetic particle events and implications for the radiation environment at the Space Station. J. Geophys. Res.: [Space Phys.], 106 (A12), 30011–30022, 2001, DOI: 10.1029/2000JA000212. [CrossRef] [Google Scholar]
 Marple, S.R., and F. Honary. A multiinstrument data analysis toolbox. Adv. Polar Upper Atmos. Res., 18, 120–130, 2004. Available online: http://polaris.nipr.ac.jp/~uap/apuar/apuar18/PUA1812.pdf. [Google Scholar]
 Morrison, M.A., W. Sun, W.A. Isaacs, and W.K. Trail. Ultrasimple calculation of very low energy momentumtransfer and rotationalexcitation cross sections: e – N_{2} scattering. Phys. Rev. A, 55 (4), 2786–2798, 1997, DOI: 10.1103/PhysRevA.55.2786. [CrossRef] [Google Scholar]
 Neal, J.J., C.J. Rodger, and J.C. Green. Empirical determination of solar proton access to the atmosphere: Impact on polar flight paths. Space Weather, 11, 420–433, 2013, DOI: 10.1002/swe.20066. [CrossRef] [Google Scholar]
 Nesse Tyssøy, H., J. Stadsnes, F. Søraas, and M. Sørbø. Variations in cutoff latitude during the January 2012 solar proton event and implication for the distribution of particle energy deposition. Geophys. Res. Lett., 40, 4149–4153, 2013, DOI: 10.1002/grl.50815. [CrossRef] [Google Scholar]
 Nielsen, E., and F. Honary. Observations of ionospheric flows and particle precipitation following a sudden commencement. Ann. Geophys., 18, 908–917, 2000, DOI: 10.1007/s005850000908y. [CrossRef] [Google Scholar]
 NOAA. Solar Proton Events Affecting the Earth Environment, 2014. Available online: http://www.swpc.noaa.gov/ftpdir/indices/SPE.txt. [Google Scholar]
 NORSTAR. CANOPUS Quiet Day Curve Generation, 2014. Available online: http://aurora.phys.ucalgary.ca/norstar/rio/doc/CANOPUS_Riometer_Baselining.pdf (accessed 18 Sep. 2014). [Google Scholar]
 Parthasarathy, R., G.M. Lerfald, and C.G. Little. Derivation of electron density profiles in the lower ionosphere using radio absorption measurements at multiple frequencies. J. Geophys. Res, 68 (12), 3581–3588, 1963, DOI: 10.1029/JZ068i012p03581. [CrossRef] [Google Scholar]
 Patterson, J.D., T.P. Armstrong, C.M. Laird, D.L. Detrick, and A.T. Weatherwax. Correlation of solar energetic protons and polar cap absorption. J. Geophys. Res., 106 (A1), 109–163, 2001, DOI: 10.1029/2000JA002006. [Google Scholar]
 Picone, J.M., A.E. Hedin, D.P. Drob, and A.C. Aikin. NRLMSISE00 empirical model of the atmosphere: statistical comparisons and scientific issues. J. Geophys. Res., 107 (A12), SIA 151–1516, 2002, DOI: 10.1029/2002JA009430. [Google Scholar]
 Porter, H.S., C.H. Jackman, and A.E.S. Green. Efficiencies for production of atomic nitrogen and oxygen by relativistic proton impact in air. J. Chem. Phys., 65, 154–167, 1976, DOI: 10.1063/1.432812. [NASA ADS] [CrossRef] [Google Scholar]
 Potemra, T.A. The empirical connection of riometer absorption to solar protons during PCA events. Radio Sci., 7 (5), 571–577, 1972, DOI: 10.1029/RS007i005p00571. [CrossRef] [Google Scholar]
 Ranta, H., H. Yamagishi, and P. Stauning. Twilight anomaly, midday recovery and cutoff latitudes during the intense polar cap absorption event of March 1991. Ann. Geophys., 13, 262–276, 1995, DOI: 10.1007/s0058599502621. [CrossRef] [Google Scholar]
 Reagan, J.B., and T.M. Watt. Simultaneous satellite and radar studies of the D region ionosphere during the intense solar particle events of August 1972. J. Geophys. Res., 81 (25), 4579–4596, 1976, DOI: 10.1029/JA081i025p04579. [CrossRef] [Google Scholar]
 Ritchie, S.E., and F. Honary. Storm sudden commencement and its effect on highlatitude HF communication links. Space Weather, 7, S06005, 2009a, DOI: 10.1029/2008SW000461. [CrossRef] [Google Scholar]
 Ritchie, S.E., and F. Honary. Observed characteristics of sudden commencement absorption. J. Atmos. Sol. Terr. Phys., 71, 609–617, 2009b, DOI: 10.1016/j.jastp.2008.11.011. [CrossRef] [Google Scholar]
 Rodger, C.J., A.J. Kavanagh, M.A. Clilverd, and S.R. Marple. Comparison between POES energetic electron precipitation observations and riometer absorptions: Implications for determining true precipitation fluxes. J. Geophys. Res.: [Space Phys.], 118, 7810–7821, 2013, DOI: 10.1002/2013JA019439. [CrossRef] [Google Scholar]
 Rodriguez, J.V., J.C. Krosschell, and J.C. Green. Intercalibration of GOES 8–15 solar proton detectors. Space Weather, 12, 92–109, 2014, DOI: 10.1002/2013SW000996. [CrossRef] [Google Scholar]
 Rostoker, G., J.C. Samson, F. Creutzberg, T.J. Hughes, D.R. McDiarmid, A.G. McNamara, A. Vallance Jones, D.D. Wallis, and L.L. Cogger. CANOPUS – a groundbased instrument array for remote sensing the high latitude ionosphere during the ISTP/GGS program. Space Sci. Rev., 71, 743–760, 1995, DOI: 10.1007/BF00751349. [CrossRef] [Google Scholar]
 Sauer, H.H., and D.C. Wilkinson. Global mapping of ionospheric HF/VHF radio wave absorption due to solar energetic protons. Space Weather, 6, S12002, 2008, DOI: 10.1029/2008SW000399. [CrossRef] [Google Scholar]
 Schumer, E.A. Improved modeling of midlatitude Dregion ionospheric absorption of high frequency radio signals during solar xray flares. Ph.D. dissertation, AFIT/DS/ENP/09J01, United States Air Force, WrightPatterson Air Force Base, Ohio, USA, p. 49, June 2009. [Google Scholar]
 Seber, G.A.F., and C.J. Wild. Nonlinear Regression, Wiley, ISBN: 9780471617600, 2003, DOI: 10.1002/0471725315. [Google Scholar]
 Sellers, B., F.A. Hanser, M.A. Stroscio, and G.K. Yates. The night and day relationships between polar cap riometer absorption and solar protons. Radio Sci., 12, 779–789, 1977, DOI: 10.1029/RS012i005p00779. [CrossRef] [Google Scholar]
 Smart, D.F., M.A. Shea, E.O. Fluckiger, A.J. Tylka, and P.R. Boberg. Changes of calculated vertical cutoff rigidities at the altitude of the International Space Station as a function of geomagnetic activity. Proc. Int. Cosmic Ray Conf., 7, 337–340, 1999. [Google Scholar]
 Smart, D.F., and M.A. Shea. Chapter 6 Galactic Cosmic Radiation and Solar Energetic Particles. In A.S., Jursa, Editor, Air Force Handbook of Geophysics and the Space Environment, Air Force Research Laboratory, USA, 1985. National Technical Information Service, (www.ntis.gov) Document accession number ADA167000. [Google Scholar]
 Smart, D.F., and M.A. Shea. A review of geomagnetic cutoff rigidities for earthorbiting spacecraft. Adv. Space Res., 36 (10), 2012–2020, 2005, DOI: 10.1016/j.asr.2004.09.015. [Google Scholar]
 SSL. GOES IM Data Book, DRL 1010801 Ed. Space Systems Loral, Revision 1, Aug. 31, 1996. Available online: http://rsd.gsfc.nasa.gov/goes/text/goes.databook.html. [Google Scholar]
 Störmer, C. The Polar Aurora, Clarendon Press, Oxford, 1955. [Google Scholar]
 Thomas, L.D., and R.K. Nesbet. Low energy electron scattering by atomic oxygen. Phys. Rev. A, 12 (4), 1729–1730, 1975, DOI: 10.1103/PhysRevA.11.170. [CrossRef] [Google Scholar]
 Tsyganenko, N.A. A magnetospheric magnetic field model with a warped tail current sheet. Planet. Space Sci., 37 (1), 5–20, 1989, DOI: 10.1016/00320633(89)900664. [Google Scholar]
 Uljev, V.A., A.V. Shirochkov, I.V. Moskvin, and J.K. Hargreaves. Midday recovery of the polar cap absorption of 19–21 March 1990: a case study. J. Atmos. Terr. Phys., 57 (8), 905, 1995, DOI: 10.1016/00219169(94)000897. [CrossRef] [Google Scholar]
 Vickrey, J.F., R.R. Vondrak, and S.J. Matthews. Energy deposition by precipitating particles and Joule dissipation in the auroral ionosphere. J. Geophys. Res., 87, 5184, 1982, DOI: 10.1029/JA087iA07p05184. [CrossRef] [Google Scholar]
Cite this article as: Rogers NC & Honary F. Assimilation of realtime riometer measurements into models of 30 MHz polar cap absorption. J. Space Weather Space Clim., 5, A8, 2015, DOI: 10.1051/swsc/2015009.
All Tables
Initial parameter estimates and trustregion bounds used in the optimisation process.
Error and correlation statistics for the Sauer & Wilkinson, 2008 PCA model (DRAP) based on 14 riometers and all 94 SPEs (1995–2005) for both original and optimised parameter sets. The Smart et al. (1999) CL model is used (with an optional equatorward shift of Δλ).
Optimal parameters and error statistics as for Table 3 (columns 1–4), substituting the Dmitriev et al. (2010) CL model.
Initial and optimised parameters of a Type 2 PCA model for the Taloyoak riometer (all 94 SPEs).
Thirteen large, multiday SPE events in the period 1995–2005 (following the selection of Akmaev 2010).
All Figures
Fig. 1. Riometers of the NORSTAR array. Contours (red) indicate corrected geomagnetic latitudes for epoch 2000. 

In the text 
Fig. 2. Altitude profiles of α_{eff} based on regression fits to measurements from three published sources. 

In the text 
Fig. 3. An example comparison of Smart et al. (1999) and Dmitriev et al. (2010) models of the 10 MeV proton cutoff latitudes for three levels of K_{p}. (Date: 22/4/2000, D_{st} = 0.) 

In the text 
Fig. 4. The proton flux spectrum, J(E, t) (interpolated between GOES energy channels at 0.1 MeV resolution) during the SPE of April 1998. The rigidity cutoff region at low energies (the black region) is determined from Dmitriev et al. (2010) at the Fort McMurray riometer. The dashed line represents the Smart et al. (1999) cutoff energies. 

In the text 
Fig. 5. 30 MHz absorption at Fort McMurray during the April 1998 SPE, plotted with superposed predictions of Sauer & Wilkinson (2008) (S&W), which uses the Smart et al. (1999) CL model, S&W with a 2° equatorward shift of the CL, and S&W with CL defined by Dmitriev et al. (2010). 

In the text 
Fig. 6. Errors and correlation coefficients of the DRAP (Sauer & Wilkinson 2008) model for riometers in the meridional “Churchill line”. Values are compared for both the DRAP CL model (Smart et al. 1999) (solid lines), and an alternative CL model (Dmitriev et al. 2010) (dashed lines). 

In the text 
Fig. 7. Errors and correlation coefficients of the DRAP (Sauer & Wilkinson 2008) model with a range of corrections to the CL based on all 94 SPEs. (a) Fort McMurray riometer (L = 5.5), (b) Island Lake riometer (L = 5.5), (c) Pinawa riometer (L = 4.3). 

In the text 
Fig. 8. Predictions vs. measurements of 30 MHz zenithal CNA at all 14 riometers using the original DRAP model (circles) and a modified generalised bestfit model (crosses) (righthand column of Table 3). The dashed line indicates equality. 

In the text 
Fig. 9. Predictions vs. measurements of 30 MHz zenithal CNA at the Taloyoak riometer using the original DRAP model (circles) and a modified generalised bestfit model (crosses) (righthand column of Table 6). The dashed line indicates equality. 

In the text 
Fig. 10. h_{Dn} and h_{Dd} optimised by regression to Taloyoak measurements and updated at 1h intervals during the SPE of April 1998. Values are shown for two values of the characteristic decay time, τ = 6 h and 24 h. 

In the text 
Fig. 11. Riometer absorption at Taloyoak (“+”) with Type 2 model predictions (solid line) using h_{Dd} and h_{Dn} optimised with τ = 6 h updated at 1h intervals. 

In the text 
Fig. 12. Absorption measured at Rankin Inlet riometer (“+”) with Type 2 model predictions (solid line) using the h_{Dd} and h_{Dn} (from Fig. 10) optimised for the Taloyoak riometer with τ = 6 h, updated at 1h intervals. The dashed line represents predictions using fixed h_{Dd} and h_{Dn} from Gledhill (1986). 

In the text 
Fig. 13. (a) RMS errors, (b) biases, (c) mean ratios and (d) correlation coefficients for Type 2 and Type 1 PCA model predictions, with original and realtime optimised parameters, averaged over 13 SPEs of Table 8. Optimisations (second and fourth bars in each group) use the Taloyoak riometer measurements only. All models apply the Dmitriev cutoff latitude model. 

In the text 
Fig. 14. Disc colours indicate absorption measured on riometers from the NORSTAR and SGO arrays on 15/7/2000 08:00 UT. Background colour represents (a) DRAP predictions, (b) DRAP modified with parameters m_{n} = 0.065, m_{d} = 0.022 found by regression to all riometers in the prior 30min period. 

In the text 
Fig. 15. As for Figure 14 but with background colours representing (a) a 6th order sum of spherical harmonics fit by regression to the riometer data and DRAP data points; (b) as for (a) but using optimised DRAP parameters m_{n} and m_{d} as in Figure 14b. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.