Scaling behaviour of lateral distribution of electrons in EAS

A number of recent studies indicate that the average shape of several distributions of electrons in very high energy cosmic ray extensive air showers (EAS), such as the energy distribution or angular distribution, primarily exhibits the so called universality [1–4]: it depends only on the stage of the longitudinal shower development in the atmosphere or equivalently on the longitudinal shower age parameters (s_||) (that represents the variation of the total number of EAS electrons with the atmospheric depth and hence describes the longitudinal shower development; details about the parameters are given in the next section) irrespective of the nature of the primary particle and energy. Such a feature was first divulged from the early work by Kamata and Nishimura [5] in the context of the development of cosmic ray cascades in the atmosphere. The experimental data also appear to substantiate this universality behaviour on an average basis [6]. The universality property is quite advantageous for the analysis of high energy cosmic ray data as it helps to parameterize the electron positron distributions, it seems useful for an accurate estimation of the muon and electromagnetic (EM) contents in an EAS [7] and also it assists to infer the primary mass composition and the nature of the first few interactions from the observed EAS data [6].

The observed lateral density distribution (LDD) of EAS electrons is, however, usually described in terms of the lateral shower age (s_⊥), which essentially describes its slope. Theoretically, the relation s_|| = s_⊥ holds for both EM showers and hadron initiated EAS [2, 5]. In most experiments, however, the estimated s_⊥ differs from s_|| for an EAS with hadrons as primary. Here note that s_|| can be estimated observationally only if the EAS experiment is equipped with Cherenkov or fluorescence detectors, whereas s_⊥ follows immediately from the lateral distribution of electrons, which is a basic measurement of any conventional EAS array consisting of particle detectors. Hence it is imperative to explore the universality of LDD of EAS electrons in terms of the lateral shower age.

A major challenge, however, is the reliable and unambiguous estimate of s_⊥ from experimentally measured electron densities. Usually the LDD of electrons in an EAS is approximated by the well known Nishimura–Kamata–Greisen (NKG) structure function [8] and the shower parameters, namely the shower size (N_e) (the total number of electrons in an EAS at an observational level) and s_⊥ are evaluated by fitting the measured densities with the NKG function. However, experimentally it is observed that the NKG function with a single s_⊥ is insufficient to describe the LDD of EAS electrons properly at all distances, which implies that the lateral age changes with the radial distance. Subsequently some modifications of the NKG structure function [9] were proposed but the radial dependency on the shower age could not be removed totally. Under the circumstances, the notion of local shower age parameters (LAP) was introduced [10] which is in essence the lateral age at a point.

In this work, from a detailed Monte Carlo (MC) simulation study, we show that the shape of the radial variation of the LAP (and hence the LDD of the electrons in an EAS) exhibits some sort of scaling (energy and mass independent) behaviour. The shape of the radial variation of the LAP is, however, found to depend on the choice of the effective Moliére radius in the NKG function. Such a scaling feature provides a better description of the radial electron distributions in EAS and should help to estimate the electron content in an EAS accurately. We also investigate the characteristics of the LAP and its correlation with other EAS observables and thereby the possible role that the parameters may play in a multi-parameter approach to studying EAS, in order to understand the nature of shower initiating particles.

The plan of this paper is the following. In the next section the shower age parameters of EAS will be introduced. In section 3 the simulation procedure adapted in this work will be described. The estimation of the shower age parameters is given in section 4. The results on the characteristics of the local shower age and its correlation with other EAS observables will be given in section 5. Finally we will conclude in Section 6.

The average longitudinal profile of an EM cascade, which is developed in a medium through a multiplicative process involving the interactions of electrons and photons when passing through it, was provided by Greisen [8] in the so-called Approximation B (i.e., considering the processes like a pair production by the photons, bremsstrahlung by the electrons and the ionization loss suffered by the electrons while neglecting the Compton scattering)

$$\begin{equation} N_{e} = \frac{0.31}{\sqrt{\ln (E_{0}/\epsilon _{0})}} \exp [t(1-1.5 \ln (s_{||})]. \end{equation} \tag{ 1 }$$

Here E₀ is the energy of the primary photon generating the cascade, ₀ is the critical energy (below which ionization losses predominate over that due to pair production) with a value of ∼82 MeV. Note that t is expressed here in cascade units (the atmospheric depth has been divided by the electron radiation length in air, taken as 37.1 gmcm⁻²). The longitudinal age s_|| is defined by the simple relation

$$\begin{equation} s_{||}=\frac{3 t}{t+2\ln (E_{0}/\epsilon _{0})} \;. \end{equation} \tag{ 2 }$$

The developmental stage of a pure EM cascade is thus characterized solely by s_||. Basically this parameter represents the change of the total number of EAS electrons with the atmospheric depth (dN/dt) [11]. Recent MC simulation studies [3, 4] exhibited a possible universality for large EAS with primaries in the EeV energy range in terms of s_||.

Advancing the preliminary work of Moliére and Bethe on the LDD of an electron near the cascade maximum, Nishimura and Kamata obtained an expression of radial distance dependence of s_|| by solving the three dimensional diffusion equations in Approximation B and taking into account the hegemony of the multiple Coulomb scattering [5]

$$\begin{equation} s_{||}(r)=\frac{3 t}{t+2\ln (E/\epsilon _{0})+ 2\ln (r/r_{m})} \end{equation} \tag{ 3 }$$

where r is the radial distance from the shower core and r_m is the Moliére radius (near 80 m at sea level) which is a characteristic constant of a medium defined as the radius of a cylinder containing on average 90% of the EMS's energy deposition.

In the same theoretical context, the LDD of cascade particles can be approximated by the NKG structure function [8], given by

$$\begin{equation} f(r)=C(s_{\bot })(r/r_{m})^{s_{\bot }-2}(1+r/r_{m})^{s_{\bot }-4.5} \;, \end{equation} \tag{ 4 }$$

where the normalization factor C(s_⊥) is given by

$$\begin{equation} C(s_{\bot }) = \frac{\Gamma (4.5-s_{\bot })}{2\pi \Gamma (s_{\bot }) \Gamma (4.5-2s_{\bot })} \;. \end{equation} \tag{ 5 }$$

The NKG formula has the advantage of normalization as it is integrable in the Euler Beta function provided s_⊥ is independent of r. The normalization of f(r) implies that ρ(r) = N_ef(r), ρ(r) being the electron density at r. The relation s_|| = s_⊥ was considered to hold for pure EMS [5]. The equations (1), (2) and (4) together provide an attractive and a complete procedure for calculating the 3D development of the EM cascade, as first pointed out by Cocconi [12].

The superposition of many such pure EM cascades build the electron component of a hadron initiated EAS. It was suggested [2, 5] that for hadron initiated EAS, both the longitudinal structure and lateral structure of soft components can be described by that of a resulting single cascade, assigning a suitable value to the age parameter.

The various approximations made in obtaining the solutions, as well as due to an over-simplification of the adopted 3D transport equations, the analytical expressions of EM cascades are of restricted applicability. It was observed that near the shower axis, the NKG predicted densities were somewhat lower than those with the original NK formula [13] for s ⩾ 1.2. On the other hand, experimentally observed densities were noted to be larger than those given by the NKG far from the axis [8], which was inferred as a possible contribution from the muon decay. A shorter value of the second momentum of the distribution than in the NKG was observed [14] and a couple of years later a steeper profile was exhibited from an MC calculation near 100 GeV [15]. An improvement of the NKG function was subsequently proposed by adopting a modulated, longitudinal age parameter s_|| dependent effective Moliére radius so that

$$\begin{equation} \rho _{el}(r) = (mr_{m})^{-2}\rho _{{\rm NKG}}(r/m) \end{equation} \tag{ 6 }$$

where m = 0.78 − −0.21s_||.⁴

On the other hand, observing that the experimental LDD of electrons in EAS was steeper than that given by the ρ_NKG and was in better agreement with the MC calculations of Hillas [15] at lower energies, Capdevielle et al [10] introduced the notion of local age. After testifying the behaviour of the LAP on experimental lateral distributions [10] and reaffirming it with the Akeno observations [16], this approach was validated by the rapporteurs of the International Cosmic Ray Conferences during the period 1981 to 1985 [17]. The whole procedure was also employed in the calculation of the radio effects of EAS [18].

If electron densities are describing any NKG-like distribution f(x), where $x={r\over r_{m}}$, for two neighbouring points, i and j, we have the (local) lateral age

$$\begin{equation} s_{ij}^{\rm local} = {{\ln \big(F_{ij} X_{ij}^{2} Y_{ij}^{4.5}\big)} \over {\ln (X_{ij} Y_{ij})}} \end{equation} \tag{ 7 }$$

where F_ij = f(r_i)/f(r_j), X_ij = r_i/r_j and Y_ij = (x_i+1)/(x_j+1). More generally, if r_i → r_j, this suggests the definition of the LAP s^local_⊥(r) at each point :

$$\begin{equation} s_{\bot }^{\rm local}(r) = {1 \over {2x+1}} \left( (x+1) {{\partial {\ln f}} \over {\partial {\ln x}}} + (2 + \beta _{0})x + 2 \right) \end{equation} \tag{ 8 }$$

If β₀ = 4.5, f_NKG(r) with s_⊥ ≡ s^local_⊥(r) can be used to fit f in the neighbourhood of r.

Typical behaviour was predicted with a characterized minimum value of s^local_⊥(r) near 50 m from the axis, followed by a general increase at a large distance [10]. The relation s^local_⊥(r) = s_ij^local for $r=\frac{r_{i} +r_{j}}{2}$ was found to be valid for the experimental distributions (taking F_ij = ρ(r_i)/ρ(r_j) as far as they were approximated by monotonic decreasing functions versus distance).

Such a prediction was substantiated by the Akeno [16], North Bengal University (NBU) EAS experiment [19] and other experiments [10]. The LAP depends mainly on the logarithmic derivative of the density versus the distance as it appears in the relation (9); however, this pure mathematical approach may not be attained in practice at any radial distance, due to the experimental uncertainties arising mainly from the use of a finite number of detectors for the density measurements, triggering conditions and errors in the determination of the shower core position. Therefore, s^local_⊥(r) is estimated via relation (7) using physical bands of distance [r_i, r_j]; for experiments with very dense grids of detectors, such distance bands may be reduced to 5–10 m, but they may have to enlarge up to about 20 m for arrays with a lower resolution, as well as in the case of individual showers with large fluctuations. For very large and giant EAS, the interval [r_i, r_j] maybe required to exceed 100 m or so. We preferred to conserve the characteristic parameters of EM cascades in relation (7) including the value of the Moliére radius to facilitate the comparison with the experimental data, which is most frequently expressed in NKG formalism.

The dependence of s_⊥(r) on r rules out a consistent integration via relation (4) casting some doubt on the accurate relation between density and size; it was shown that such a dependence of s_⊥(r) on r is mainly a basic feature of pure EM cascades [10, 20].

For generating EAS events, the air shower simulation program CORSIKA [21] is exploited here. Here our discussions are mainly restricted to cosmic rays in the knee region of the primary spectrum. In this work, the high energy (above 80 GeV/n) hadronic interaction model QGSJET 01 version 1c [22] was used in combination with the low energy (below 80 GeV/n) hadronic interaction model GHEISHA (version 2002d) [23] or FLUKA [24], depending on the primary energy in the framework of the CORSIKA MC program version 6.600/6.735 [21] to generate EAS events. Note that the low energy interaction model GHEISHA exhibits a few shortcomings [25, 26] but the LDD of EAS electrons does not depend much on the low energy hadronic models, except at large distances [25]. Hence for very high energy events involving large radial distances we employed FLUKA [24]. A relatively smaller sample was also generated using the high-energy interaction model SIBYLL (version 2.1) [27] to judge the influence of the hadronic interaction models on the results.

The CORSIKA program allows one to choose either of the two options, the EGS4 (electron gamma shower system version 4) [28] and the NKG for obtaining a lateral distribution of the charge particles. The former option facilitates a detailed MC simulation of the EM component of a shower that incorporates all the major interactions of electrons and photons (see [11]), whereas the NKG option relies on an analytical approach rather than a full MC simulation. In the NKG option, the electron density of an EM sub-shower is calculated straightaway using the NKG function with a reduced Moliére radius [9, 10]. One gets better accuracy and more detailed information about the EM component with the EGS4 option, at the expense of long computing time. We underline here that the NKG option (subroutine NKG inside CORSIKA)⁵ is dealing mainly with the relations ((6), (7)) and not directly with equation (4).

We have considered the US-standard atmospheric model [31] with a planar approximation. The maximum primary zenith angle was restricted to 50°. The EAS events were generated mainly for proton and iron nuclei as primaries. A few events were also generated for γ as primary. Irrespective of the nature of the primaries, the slope of primary power law spectra was taken as −2.7 below the knee (3 × 10¹⁵ eV) and as −3.0 above. The EAS events were simulated at different geographical positions corresponding to the experimental sites of AKENO [16], KASCADE [30] and NBU [32]. The magnetic fields are provided accordingly. On the observational level, the kinetic energy thresholds were chosen as 3 MeV for electrons (e⁺ and e⁻) irrespective of the primary species and energies.

3.1. Generation of the EAS Monte Carlo library

The simulated shower library consists of more than 30 000 EAS events with the EGS4 option and more than 180 000 events with the NKG option in the primary energy interval of 10¹⁴ eV to 3 × 10¹⁶ eV. In order to appreciate the asymptotic tendencies at ultra high energies, our library has also been enriched by about 1000 events simulated at E_o = 5 × 10¹⁷ and 10¹⁸ eV for proton and iron primaries: apart from the thinning factor which is taken as 10⁻⁶ with optimum weight limitation [33] (i.e., all particles are followed up to an energy E_th, where E_th/E_o = 10⁻⁶, after which only one of those particles is tracked giving appropriate weight to it) the simulation conditions are here identical concerning the hadronic interaction models, the zenith angle range, both the EGS4 and the NKG options. Here we would like to specify that the optimized thinning factor 10⁻⁶ with the optimum weight limitation for the CORSIKA version used here is considered as the best compromise between the computing time and the accuracy at ultra high energies [34].

In all cases involving the EGS4 option, the longitudinal development is restored numerically and that s_|| is computed [35], instead of relation (2), by:

$$\begin{equation} s = \exp \left[ { {{2}\over {3}} \times \left\lbrace {1+{{\alpha }\over {t}} -\tau } \right\rbrace } \right] \hspace{7.11317pt} \ \ \ {\rm with} \ \tau = {{t_{\rm max}\over {t}}}, \ \ \alpha = \ln {{N_{\rm max}\over {N_{e}}}} \end{equation} \tag{ 9 }$$

where t_max, N_max are respectively depth and size read at the cascade maximum.

3.2. The NKG and the EGS options

We used both the EGS4 and the NKG options simultaneously for about 30 000 events. In figure 1(a) we compare the LDD of EAS electrons obtained with the stated two options for proton, iron and gamma primaries. It is noted that the NKG option gives a higher density with a steeper radial distribution compared to the EGS option. A small density excess appears for the pure EM cascades near the axis for the NKG option; such an excess presents also in the proton initiated air showers. However, for the proton showers, a tolerable agreement between the output of the two options was noted over a large band of densities between radial distances of 10–100 m from the shower axis; it reconfirms that for proton and photon initiated showers the NKG option is quite useful to calculate a large number of cascades in a short time.

Zoom In Zoom Out Reset image size

For Fe primaries, both the options indicate an older profile near the axis. The NKG option was found to give an excess density between 2–10 m distance. The average energy of the positrons was quite a bit lower in the case of iron initiated showers and the cross section of positron annihilation becomes more important for the lower part of the cascade. This effect is probably enhanced by the longer path of the electrons in the geomagnetic field and the larger energy loss by ionization. The NKG option is, therefore, not so accurate for the simulation of heavy nuclei initiated showers with large zenith angles after the shower maximum. For vertical showers, the output of the NKG option is nevertheless acceptable. At larger distances a slight deficit in the densities appears with the NKG option; this probably comes from the different treatment of the multiple Coulomb scattering in the NKG option than in the EGS. Furthermore, Bhaba and Moller scattering are treated in complement with the separated MC procedures. On the other hand, the geomagnetic field of the earth enhances the path of the muons. Consequently their loss by ionization and their decay give more electrons, which is not incorporated in the NKG option. Besides, the NKG option does not accommodate photo-production inside the EM sub-cascades.

The simulated data were analysed using the reconstruction algorithms developed to obtain shower size and shower age. We adopted two different methods. First following the traditional approach, we estimated the basic shower parameters by fitting the density data to the NKG function by the method of chi-square minimization through an iterative procedure based on the method of steepest decent. As for example, in figure 1(b) the simulated particle densities at different radial distances are plotted along with the fitted curve obtained with the NKG function. Here it is noted that the majority of the EAS groups traditionally estimate basic shower parameters based on the NKG function. The error in estimating lateral shower age in the shower size interval 10³–10⁵ particles (corresponding to the primary energy range 10¹⁴–3 × 10¹⁵ eV) was found to be ±0.03 for the QGSJET model and ±0.05 for the SIBYLL. The larger error for the SIBYLL model seems solely statistical, due to the generation of a relatively fewer number of EAS events using the later model.

The local age for EAS charged particles was computed for each individual event straightaway, applying equation (7). When estimating the LAP, the main sources of error are the fluctuations in particle density and the uncertainties in radial distance estimation. In simulated data the radial distance of each particle is known with a high accuracy. In this work the error in the LAP, due to uncertainties in radial distance estimation, was kept small by taking small radial bins. For minimizing the statistical fluctuations in particle density at different radial bins, a large number of events need to be considered. In this analysis, the error of the LAP for EAS with the primary energy in the PeV range remains within 0.05 for 10 m < r < 250 m, whereas for r < 10 m or when r > 300 m the error of the LAP is found to be higher, about 0.1. At higher primary energies (5 × 10¹⁷–10¹⁸ eV) the error of the LAP is found at about 0.12 near the core, which decreases to about the 0.07 level when 20 m < r < 300 m but increases again with the radial distance and reaches to about 0.15 when r approaches 1000 m.

The variation of the LAP with the radial distance from the shower core is shown in figure 2. It is known from previous studies [10] that with an increase of the radial distance, the LAP initially decreases, reaching a minimum at around 50 m and then increases, as was also noted from the experimental data [16, 19,10]. Here we noticed two other interesting features (figure 2(a)–(c)): the local age again starts to decrease at around 300–400 m. To examine whether the experimental data also demonstrates a fall in the local age at large radial distances, we compute the local age from the LDD data of total charged particles, as measured by the AGASA experiment [36] for primary energy 2 × 10¹⁸ eV and compared these values with our simulation results in figure 2(c). The experimental data clearly support the trend predicted by the simulation results at larger radial distances. The characteristic high–low–high kind of radial variation in the local age at relatively smaller distances (within 300 m or so) could not be substantiated by the AGASA data, due to the large separation of the detectors of the array. It is worthwhile mentioning that there was an indication for such a decrease of s^local_⊥(r) at around 300 m in the experimental results obtained by Akeno [19] and KASCADE-Grande [37]. Such behaviour is also depicted in figure 3 of the Kascade report [30], where one may notice a maximal deficit at 50–80 m in the ratio of the measured and the fitted electron densities, as well as an excess at large distances when fitting with the NKG formula. The later feature, however, has not been thoroughly investigated.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

These findings are important for an analysis of very large air showers observed/to be observed by the KASCADE-Grande, AGASA, AUGER, Yakutsk, Telescope Array involving large radial distances. The large EAS experiments often treat charged particle densities at large radial distances, such as at 500, 600, 1000 m from the shower core as an estimator of the primary particle energy, though such a technique involves several uncertainties [38]. These findings, of a rapid change in the slope of the radial distribution of electrons at large radial distances, suggest that more controls should be adopted in the estimation of the primary energy of large showers, for instance by taking particle densities at more than one radial distance.

Another important observation is that in general the nature of the variation of the local age with the radial distance appears nearly the same for all of the primary energies, i.e. the nature of the variation is practically independent of the energy of the shower initiating particles, which implies that the local age (or the lateral distribution of electrons in EAS) exhibits some sort of scaling behaviour in respect to the radial dependence from the shower core.

To examine systematically the influence of the effective Moliére radius on the shape of the lateral distribution of charged particles in EAS, we study the radial variation of the LAP for different effective Moliére radii, which is shown in figure 3(a) for a proton primary with a primary energy 5 × 10¹⁷ eV.

A string-like feature emerged with two nodes, as seen from the figures, one close to the shower core while the other at around r ∼ 400 m that increases slowly with the effective Moliére radius; the effective Moliére radius behaves somewhat like the tension in a piece of string. Beyond the second node, however, the LAP is found to decrease monotonically with an increase of the radial distance, irrespective of the choice of the effective Moliére radius.

In order to explore the inherent cause of such a feature of the LDD of electrons in hadron initiated EAS, we studied the radial variation of the local age for γ ray initiated showers, and one such plot at the primary energy 10¹⁵ eV is shown in figure 3(b). We found the similar nature of radial variation of the LAP as in the hadron initiated showers. As ascertained in previous simulations [10], the behaviour of s^local_⊥(r) comes mainly from the discrepancies between the EGS output and the NKG function, i.e. between the rigorous descriptions of the EM cascade adopting the basic EM processes as well as the Moller, Bhabha scattering, positron annihilation, dependence of the cross section on energy on the one hand, Approximation B combined with Landau and small angle approximations in a single description of the multiple Coulomb scattering on the other hand. When the experimental data [30, 35] are superimposed on figure 3, we understand that a reduced Moliere radius (between 20–50 m) is favoured for all primary energies, implying a dramatic reduction in the mean scattering angle connected with the scattering energy of 21 MeV.

To explore the physical nature associated with the lateral shower parameter, if any, we studied the details of the characteristics of the shower age. For the local age, we considered two different parameters: the minimum value corresponds to the local age at the radial distance, about 40 m, and an average value between 40 to 300 m.

5.1. Distribution of shower age and its fluctuation

The distributions of the LAP and the lateral shower age were studied for the primary energy range 3 × 10¹⁴ to 3 × 10¹⁶ eV. The fluctuations in the shower age were found to be larger for proton initiated showers compared to those initiated by a heavier primary. If we consider a small primary energy bin instead of a wide one, for instance by selecting the showers inside a small muon size band, we observed that both the distributions of p and Fe could be separated, which in fact becomes very contrasted as shown in figure 4(a); this approach, if adopted with the experimental data, may yield important information on the primary composition around the knee region.

Zoom In Zoom Out Reset image size

The fluctuations (σ) (variance) in the LAP in different shower age bins are estimated and as a function of the shower size are drawn in figure 4(b) for proton and iron primaries for the interaction model QGSJET. In accordance with expectations [39], the fluctuations in the LAP were found to be larger for the proton initiated showers in comparison to those initiated by the primary Fe, except at lower energies.

5.2. Longitudinal shower age versus lateral shower age

For each simulated event, the longitudinal shower age was also estimated using the relation (9) and the difference between the two age parameters, longitudinal age and average local age, was obtained. The frequency distributions of the differences between the two age parameters for proton and Fe primaries are given in figure 5. The frequency distribution for where the proton primary peaks, at around 0.2, is consistent with the early observations [40, 10], whereas for Fe initiated showers the peak difference is much lower, at about 0.07. However, for a non-negligible fraction of events, the differences between the two shower age parameters were found to be substantial.

Zoom In Zoom Out Reset image size

5.3. Variation of the shower age with electron size

In figure 6, we plot the local shower age at a radial distance of about 40 m (minimum value) versus an average shower size, obtained from the simulation results for both proton and iron primaries at the Akeno and KASCADE locations. The corresponding observational results extracted from the Akeno and KASCADE experimental data are also given in figure 6. For the Akeno experimental data, we extracted a minimum local age for the different shower sizes from reference [16], whereas for KASCADE we estimated it from their measured lateral distribution [30, 37,42].

Zoom In Zoom Out Reset image size

The comparison of the simulated results with the experimental observations from both the Akeno and KASCADE EAS experiments (figure 6) indicate a need for a change in the primary composition towards a heavier primary, as the energy increases across the knee of the primary energy spectrum. The KASCADE group also reached a similar conclusion using the shape parameter instead of the shower age [42], as well as from the study of the muon content in EAS [43]. The present data of the LHC, especially the pseudo-rapidity density distributions, suggest larger multiplicities and inelasticities than in the models used in the CORSIKA simulations [44]. However, up to an energy of 2.6 × 10⁷ GeV, this could result in a very small reduction of the reported enhancement of the primary mass with energy.

From the present analysis we conclude the following.

• (1)  
  The lateral distribution of electrons in EAS exhibits some sort of scaling (energy independent) behaviour in terms of the local age. The characteristic feature of the local age versus the radial distance curve is that with an increase of the radial distance, the local age decreases initially and reaches a minimum, at around 40 m, then it starts increasing, attaining a local maxima, at around 300–400 m, and then starts decreasing again. Such a feature appears to be independent of the energy of the EAS initiating particle, at least there is no strong dependence on the primary energy. Such a characteristic radial variation in the local age is found as a generic feature of EM cascades.
• (2)  
  The local age offers a good solution towards an unambiguous estimation of the shower age. Since the shower age varies with the radial distance, even for the modified NKG functions, a comparison of the lateral shower age of different EAS experiments is not meaningful, as the radii of the shower discs naturally differ from experiment to experiment, depending on the experimental set up. Even in a single EAS experiment different events have different radial extensions and thus a lateral age obtained through fitting with the NKG function seems ambiguous. The local age at a particular distance (say at about 40 m where it takes the minimum value) is, however, not always practical owing to the large fluctuation in the electron density data in a real measurement. A rational idea could be to take some sort of average local age between the first minimum, at around a core distance of 40 m, and the subsequent local maximum, at around 300 m. Experimentally the radial variation of the LAP can be checked properly with a full coverage detector array, like the ARGO-YBJ [45].
• (3)  
  The use of a reduced effective Moliére radius in the NKG function leads to a roughly constant age over a radial distance up to about 200 m [41], but the shower age estimated in such a manner is found take quite a higher value than that of the longitudinal age. More importantly, even with a reduced effective Moliére radius, the local age is still found to vary after a radial distance of about 200 m.
• (4)  
  The lateral age offers a good estimator of the longitudinal development of an EAS cascade, as already noted in some earlier works [46, 40,10]. However, the parameter correlates with the stage of the shower development on a statistical basis; the average of this parameter increases as air showers traverse an increased thickness of atmosphere. The experimental observations [20] also substantiate such behaviour. The distribution of the differences between the local age and the longitudinal age also indicate the strong correlation between the two ages. Such a feature has been noted for two different hadronic interaction models, the QGSJET and the SIBYLL, and hence appears robust. It is imperative to examine such correlations using EPOS [47], the only model that seems to be providing quite a consistent description of the longitudinal and lateral EAS profiles [48], which would need to be performed in future work.
• (5)  
  The fluctuation of the LAP was found to be sensitive to the nature of the primary particle. However, the level of uncertainty in determining the lateral shower age from the experimental data is comparable with the magnitude of fluctuation and hence deriving any firm conclusion on the nature of the primary only from the shower age fluctuation is difficult. If showers within a small bin of the primary energy could be selected, for instance by considering shower events in a small muon size interval, the distribution of shower age was found to be quite sensitive on primary composition; this property may be useful (in conjunction with other primary sensitive parameters) to extract information about primaries.

The comparison of the simulation results with the Akeno and the KASCADE observations in respect to a variation of the shower age with shower size around the knee indicates a change in the primary composition towards heavier primaries across the knee. This finding supports the results obtained from the study of the muon content in EAS. The 3D plot of the shower size versus muon size and shower age may improve accuracy over the conventional approach: deducing the composition through the shower size versus muon size curve. It would be an interesting task to apply such a method for a composition study using observed EAS data from a closely packed air shower array with the facility of concurrent muon measurements, such as the GRAPES experiment at Ooty [49].

The authors sincerely thank two anonymous referees and a board member for critical comments and helpful suggestions. AB would like to thank the Department of Science and Technology (Govt. of India) for support under the grant no. SR/S2/HEP-14/2007.