GLE definition
The Ground Level Enhancements (GLE) events are caused by the relativistic solar cosmic ray effect on the ground based detectors (mainly neutron monitors). The solar cosmic rays consist basically of protons. And the particles responsible for the GLE are the relativistic solar protons (RSP). Their energy usually is of order of 1-10 GeV. In rare occasions it can exceed 20 GeV.

Our real-time GLE modeling technique
The parameters of RSP can be determined from the data of the ground based neutron monitor network, with the help of a GLE modeling technique. Such a modeling includes the next steps:

1. Definition of asymptotic viewing cones of the NM stations under study by the particle trajectory computations in a model magnetosphere.
2. Calculation of the NM responses at variable primary solar proton flux parameters.
3. Application of a least square procedure for determining primary solar proton parameters (namely, energy spectrum, anisotropy axis direction, pitch-angle distribution) outside the magnetosphere by comparison of computed ground based detector responses with observations.

The GLE modeling technique was presented in [1]. The improved version of [1] was presented in [2], which took into account the contribution in the NM response not only vertical but oblique incident particles. The account of oblique incident particles is made in our version of the GLE modeling technique [3] too. It is necessary to note, that for real-time of definition of RSP parameters it is required as much as possible simplified technique, which, nevertheless, would give exact enough results. Taking into account such limiting factors, as the limited number of available stations at the moment (~20 versus 35-40 in case of [1-3]), limiting time of calculations, etc., we suggest a simplifying version of our modeling technique [3] which can be used for real-time conditions. First of all, we are limiting computations with accounts only vertically incident particles. And the asymptotic directions of viewing for them are calculated in real time with the Tsyganenko 1989 magnetosphere model with the Kp index of a current geomagnetic activity as a parameter. The time of one set of 25 NMBD stations asymptotic cones (AC) computations takes about 3 min. The file containing asymptotic directions of viewing for available NMDB stations is forming as one of the input files of the model.

The task to determine parameters of the RSP flux outside the Earth's magnetosphere using data of the ground based network of NMs is an 'inverse problem'. The corresponding 'direct problem' (to determine the GLE increase on a given NM provided primary RSP flux parameters are known) can be expressed by the relation:

  Direct problem, integral (1)

where ΔN/N is a relative increase of the NM count rate; R is a particles' rigidity; R C is a cutoff rigidity; J(R) is a rigidity spectrum; F(Θ) is a pitch-angle distribution. Pitch-angle Θ depends on R because given NM accepts those particles with given R which had specific Θ outside the magnetosphere. Finally, S(R), Specific Yield Function (SYF, taken from [4]) is a function determining a NM response to RSP flux at a given R.

For computer calculations we substitute the formulum (1) with its discrete analogue:

  Direct problem, discrete (2)

Here a factor A(R) is introduced which is a discrete function A(R)=1 for allowed trajectory (proton with such rigidity can reach the station) and A(R)=0 for forbidden trajectory. Parameter A is determined at the asymptotic cone calculations.

In our model we set:

RC=1 GV (atmospheric cutoff rigidity), RUP=10 GV is the arbitrary upper limit of solar cosmic ray spectrum.
J(R) = J0R-γ* is a modified power rigidity spectrum [2] with variable slope γ* = γ + Δγ ∙ (R-1); γ is a power-law spectral exponent at R = 1 GV, Δγ is a rate of γ increase per 1 GV.
F(θ(R)) = exp(-θ2/2σ2) is a Gaussian pitch-angle distribution. Anisotropy axis is defined by two angles φ and λ which are GSE latitude and longitude.

So the model has 6 parameters to be determined:
J0 intensity of the flux, (m2 s ster GV)-1
γ power-law spectral exponent at R = 1 GV
Δγ rate of γ increase per 1 GV
σ Gaussian parameter of the pitch-angle distribution, deg
φ, λGSE latitude and longitude of anisotropy axis, deg
The inverse problem solution (optimization)
To determine the above parameters we have to minimize the function
  Reverse problem (3)

where P is vector of parameters {J0, γ, Δγ, σ, φ, λ}; (ΔN/N)i mod is a model relative increase at a given NM station i; (ΔN/N)iobs is an observed increase at the station. The summation is taken over all NM stations whose data are available.
The optimization quality may be characterized by the value Discrepancy
As our practice shows good convergence of optimization process is at DΣ ≤ 5% [4,5].

NM data preparation for the model
1. After receiving the GLE alert signal the 1-min count rate data together with barometric pressure from all the NMDB stations start to enter the specified file. The hourly count rate and pressure data of the same stations for the hour previous to one when GLE has started are also fill in the file.
2. The recorded 1-min count rate data are corrected for the barometric effect by the two attenuation method [5,6]. The pre GLE hourly count rate data are corrected for pressure with the one barometric coefficient.
3. For the majority of events operating with 1- min data result in the large errors (unless an increase is of hundreds percents). The 1-min data are summed in 3 or 5-minute meanings, which then are used in the solving of the least squares problem.
The code of solving the least square problem is created in the Mathcad 2001. The code use as input parameters the corrected for pressure data of NMDB stations and calculated asymptotic directions for these stations.

Express and full-scale optimization
For definition of RSP flux parameters in real time the cut down version of a complete technique was developed. Restrictions, imposed on accounts in real time are:
a) less number of NM stations accessible in real time: up to 25 stations NMDB in comparison with 30-35 in a complete technique;
b) limited time for asymptotic cones computations. Thus trajectories of vertically incident particles are computed only. In a complete technique the contribution of oblique incident on NM particles was taken into account;
c) larger step on rigidity (ΔR): 0.01 GV against 0.001 GV in a complete technique.
The reduction of input parameters amount also promoted reduction of operation time of the program solving a least square problem.

[1] M.A. Shea and D.F. Smart, Space Sci.Rev., V.32, P.251-271,1982.
[2] J.L. Cramp, M.L. Duldig, E.O. Flueckiger, J.E. Humble, M.A.Shea, and D.F. Smart, J.Geophys.Res., V.102, P.4237-24248, 1997.
[3] E.V. Vashenyuk, Yu.V. Balabin, B.B. Gvozdevsky, Relativistic solar cosmic ray dynamics in large ground level events, Proc. 21ECRS, Institute Exper. Physics, Kosice, Slovakia, 264, 2009.
[4] H. Debrunner, E. Flueckiger, and J.A. Lockwood. Response of Neutron Monitors to Solar Cosmic Ray Events: 8th European Cosmic Ray Symposium, Rome, 1984, Book of abstracts.
[5] K.G. McCracken, The cosmic ray flare effect 1. Some new methods and analysis, J.Geophys.Res. V.67, P.423-434, 1962.
[6] N.S. Kaminer. About accounts of the barometrical coefficient during cosmic ray flares, Geomagnetism and Aeronomia, V.7, P.806-809, 1967.