text:
# The information from muon arrival time distributions of high-energy EAS as measured with the KASCADE detector

## 1 Introduction

The time delay of the particles inside the front of Extensive Air Showers (EAS) and the temporal structure of different EAS components are a subject of longstanding interest of cosmic ray research. In fact the first experimental studies have been performed in 1953 by Bassi, Clark and Rossi  and Jelly and Whitehouse , followed by many others -. The renewed interest arises from recent measurements using advanced detector facilities like the COVER PLASTEX detector within the GREX array -, the EAS TOP array  or the facilities of the KASCADE experiment . Investigations of the structure of the EAS disc by the KASCADE experiment are focused in particular on the muon component -. Due to the reduced influence of multiple Coulomb scattering of GeV muons and the absence of absorption, muons travel nearly in straight lines from the locus of production to the observing detector. Thus muon arrival time distributions, observed at large distances from the shower axis, are expected to map the longitudinal EAS development and to reflect the distribution of production heights via the time-of-flight of the muons. With some simplifications, using a triangulation procedure, the distribution of the production heights could be estimated from the time delay of the observed muons relative to the arrival time of the shower centre . Basically the same information could be alternatively deduced from the distribution of angles of muon incidence relative to the shower axis . Corresponding measurements are also a current subject of KASCADE investigations using the Muon Tracking Device . Angle-of-incidence observations or the combination with arrival time measurements (as discussed with the term Time Track Complementary ) do not reveal additional basic information  though they provide an interesting practical alternative.
The basic sensitivity of muon arrival time distributions (whose phenomenological features seen with the KASCADE experiment are reported in Ref. ) to longitudinal EAS development and to the elongation rate  puts the question which particular parameters of the observed distributions provide useful signatures of the mass of the EAS primary and under which conditions they are helpful for the determination of the mass composition of the primary cosmic rays. The present paper addresses experimentally the question of the sensitivity to the primary mass on the basis of results from the KASCADE experiment, in particular of EAS event-by-event measurements of the temporal EAS structure (shower thickness) observed at relatively small distances to the shower axis $`R_\mu <100`$ m for the primary energy range around the knee. Methodically non-parametric statistical analysis techniques of multivariate distributions  are applied for the mass classification of the observed EAS by their time parameters, derived from the muon arrival time distributions, and by their correlations with other EAS observables (like the shower size $`N_e`$ and the muon content $`N_\mu `$). For muon arrival time distributions the methods are outlined in Refs. . These procedures avoid the bias of pre-chosen functional forms and compare to reference patterns, derived from extensive Monte Carlo simulations, by use of the EAS Monte Carlo simulation code CORSIKA , thus including all natural EAS fluctuations. By applying additional detailed detector simulations based on GEANT , the distortions arising from the experimental conditions  are taken into account.
The analysis of the EAS observations in terms of the mass of the primary introduces necessarily some model dependence by the high-energy hadronic interaction models invoked as generators of the Monte Carlo simulations. Therefore the results are subject of the uncertainities of the hadronic interaction model, used for the analysis. A possible way to evaluate the quality of a model is to derive the primary mass composition from the analysis of different observables and to consider the agreement or disagreement of the results. This idea has been worked out in a consistent and efficient manner  with an extended set of observables of the KASCADE experiment, showing that there exist systematic differences in the estimate of $`lnA`$, e.g., derived from different combinations of (correlated) observables , revealing the limitations of popular models like QGSJET , VENUS  and SIBYLL . This indicates that improvements of the models are urgently needed on basis of experimental findings. The QGSJET model is actually considered to be that one of the best internal consistency . Hence we refer furtheron to the QGSJET model.
The muon arrival time measurements imply a particular selection of EAS, with a distortion of the original mass composition in the subset of classified EAS. This is due to the energy threshold $`E_{th}`$ (2.4 GeV for KASCADE) of the muon detection and the condition of a minimum multiplicity $`n_{th}`$ ($`=3`$ muons detected by the KASCADE Central Detector) for reconstructing an arrival time distribution. The distortion is dependent on the particular values of $`E_{th}`$, $`n_{th}`$, on the mass $`A`$ of the primary and the muon content $`N_\mu `$, but also on the radial distance $`R_\mu `$ of the arrival time measuring detector from the shower axis. In order to restore the original mass composition corresponding efficiency (acceptance) correction factors have to be applied to the identified mass groups. These factors can be only determined by Monte Carlo simulations, studying the cuts implied by the experiment and their effect on the efficiency. An efficient and sensitive test of the interaction model used for the simulations and of the particle tracking procedures is enabled by applying the calculated correction factors to the results found at different $`R_\mu `$ and to look for the consistency of the resulting mass composition. Such a test of the Monte Carlo simulations and of the QGSJET model, especially in view of predictions of the muon energy spectra and muon densities and their fluctuations in event-by-event observations, is an important aspect of muon arrival time studies of this paper. The results add to the conclusions about the difficulties of the model in interpreting consistently muon density measurements with different muon energy thresholds .

## 2 Muon arrival time measurements

Measured arrival times $`\tau _{1\mu }<\tau _{2\mu }<\tau _{3\mu }<\mathrm{}`$ of muons registered by timing detectors at a distance $`R_\mu `$ from the shower axis, corrected by $`\pm R_\mu \mathrm{tan}\theta /c`$ ($`c`$ \- speed of light) in order to eliminate the distortions due to the shower inclination, refer to an experimentally provided zero-time. Depending on the choice of the kind of zero-time reference, there are two different types of muon arrival time distributions. By the use of the arrival time $`\tau _c`$ of the shower core as reference global arrival times are observed :

$`\mathrm{\Delta }\tau _i^{glob}(R_\mu )=\tau _{i\mu }(R_\mu )\tau _c,`$ (1)
e.g. for the foremost muon registered at $`R_\mu `$:

$`\mathrm{\Delta }\tau _1^{glob}(R_\mu )=\tau _{1\mu }(R_\mu )\tau _c.`$ (2)
This type of time distributions informs about the curvature of the shower front as well as about the structure of the disc. The arrival time of the shower core is difficult to determine with sufficient precision. Therefore the analysis has been restricted to ”local” times, which refer to the arrival of the foremost muon locally registered by the detector:

$`\mathrm{\Delta }\tau _i^{loc}(R_\mu )=\tau _{i\mu }(R_\mu )\tau _{1\mu }(R_\mu ).`$ (3)
(with omitting further the label loc), informing only about the thickness and the structure of the muon disc. Implications of observations of local muon arrival time distributions due to the fluctuations of the arrival of the first registered muon have been discussed in Ref. . For event-by-event observations with a fluctuating number of muons (multiplicity), the individual relative arrival time distributions can be characterised by the mean values $`\mathrm{\Delta }\tau _{mean}`$, and by various quantiles $`\mathrm{\Delta }\tau _q`$, like the median $`\mathrm{\Delta }\tau _{0.50}`$, the first quartile $`\mathrm{\Delta }\tau _{0.25}`$ and the third quartile $`\mathrm{\Delta }\tau _{0.75}`$ (for details of the definition see Refs. ). Their mean values and dispersion (standard deviations) represent the time profile of the EAS muon component.
KASCADE, whose layout is described in more detail in Ref. , is a multidetector system, installed in Forschungszentrum Karlsruhe (110 m a.s.l.), Germany, for the observation of extensive air showers in the primary energy range around the knee. One part is an array of 252 detector stations, distributed over an area of $`200\times 200`$ m2 on a grid of 13 m spacing for measuring the electron-photon component and the muon component with a threshold of 5 MeV and 230 MeV, respectively and providing the basic information about arrival direction, core position, electron and muon sizes of the observed EAS. In particular, from the data of the field array the so-called truncated muon number $`N_\mu ^{tr}`$, i.e. the muon density integrated between 40 m and 200 m, is derived and used in the KASCADE case (between $`10^{14}`$ and $`10^{16}`$ eV) as an approximate mass independent energy identifier . The location of the shower core can be determined (inside the fiducial area) with an uncertainty of less than 3 m. The arrival direction of the shower is reconstructed with an uncertainty better than $`0.5^{}`$. Details of the reconstruction procedures are described in Ref. . As additional muon detector Limited Streamer Tubes tracking detectors have been installed in an underground tunnel for measurements of the lateral distribution ($`E_{th}=0.8`$ GeV) and of muon angles-of-incidence distributions .
The muon arrival time measurements use, in particular, the facilities of the Central Detector  of KASCADE. It is basically an iron sampling calorimeter (with an area of $`16\times 20`$ m2, set up with liquid TMS and TMP ionisation chambers) for the identification and energy measurement of hadrons. In the basement of the set-up, below 3800 t of iron and concrete, large-area position-sensitive multiwire proportional chambers (MWPC)  are operated for the identification of muons with 2.4 GeV energy threshold. The performance of the detection system of the MWPC is improved by a layer of streamer tubes . The trigger plane of the calorimeter is a system of 456 plastic scintillation detector elements ($`47.5\times 47.5\times 3`$ cm3 in size, each separated by a wavelength shifter) for providing a fast trigger signal (in addition to the trigger from the field array) for the MWPC, and for the timing measurements. Fast electronics  records low (muons) and high energy deposits (cascading hadrons). In order to remove the signals from cascading hadrons an upper limit of the energy deposit of 30 MeV in each of the 456 scintillation counters is imposed. An amount of 24 millions of KASCADE events are used for the analysis . Muon arrival time distributions have been reconstructed for muons with $`E_{th}=2.4`$ GeV and for events with a number $`n3`$ of registered muons; after applying also some general cuts concerning the core position (within 100 m from the array centre), the angle of EAS incidence ($`<40^{}`$) and $`log_{10}N_\mu ^{tr}`$ ($`>3.6`$) the sample shrinked to approx. 240 000 showers. The phenomenological features of the observed muon arrival time distributions have already been reported in Ref. , where further experimental details are communicated.

## 3 EAS simulations

The interpretation of the measured muon arrival time distributions and their correlations with other shower parameters need a-priori knowledge to be deduced from Monte Carlo simulations of the EAS development. The present analysis is based on simulations with the code CORSIKA (version 5.62) with a full and detailed simulation of the detector response. The simulations use the QGSJET (version 1998) model  as generator for high-energy interactions and GHEISHA  for interactions below $`E_{lab}=80`$ GeV. The electromagnetic part is treated by the EGS4 program . Earth magnetic field, observation level and particle detection thresholds have been chosen in accordance to the experimental situation. The U.S. standard atmosphere  has been adopted for the simulations. The simulations have been performed for three different classes of primaries: protons (H) for the light group, oxygen (O) for the CNO-group and iron (Fe) for the heavy group. The energy range covered by the simulations extends from $`5.010^{14}`$ eV to $`3.0610^{16}`$ eV and the zenith angles comprise the range of $`0^{}\theta 40^{}`$. The centres of the simulated showers have been positioned inside a quadratic area of $`210\times 210`$ m2, slightly extending the KASCADE area. For each primary type approximately 90000 showers have been simulated, with decreasing numbers of shower events for higher primary energies and the larger zenith angles, due to restrictions in the computing time. Finally a sample for $`0^{}\theta 24^{}`$ (corresponding to a third of the considered $`\mathrm{sec}\theta `$ \- range of $`(1,\mathrm{sec}40^{})`$) was used in the analysis. The simulated primaries (separately for each of the cases H, O, Fe) have been weighted according to:

$`dN=const.E_0^\gamma \mathrm{sin}\theta \mathrm{cos}\theta D_{core}dE_0d\theta dD_{core}`$ (4)
where $`dN`$ is the number of primaries with energy between $`E_0,E_0+dE_0`$, zenith angle between $`\theta ,\theta +d\theta `$, intersecting the plane of the detector array between $`D_{core},D_{core}+dD_{core}`$. A fixed spectral index $`\gamma =2.7`$ over the whole primary energy

range is adopted for all primaries (H, O, Fe). More details about the simulations are given in Refs.  and . As input for the further analysis of the muon arrival time distributions and their comparison with model predictions the following shower parameters and their correlations are regarded:

* the shower size $`N_e`$;
* the truncated muon number $`N_\mu ^{tr}`$;
* the multiplicity $`n`$ of muons ($`E_{th}=2.4`$ GeV) detected in the facilities of the Central Detector;
* the quantiles $`\mathrm{\Delta }\tau _q`$ of the local muon ($`E_{th}=2.4`$ GeV) arrival time distributions at various distances $`R_\mu `$ from the shower centre;
* the reduced quantiles $`\mathrm{\Delta }\tau _q^{}=\mathrm{\Delta }\tau _q/\rho _\mu `$ of the local muon arrival time distributions ($`E_{th}=2.4`$ GeV), i.e. the quantiles divided by the density $`\rho _\mu `$, where $`\rho _\mu `$ is estimated from the observed multiplicity and the effective area of the muon detector set-up.

Fig. 1 displays frequency distributions of some shower parameters, prepared by the simulation for different primaries and compared to the actual experimental observations for particular $`log_{10}N_\mu ^{tr}`$ and $`R_\mu `$ bins. Each distribution is normalised (only for the presentation of $`N_e`$ the distribution of oxygen-induced EAS is additionally shown). There are some discrepancies in the multiplicity distribution between simulations and experimental observations. They may arise from an imperfect adoption of the primary energy spectrum and from efficiency effects not fully corrected for. Since they also affect the fluctuations of the local time parameters originating from the multiplicity dependence (discussed in detail in Ref. ), this analysis uses preferentially the reduced quantities $`\mathrm{\Delta }\tau _q^{}`$. The fluctuations largely cancel in the reduced parameters. Actually it has been argued that the reduced quantiles $`\mathrm{\Delta }\tau _q^{}`$ could exhibit enhanced mass sensitivity since they include approximately the $`\mathrm{\Delta }\tau _q(R_\mu )\rho _\mu (R_\mu )`$ correlation . However this sensitivity is obscured by the limitations of the experimental response. CORSIKA simulations show  that the age parameter is carrying relevant information about the longitudinal development. This aspect has not been worked out by the present analysis.

## 4 Non-parametric multivariate analyses

Non-parametric statistical methods are applied in studies of multidimensional distributions of observables allocating the single observed events to different classes (in our case to proton, oxygen, or iron primaries) by comparing the observed events with the simulated distributions without using a pre-chosen parameterisation. The methods of decision making and procedures of applications to cosmic ray data analyses are extensively described in Refs.  and outlined for analyses of muon arrival time distributions in Refs. . The procedures take into account the effects of the EAS fluctuations in a quite natural way and are able to specify the uncertainties, by an estimate of the true-classification and misclassification probabilities. The classification probabilities are determined by the extent to which the likelihood functions of the individual classes, derived from the simulations, are overlapping. Basically such pattern recognition methods, using trained neural networks or Bayes decision rules, depend on the hadronic model generating the reference patterns for the experimental observables. The concept of the present analysis is to determine in a first step the sensitivity of various observable correlations to the mass of the primary, by applying the one-leave-out test to the sample of distributions prepared by Monte Carlo simulations (QGSJET model). The one-leave-out test determines the probability that a multidimensional event, taken from the considered (simulated) distribution, will be correctly (”true”) or incorrectly (”false”) classified by the pattern recognition procedure (see appendix of Ref. ). With respect to the experimental data sample the studies are performed not only for six different $`log_{10}N_\mu ^{tr}`$ ranges between 3.60 and 5.00 (covering the knee-range) but also for four different $`R_\mu `$-ranges (45-65 m, 65-80 m, 80-90 m, 90-100 m). This allows to explore trends of the classification and misclassification probabilities and other quantities with varying $`R_\mu `$. Using the determined classification probabilities the mass composition is reconstructed from the data samples measured at different $`R_\mu `$-ranges. The resulting mass composition generally differs from the primary mass composition searched for, since the analysed data samples represent selections in the multidimensional space of all EAS observables. Hence efficiency correction factors have to be determined, whose application should lead to a primary mass composition consistent for all $`R_\mu `$-ranges. The determination of the efficiency correction factors invokes again Monte Carlo simulations, and the extent of agreement of the primary mass composition resulting from different $`R_\mu `$ ranges provides a test for the procedures and the generator of the Monte Carlo simulations.

### 4.1 The classification and misclassification probabilities

As example in Fig. 3 the true-classification and misclassification probabilities, inferred with the one-leave-out test from the simulation sample, are displayed for H, O, Fe primaries and for different combinations of observables of EAS events. The events have been registered with the conditions of muon arrival time measurements in a particular $`\{log_{10}N_\mu ^{tr},R_\mu \}`$ \- range. With little surprise it is noted that the mass discrimination is dominated by the $`\{N_\mu ^{tr},N_e\}`$ correlation (see Fig. 3), and it is also evident that local time parameters (shown for the median $`\mathrm{\Delta }\tau _{0.50}`$ or the third quartile $`\mathrm{\Delta }\tau _{0.75}`$ and their reduced values $`\mathrm{\Delta }\tau _{0.50}^{}`$, $`\mathrm{\Delta }\tau _{0.75}^{}`$) have minor influence on the discrimination. We emphasise that this statement holds for local arrival time distributions in the studied range of relatively small distances $`R_\mu `$ ($`<`$ 100 m) from the shower centre. Theoretical studies

 indicate more pronounced differences in the temporal structure of the shower disc with increasing $`R_\mu `$ and $`log_{10}N_\mu ^{tr}`$. In contrast to results ignoring the detector response  the use of reduced parameters $`\mathrm{\Delta }\tau _{0.50}^{}`$, $`\mathrm{\Delta }\tau _{0.75}^{}`$ does not significantly improve the mass discrimination. Fig. 3 shows the true classification probabilities (and uncertainties) averaged over all classes with increasing $`log_{10}N_\mu ^{tr}`$ indicating a marginal systematic improvement if the time information is added. The uncertainties are determined by the bootstrap method , which basically consists of applying the classification procedure for a test sample several times, thus deriving an average value and the variance.

### 4.2 Reconstruction of the mass composition from the observed data samples

The true-classification $`P_{ii}`$ and misclassification

$`P_{ij}`$ probabilities, deduced for all studied $`\{log_{10}N_\mu ^{tr},R_\mu \}`$ ranges, are used for the reconstruction  of the mass composition of the samples of registered events by inverting the system of linear equations:

$`\begin{array}{c}n_H^{}=n_HP_{HH}+n_OP_{OH}+n_{Fe}P_{FeH}\hfill \\ n_O^{}=n_HP_{HO}+n_OP_{OO}+n_{Fe}P_{FeO}\hfill \\ n_{Fe}^{}=n_HP_{HFe}+n_OP_{OFe}+n_{Fe}P_{FeFe}\hfill \end{array}`$ (8)
where $`n_H`$, $`n_O`$, $`n_{Fe}`$ are the true numbers, defining the mass composition in the sample $`N=n_H+n_O+n_{Fe}`$, getting altered to $`n_H^{}`$, $`n_O^{}`$, $`n_{Fe}^{}`$ because of the misclassifications. In Fig. 4 the reconstruction of an arbitrarily adopted mass composition of the event samples (displayed by the corresponding symbols), is shown with the resulting uncertainties as varying with $`log_{10}N_\mu ^{tr}`$. The result of the application of the reconstruction procedures to experimental samples measured with KASCADE, observed at different distances $`R_\mu `$ from the shower centre, are shown in the upper panels of Fig. 6. The uncertainties are influenced by the limited statistical accuracy of the data samples, in particular for the bins with larger $`log_{10}N_\mu ^{tr}`$ values i.e. primary energies. Nevertheless it is obvious that mass compositions of measured KASCADE samples differ for different $`R_\mu `$ bins, since the observation conditions lead to mass dependent differences in the observation efficiency at different $`R_\mu `$, leading to distortions of the deduced primary mass composition. This feature should be removed after applying a correct efficiency correction procedure. That necessarily implies again the use of the particular hadronic interaction model.

### 4.3 Reconstruction of the primary mass composition

Efficiency correction factors ($`C_H`$, $`C_O`$, $`C_{Fe}`$) have been calculated in order to adjust the mass composition of measured KASCADE samples ($`P_H`$, $`P_O`$, $`P_{Fe}`$) to the primary mass composition ($`P_H^{}`$, $`P_O^{}`$, $`P_{Fe}^{}`$) according to the relation :

$`P_H^{}:P_O^{}:P_{Fe}^{}={\displaystyle \frac{P_H}{C_H}}:{\displaystyle \frac{P_O}{C_O}}:{\displaystyle \frac{P_{Fe}}{C_{Fe}}}.`$ (9)
As a first step, the simulated spectra given by eq. 4 have been normalised to the same value (=1000) for each type of primaries $`A`$ (H, O, Fe), and for all simulated cases of $`E_0`$, $`\theta `$ and core position ranges. For a given primary $`A`$, the detected spectra will appear distorted at the ground level due to the absorption in the atmosphere and the selection cuts. Only a fraction of the original events are detected by the KASCADE detector and reconstructed successfully. The major influence on the values of the efficiency correction factors arises from the applied cuts on the reconstructed shower events. The correction factors depend on the $`R_\mu `$, $`\theta `$ and $`log_{10}N_\mu ^{tr}`$ ranges, on the multiplicity threshold for the 2.4 GeV muons and, of course, on the primary type. For a given primary $`A`$, in a certain $`\{R_\mu ,R_\mu +\mathrm{\Delta }R_\mu ;\theta ,\theta +\mathrm{\Delta }\theta ;log_{10}N_\mu ^{tr},log_{10}N_\mu ^{tr}+\mathrm{\Delta }log_{10}N_\mu ^{tr};n,n+\mathrm{\Delta }n\}`$ (multidimensional) bin, the correction factor ($`C_A`$) is given by the sum of the weights ($`w_i`$) of the $`p`$ simulated events which are accepted showers in the multidimensional bin:

$`C_A={\displaystyle \underset{i=1}{\overset{p}{}}}w_i.`$ (10)
The weights $`w_i`$ are determined via the normalised mass spectra (eq. 4) which get filtered by the observation conditions. Considering the equivalent number of events

$`q`$  given by:

$`q={\displaystyle \frac{\left(\underset{i=1}{\overset{p}{}}w_i\right)^2}{\underset{i=1}{\overset{p}{}}w_i^2}},`$ (11)
the statistical uncertainty (error of the mean value) has been calculated for $`C_A`$ by:

$`\delta C_A=C_A{\displaystyle \frac{\sqrt{q}}{q}}={\displaystyle \frac{C_A}{\sqrt{q}}}.`$ (12)
Fig. 6 displays for different $`R_\mu `$-ranges the

$`log_{10}N_\mu ^{tr}`$-variation of the mass composition of the measured KASCADE samples extracted from the $`\{N_\mu ^{tr},N_e,`$ $`\mathrm{\Delta }\tau _{0.50}^{},\mathrm{\Delta }\tau _{0.75}^{}\}`$ correlation and the primary mass composition resulting after the correction concerning the biased acceptance by the specific observational conditions. In the lower part of Fig. 6 the results are shown as variation of $`lnA`$ with $`R_\mu `$ for different $`log_{10}N_\mu ^{tr}`$ ranges. The uncertainties of the correction factors (eq. 12) have been included by propagation via eq. 9. It should be stressed that the displayed error bands do not include the uncertainty of the adopted spectral indices $`\gamma `$. Within the overall uncertainties the results obtained for different $`R_\mu `$ ranges are in fair agreement. This can be considered as a proof of the consistency of the Monte Carlo simulations invoked for the analysis of the data. The correction factors $`C_A`$ give the fraction of the primary energy spectrum of the considered primary mass which contributes to the $`\{N_\mu ^{tr},N_e,\mathrm{}\}`$ sample observed with conditions specified by the energy threshold of the registered muons $`E_{th}`$, the distance from the shower centre and the multiplicity threshold $`n_{th}`$. In order to give an impression on these factors, they are shown for various cases with their variation with $`n_{th}`$ in Fig. 7. The variation of $`C_A`$ is more pronounced at smaller $`N_\mu ^{tr}`$ values (i.e. primary energies), and as expected from the muon lateral distributions the acceptance of Fe events in the observation sample is higher than for proton events.

## 5 Test of the consistency of Monte Carlo simulations

As emphasised above, the event samples registered for these measurements are particular selections of all EAS events, determined by the observation conditions, especially by the energy threshold of the detected muons and the multiplicity threshold. Thus the event selections are affected by the lateral variation and mass dependence of the muon energy spectrum of EAS. The mass composition of these samples, as inferred from the correlations of the observed EAS observables, are altered as compared to the primary composition, in a way varying with the distance $`R_\mu `$ from shower centre. The acceptance or efficiency factors $`C_A`$ for correcting this effect depend on the Monte Carlo simulations and their ingredients like the model descriptions of the hadronic interaction. Since the measurements are performed at different distances from the EAS centre, the variation of the results after applying the $`C_A`$ factors implies a consistency test of the Monte Carlo simulations since the final result should be independent from $`R_\mu `$ within the given uncertainties. This test can be refined and be more stringent when not only the variation with $`R_\mu `$, but also the variation with the multiplicity threshold $`n_{th}`$ is scrutinised. In Fig. 8 the results of such a refined test are shown, based on the dominant $`\{N_\mu ^{tr},N_e\}`$ correlation within the observed samples. Observables describing the (local) EAS time structure (which are shown to be of minor sensitivity) have not been included in the classification procedure though their measurement defines the test sample. The results have been obtained by repeating the full data analysis, described in the previous sections, for different multiplicity thresholds. In view of the rather specific observation conditions, sensitive to variations of the EAS muon energy spectrum, we consider the results, even if displaying larger uncertainties, as a remarkable confirmation of the consistency of the Monte Carlo simulations performed with the CORSIKA code using the QGSJET model. The test could be improved when samples resulting from specific cuts would be analysed by simultaneously classifying the primary energy and the mass, using efficient pattern recognition methods as applied in other KASCADE studies . The mass composition resulting from the present studies

is in fair agreement with the results of those investigations. In particular it corroborates the finding of an increase of $`lnA`$ beyond the knee, i.e. for $`log_{10}N_\mu ^{tr}>`$ ca. 4.2.

## 6 Conclusions

The present investigation considers experimental data of the KASCADE experiment which are based on observations of local muon arrival time distributions, representing the variation of the EAS disc thickness with the distance from the shower axis. The first focus of the studies is the question to which extent these time quantities and their correlations with other EAS observables can help to improve the mass discrimination in the event samples observed with realistic experimental, i.e. with KASCADE conditions. Methodically the event-by-event analyses apply a non-parametric approach with reference patterns provided by CORSIKA simulations of the EAS development, using as generators the QGSJET model for the high energy interaction part and GHEISHA for the low energy part. It turns out that similarly to other EAS observables, the local time quantities provide only a marginal contribution to the mass discrimination as compared to the dominant $`\{N_\mu ^{tr},N_e\}`$ correlation, at least in the relatively limited range of $`R_\mu <100`$ m. It should be emphasised that this statement does by no means disqualify the sensitivity of global muon arrival time distributions measured relative to the arrival time of the shower centre. In fact, studies including the curvature of the shower front indicate a considerable enhancement of the mass discrimination, in particular at larger distances from the shower core and higher primary energies .
The event sample collected with the observation of muon arrival time distributions in the KASCADE experiment is a subset of all EAS events measured with $`100\%`$ efficiency for $`log_{10}N_\mu ^{tr}`$ approx. 3.6. As consequence of the muon energy threshold $`E_{th}=2.4`$ GeV and the condition of a registered multiplicity $`n3`$ the original primary mass composition in the registered samples gets distorted in a way which is sensitive to lateral variations of the integral EAS muon energy spectrum. In order to determine the fraction of all shower events which are accepted in the specific subset, specified by the measuring conditions, Monte Carlo simulations have to be invoked for the calculation of corresponding acceptance or efficiency factors. The variation of such corrections with the distance of observation $`R_\mu `$ from the shower axis and with the multiplicity threshold $`n_{th}`$ provides the possibility to test the consistency of the Monte Carlo simulations with the data. Even admitting the large uncertainties due to the limited number of events observed and simulated, our results indicate a remarkable consistency of the performed Monte Carlo simulations using the QGSJET model as generator.

We acknowledge various useful discussions with Prof. Dr. M. Petrovici about the aspects of the results elaborated in this paper. The authors would like to thank the members of the engineering and technical staff of the KASCADE collaboration who considerably contributed with enthusiasm and expert engagement to the success of the reported measurements. The work has been supported by the Ministry of Research of the Federal Republic of Germany. The Romanian Ministry of Education and Research provided a grant for supporting the collaborating group of the National Institute for Physics and Nuclear Engineering of Bucharest. The Polish collaborating group of the Cosmic Ray Division of the Soltan Institute for Nuclear Studies, Lodz, is supported by the Polish State Committee for Scientific Research (grant No. 5 P03B 133 20). The KASCADE collaboration work is embedded in the frame of scientific-technical co-operation (WTZ) projects between Germany and Romania (ROM 99/005), Poland (POL 99/005) and Armenia (ARM 98/002).


format: text

text:
# Anderson transition in low-dimensional disordered systems driven by nonrandom long-range hopping

thanks: On leave from “S.I. Vavilov State Optical Institute”, Saint-Petersburg, Russia.

## Abstract

The single-parameter scaling hypothesis predicts the absence of delocalized states for noninteracting quasiparticles in low-dimensional disordered systems. We show analytically and numerically that extended states may occur in the one- and two-dimensional Anderson model with a nonrandom hopping falling off as some power of the distance between sites. The different size scaling of the bare level spacing and the renormalized magnitude of the disorder seen by the quasiparticles finally results in the delocalization of states at one of the band edges of the quasiparticle energy spectrum. The delocalized nature of those eigenstates is investigated by numerical diagonalization of the Hamiltonian and by the supersymmetric method for disorder averaging, combined with a renormalization group analysis.

Localization of noninteracting quasiparticles in random media with time-reversal symmetry and finite-range hopping have been extensively studied since the seminal paper by Anderson Anderson58 . The hypothesis of single-parameter scaling, introduced in Ref. Abrahams79, , led to the general belief that all eigenstates of noninteracting quasiparticles were exponentially localized in one (1D) and two (2D) dimensions (see Refs. Lee85 ; Kramer93  for a comprehensive review) and that localization-delocalization transitions no longer exist in the thermodynamic limit. Despite models with finite-range hopping work nicely in describing a variety of materials, long-range hopping are often found in different physical systems (e.g. Frenkel excitons). Random long-range hopping was found to give rise to delocalization of states not only in three dimensional systems Anderson58  but in any dimension Logan87 ; Levitov89 ; Mirlin96 ; Parshin98 . Recent studies Deych00  revised the validity of the single-parameter scaling hypothesis even within the original 1D Anderson model with nearest-neighbor coupling, although did not question the statement that all eigenstates in 1D random systems are localized.

In this Letter we present analytical and numerical proofs that a localization-delocalization transition may occur in 1D and 2D systems with diagonal disorder and nonrandom intersite coupling which falls off according to a power-like law. Apart from the importance of this finding from a general point of view, it may be relevant for several physical systems. As an example, let us mention dipolar Frenkel excitons on 2D regular lattices where molecules are subjected to randomness due to a disordered environment Nabetani95 . Biological light-harvesting antenna systems represent a realization of the model we are dealing with Kopelman97 ; dendrimers . Magnons in 1D and 2D disordered spin systems provide one more example of interest.

We consider the Anderson Hamiltonian on a $`d`$-dimensional ($`d=1,2`$) simple lattice with $`𝒩=N^d`$ sites:

$$=\underset{𝐧}{}\epsilon _𝐧|𝐧𝐧|+\underset{\mathrm{𝐧𝐦}}{}J_{\mathrm{𝐧𝐦}}|𝐧𝐦|,$$
(1)
where $`|𝐧`$ is the ket-vector of the state localized at site $`𝐧`$, and $`\{\epsilon _𝐧\}`$ are random site energies, assumed to be uncorrelated for different sites and distributed uniformly within an interval $`[\mathrm{\Delta }/2,\mathrm{\Delta }/2]`$, thus having zero mean and standard deviation $`\sigma =\mathrm{\Delta }/\sqrt{12}`$. The hopping integrals between lattice sites $`𝐦`$ and $`𝐧`$ will be taken in the form $`J_{\mathrm{𝐦𝐧}}=J/|𝐦𝐧|^\mu (J_{\mathrm{𝐦𝐦}}0)`$, where $`J>0`$ without loss of generality and the lattice constant is set to unity. We stress that hopping integrals do not fluctuate.

For our qualitative reasoning we rewrite the Hamiltonian (1) in the Bloch wave representation, $`|𝐤=𝒩^{1/2}_𝐧\mathrm{exp}(i\mathrm{𝐤𝐧})|𝐧`$, with periodic boundary conditions. It then reads

$$=\underset{𝐤}{}E_𝐤|𝐤𝐤|+\underset{\mathrm{𝐤𝐤}^{}}{}(\delta )_{\mathrm{𝐤𝐤}^{}}|𝐤𝐤^{}|,$$
(2a)
$$E_𝐤=J\underset{𝐧\mathrm{𝟎}}{}\frac{e^{i\mathrm{𝐤𝐧}}}{|𝐧|^\mu },$$
(2b)
$$(\delta )_{\mathrm{𝐤𝐤}^{}}=\frac{1}{𝒩}\underset{𝐧}{}\epsilon _𝐧e^{i(𝐤𝐤^{})𝐧},$$
(2c)
where the wavenumbers $`𝐤`$ and $`𝐤^{}`$ run over the first Brillouin zone. Notice that $`\mu >d`$ to ensure the convergence of (2b) in the thermodynamic limit.

The key point of our qualitative arguments is as follows Rodriguez00 . We compare the size scaling of the typical magnitude of the scattering matrix $`(\delta )_{\mathrm{𝐤𝐤}^{}}`$ with the size scaling of the level spacing $`\delta E`$ in the bare quasiparticle spectrum of the homogeneous Hamiltonian ($`\mathrm{\Delta }=0`$). In particular, we focus our attention on those eigenstates laying close to the band edges $`𝐤=0`$ (top) and $`𝐤=𝝅`$ (bottom), where $`𝝅=\pi `$ and $`(\pi ,\pi )`$ for 1D and 2D systems, respectively. The typical fluctuation of the scattering matrix $`(\delta )_{\mathrm{𝐤𝐤}^{}}`$ is $`\sigma _{\mathrm{eff}}=\sigma /N^{d/2}`$. Thus, in spite of the fact that the magnitude of the disorder is $`\sigma `$, the quasiparticle sees an effectively reduced value $`\sigma _{\mathrm{eff}}`$. It is important that $`\sigma _{\mathrm{eff}}`$ scales inversely proportional to $`N^{d/2}`$. Straightforward calculations of the bare energy spectrum (2b) close to the band edges give the following results:

$$E_𝐤E_0JA_d(\mu )|𝐤|^{\mu d},|𝐤|0,$$
(3a)
$$E_𝐤E_𝝅+JB_d(\mu )|𝐤𝝅|^2,𝐤𝝅,$$
(3b)
where $`E_0`$ and $`E_𝝅`$ are the band-edge energies and $`A_d(\mu )`$ and $`B_d(\mu )`$ are dimensionless constants. From (3) it follows that the level spacing scales as $`\delta EN^{\mu +d}`$ at the top of the band, while at the bottom one gets $`\delta EN^2`$.

The matrix $`(\delta )_{\mathrm{𝐤𝐤}^{}}`$ couples the bare (extended) quasiparticle states to each other and may result in their localization within a region of size smaller than the system size. It seems reasonable to assume that the states will be weakly coupled and consequently will be delocalized over the whole system provided the inequality $`\sigma _{\mathrm{eff}}\delta E`$ holds. It is remarkable that for $`\mu <3d/2`$, the level spacing $`\delta E`$ at the top of the band diminishes upon increasing $`N`$ slower than the effective magnitude of disorder $`\sigma _{\mathrm{eff}}`$. Therefore, if the coupling between bare states is weak for some finite $`N`$ ($`\sigma _{\mathrm{eff}}\delta E`$) then it will become even weaker upon increasing $`N`$. Consequently, one may expect that the state will remain extended in the thermodynamic limit $`N\mathrm{}`$. It is also reasonable to assume that disorder of magnitude larger than the bare bandwidth will localize all the states. From the above arguments we conjectured the existence of an Anderson transition in 1D and 2D systems with diagonal disorder and nonrandom long-range hopping as long as $`\mu <3d/2`$. Below we provide analytical and numerical confirmations of this conjecture.

Concerning the parabolic range of the energy spectrum (close to the bottom of the band), we notice that the level spacing diminishes as $`N^2`$ upon increasing the lattice size, i.e., faster than the effective magnitude of disorder $`\sigma _{\mathrm{eff}}`$. Now, even if $`\sigma _{\mathrm{eff}}\delta E`$ for a small lattice size and the states are delocalized, the above inequality will be reverted for larger $`N`$, resulting finally in the localization of those states. The same conclusion holds for both band edges within the nearest-neighbor approximation, where the level spacing is always $`\delta EN^2`$.

A supersymmetric method for disorder averaging efetov ; gurus , combined with a renormalization group (RG) analysis, provide support to the above arguments. In short (the details will be published elsewhere), the sequence of our steps is as follows. As a first step, we consider the one particle Green’s function with the fermionic partition function $`Z_0`$ and the bare action $`S_0`$ defined as

$$Z_0=\underset{𝐧}{}d\psi _𝐧d\overline{\psi }_𝐧e^{S_0}=\mathrm{det}(I),$$
(4a)
$$S_0=\underset{\mathrm{𝐧𝐦}}{}i\overline{\psi }_𝐧(I)_{\mathrm{𝐧𝐦}}\psi _𝐦,$$
(4b)
where $`=E+i0^+`$ and $`I`$ is the identity matrix. Introducing bosonic ghosts $`\beta ,\overline{\beta }`$ and expressing $`1/Z_0`$ as a path integral gurus , we then average the one particle Green function using the Gaussian probability distribution of site energies, $`P(\epsilon _𝐧)=(1/\pi g)^{1/2}\mathrm{exp}(\epsilon _𝐧^2/g)`$, instead of the box distribution introduced in the beginning. This allows us to perform the integration over site energies explicitly. The effective action

$`S_{\mathrm{eff}}`$ $`=`$ $`i{\displaystyle \underset{\mathrm{𝐧𝐦}}{}}[\overline{\psi }_𝐧(J_{\mathrm{𝐧𝐦}}\delta _{\mathrm{𝐧𝐦}})\psi _𝐦`$ (5)
$`+`$ $`\overline{\beta }_𝐧(J_{\mathrm{𝐧𝐦}}\delta _{\mathrm{𝐧𝐦}})\beta _𝐦]`$
$`+`$ $`{\displaystyle \frac{g}{4}}{\displaystyle \underset{𝐧}{}}(\overline{\psi }_𝐧\psi _𝐧+\overline{\beta }_𝐧\beta _𝐧)^2,`$
which appears after averaging, will be the main object of our RG analysis. For doing this, it is convenient to rewrite the action in the $`𝐤`$-representation and to regroup the terms as follows $`S_{\mathrm{eff}}=S_{\mathrm{kin}}+S_{}+S_g`$, where

$`S_{\mathrm{kin}}=`$ $``$ $`iJA_d{\displaystyle }d^d𝐤|𝐤|^{\mu d}[\overline{\psi }(𝐤)\psi (𝐤)`$ (6a)
$`+`$ $`\overline{\beta }(𝐤)\beta (𝐤)],`$
$$S_{}=id^d𝐤[\overline{\psi }(𝐤)\psi (𝐤)+\overline{\beta }(𝐤)\beta (𝐤)],$$
(6b)
$`S_g`$ $`=`$ $`{\displaystyle \frac{g}{8\pi }}{\displaystyle \underset{i=1}{\overset{4}{}}d^d𝐤_𝐢\delta (𝐤_\mathrm{𝟏}+𝐤_\mathrm{𝟐}𝐤_\mathrm{𝟑}𝐤_\mathrm{𝟒})}`$ (6c)
$`\times [\overline{\psi }(𝐤_\mathrm{𝟏})\psi (𝐤_\mathrm{𝟐})+\overline{\beta }(𝐤_\mathrm{𝟏})\beta (𝐤_\mathrm{𝟐})]`$
$`\times [\overline{\psi }(𝐤_\mathrm{𝟑})\psi (𝐤_\mathrm{𝟒})+\overline{\beta }(𝐤_\mathrm{𝟑})\beta (𝐤_\mathrm{𝟒})].`$
We have absorbed the constant $`E_0`$ into $``$ in (6b). Note as well that the integration over momenta are restricted to a $`d`$ dimensional sphere of radius $`\mathrm{\Lambda }`$, i.e. $`|𝐤|<\mathrm{\Lambda }`$, with $`\mathrm{\Lambda }`$ being an ultraviolet cutoff.

The action (6) is the starting point of our RG analysis, which is inspired by Shankar’s approach shankar  to fermionic condensed matter systems. The key observation is that the kinetic part of the action, $`S_{\mathrm{kin}}`$, is invariant under the following scaling transformation of the cutoff $`\mathrm{\Lambda }`$, the momenta $`𝐤`$ and the fields:

$$\mathrm{\Lambda }\mathrm{\Lambda }^{}=\mathrm{\Lambda }/b,𝐤𝐤^{}=b𝐤,b>1$$
(7a)
$$\varphi (𝐤^{}/b)=b^{\mu /2}\varphi (𝐤^{}),\varphi =\psi ,\overline{\psi },\beta ,\overline{\beta }.$$
(7b)
For generic values of $`\mu >d`$, Eq. (7a) is a non standard scaling law which emerges from the unusual kinetic term (6a).

Driven by (7), the mass term $`S_{}`$ also transforms into itself with a new coupling constant $`^{}`$ given by $`^{}=b^{\mu d}`$. Thus $``$ is a relevant perturbation of the free action $`S_{\mathrm{kin}}`$, as it is always the case of mass terms shankar .

The term $`S_g`$ also transforms onto itself under the following RG transformation $`g^{}=b^{2\mu 3d}g`$. This equation implies that $`g`$ eventually goes to zero upon increasing $`b`$ provided that $`\mu <3d/2`$. Hence, randomness vanishes in the low energy effective theory for $`\mu <3d/2`$. On the contrary, $`g`$ runs to stronger coupling whenever $`\mu >3d/2`$.

Finally, the coupling $`g`$ is marginal at tree level for $`\mu =3d/2`$, and one has to consider the one loop effects to see its fate. This can be done using the techniques developed in Ref. shankar . Let us present our main results. The RG flows of the constants $`g`$ and $``$, up to one loop, are given by

$`{\displaystyle \frac{d\overline{g}}{ds}}`$ $`=`$ $`(2\mu 3d)\overline{g}+\overline{g}^2,`$ (8a)
$`{\displaystyle \frac{d\overline{}}{ds}}`$ $`=`$ $`(\mu d)\overline{}\overline{g},`$ (8b)
where $`s`$ is the RG parameter defined as $`b=\mathrm{exp}(s)`$, and $`\overline{g}`$ and $`\overline{}`$ are related to $`g`$ and $``$ as follows

$$g=\frac{\pi }{\mathrm{\Omega }_d}(JA_d)^2\mathrm{\Lambda }^{2\mu 3d}\overline{g},=\frac{1}{4}JA_d\mathrm{\Lambda }^{\mu d}\overline{},$$
(9)
with $`\mathrm{\Omega }_d`$ being the volume of the $`d`$-dimensional sphere.

Equation (8a) has an unstable fixed point $`\overline{g}_{}=3d2\mu `$ provided that $`\mu <3d/2`$ (see Fig. 1a). Below this point, the coupling $`g`$ goes to zero asymptotically, while above $`g`$ grows. For $`\mu >3d/2`$ the critical point disappears and the system always flows to strong coupling (see Fig. 1b). This signals about some changes in the density of states and the different nature of the eigenfunctions when passing from $`\mu <3d/2`$ to $`\mu >3d/2`$. We would like to remark the fact that the critical value $`\mu =3d/2`$ appeared in the present RG analysis coincides with that found on the basis of our qualitative arguments.

It is usually argued that the averaged one-particle Green function does not carry information about the spatial extent of the eigenfunctions. Therefore, we chose the generalized inverse participation ratios (GIPR) as the relevant quantities to characterize localization properties of the states. The standard definition reads $`I_\nu ^{(q)}=_𝐧|\mathrm{\Psi }_{\nu 𝐧}|^{2q}`$, where $`\mathrm{\Psi }_{\nu 𝐧}`$ is the probability amplitude of the normalized eigenstate $`\nu `$ at site $`𝐧`$, and the brackets denote disorder averaging. As is well known, the GIPR scales as $`I_\nu ^{(q)}𝒩^{q+1}`$ for delocalized states ($`q1`$), while shows no size scaling for localized ones.

Our main steps in calculating the GIPR are as follows. We estimated the uppermost eigenfunction, $`\mathrm{\Psi }_{0𝐧}`$, using perturbation theory with respect to the random term $`_{\mathrm{𝐤𝐤}^{}}(\delta )_{\mathrm{𝐤𝐤}^{}}|𝐤𝐤^{}|`$ in (2a), namely considering the on-site energies small compared to the hopping parameter $`J`$, $`|\epsilon _𝐧|J`$. Then, $`\mathrm{\Psi }_{0𝐧}`$ can be written in real space as

$`\mathrm{\Psi }_{0𝐧}`$ $`=`$ $`{\displaystyle \frac{e^{\varphi _𝐧}}{(_𝐧e^{2\varphi _𝐧})^{1/2}}},`$ (10a)
$`\varphi _𝐧`$ $``$ $`{\displaystyle \underset{𝐦}{}}S_{\mathrm{𝐧𝐦}}\epsilon _𝐦,`$ (10b)
$`S_{\mathrm{𝐧𝐦}}`$ $``$ $`{\displaystyle \frac{1}{𝒩}}{\displaystyle \underset{𝐤\mathrm{𝟎}}{}}{\displaystyle \frac{e^{i𝐤(𝐦𝐧)}}{E_0E_𝐤}}.`$ (10c)
In order to compute the GIPR we again made use of the the supersymmetric method for disorder averaging (as described above) with the bare action given by Eq. (5), as well as the replica trick introduced in Ref. ludwig, . The latter reads

$$I_0^{(q)}=\underset{r0}{lim}\frac{1}{N^d}\underset{𝐧,𝐧_1,\mathrm{},𝐧_{rq}}{}e^{2q\varphi _𝐧}\underset{j=1}{\overset{rq}{}}e^{2\varphi _{𝐧_j}}.$$
(11)
In doing so, we found that $`I_0^{(q)}N^{(q1)d}`$ provided when $`d<\mu <3d/2`$; in other words, the generalized dimension equals the space dimension so that the uppermost state is delocalized, in full agreement with our qualitative picture. In particular, notice that the so called inverse participation ratio (IPR) scales as $`I_0^{(2)}𝒩^1`$.

Since the previous analytical study of the GIPR was perturbative, we have also carried out a numerical study of the model to support the validity of our conclusions. We took advantage of the Lanczos method Golub96  as well as the density matrix renormalization group approach Delgado99 , allowing one to calculate some few eigenstates of the Hamiltonian (1) for rather large system size. In Fig. 2 we plotted the IPR ($`q=2`$) of the uppermost state as a function of the system size $`𝒩=N^d`$ for different degrees of disorder $`\mathrm{\Delta }`$. The behaviour of the other top states is similar to that which manifests the uppermost state. Observing Fig. 2 we conclude that the uppermost state is delocalized even for a moderately high value of the degree of disorder ($`\mathrm{\Delta }=8J`$ in 1D systems and $`\mathrm{\Delta }=40J`$ in 2D system), provided $`d<\mu <3d/2`$. For comparison, the 1D (2D) bandwidth for $`\mu =5/4`$ ($`\mu =9/4`$) in the absence of disorder is of the order of $`10.5J`$ ($`28J`$). However, for large degree of disorder the IPR remains constant on increasing the system size (see Fig. 2 for $`\mu =5/4`$, $`d=1`$ and $`\mathrm{\Delta }=40J`$), indicating that the uppermost eigenstate is localized. Therefore, the top eigenstates undergo the Anderson transition on increasing $`\mathrm{\Delta }`$ whenever $`d<\mu <3d/2`$. It is to be noticed the absence of scaling of the IPR and the subsequent localization for $`d=1`$, $`\mu =3`$, and $`\mathrm{\Delta }=8J`$. This result is in full correspondence with the analytical analysis stating that no transition is expected in this case.

In summary, we have shown analytically and numerically that a power-like nonrandom intersite hopping, $`J_{\mathrm{𝐧𝐦}}=J/|𝐧𝐦|^\mu `$, may act towards delocalization of quasiparticle states in low dimensional systems. In particular, the states of the top of the band may be delocalized at rather high magnitudes of disorder ($`\mathrm{\Delta }J`$) and undergo a localization-delocalization transition as the magnitude of disorder increases. Scaling arguments provide a clear physical picture of the underlying mechanism responsible for the delocalization of the states, in spite of the low-dimensional ($`d2`$) geometry of the system. The different size scaling of the bare level spacing, $`\delta EN^{\mu +d}`$, and the magnitude of disorder seen by a quasiparticle, $`\sigma _{\mathrm{eff}}\mathrm{\Delta }/N^{d/2}`$, is the feature of the model responsible for this unusual behavior. $`\delta E`$ decreases slower than $`\sigma _{\mathrm{eff}}`$ upon increasing the system size as long as $`d<\mu <3d/2`$, resulting in the delocalization of the corresponding quasiparticle states in the thermodynamic limit. We stress that the main finding of our study, namely the existence of the Anderson transition in a physically relevant model, has been concluded on the basis of three different approaches and the conclusions obtained are self-consistent. Most important, the validity of the scaling analysis is not limited to the present model. Indeed, it is established on solid grounds that the standard, three dimensional Anderson model manifests the localization-delocalization transition at the band center. Within this model, the bare level spacing at the band center diminishes proportionally to $`N^1`$, while the magnitude of effective disorder goes down faster, $`N^{3/2}`$, thus being unable to localize the states at the band center for a moderate disorder. A strong disorder (large compared to the band width) localizes the states, giving rise a localization-delocalization transition. We can then be confident that this kind of scaling arguments may provide physical insight in several localization problems.

###### Acknowledgements.

V. A. M. acknowledges support from MECyD (Project SAB2000-0103). A. R. and F. D-A. were supported by DGI-MCyT (Project MAT2000-0734) and CAM (Project 07N/0075/2001). G. S. and M. A. M-D. acknowledge support from PGC (Project BFM2000-1320-C02-01).


text:
# Deep optical observations of G 65.3+5.7

## 1 Introduction

Gull et al. (gul77 (1977)) reported the detection of a new supernova remnant in Cygnus during an optical survey with narrow band filters. In particular, they obtained an $`[`$O iii$`]`$ photograph showing the existence of a $``$3°.3$`\times `$4°.0 filamentary shell. Due to its large extent, high resolution imagery has been performed only in selected regions of the remnant (Fesen et al. fes83 (1983) and references therein). Subsequent, optical spectra obtained by Sabbadin & D’Odorico (sab76 (1976)) and Fesen et al. (fes85 (1985)) suggested emission from shock heated gas. A stratification among lower and medium ionization lines has been reported by Sitnik et al. (sit83 (1983)) for another area of G 65.3+5.7. Fesen et al. (fes85 (1985)) obtained long–slit spectra at two different positions of G 65.3+5.7 showing strong sulfur line emission relative to H$`\alpha `$ and especially one of the spectra was characterized by very strong $`[`$O iii$`]`$ emission compared to H$`\beta `$ ($``$40). High $`[`$O iii$`]`$/H$`\beta `$ ratios have been observed in other remnants as well and imply the presence of incomplete recombination zones (e.g. Raymond et al. ray88 (1988)). Radio observations at 1420 MHz (Reich et al. rei79 (1979)) confirmed the non-thermal nature of the radio emission from this shell and found a spectral index of $``$0.61. Current estimates of the linear diameter and the distance to the remnant are $``$70 pc and 0.9–1.0 kpc, respectively (Gull et al. gul77 (1977); Reich et al. rei79 (1979); Rosado ros81 (1981); Sitnik et al. sit83 (1983)). Soft X–ray observations show that G 65.3+5.7 is a weak source of X–ray emission and imply a shock velocity in the range of a few hundred km s-1 and an age of $``$ 20000 yrs (Snyder et al. sny78 (1978); Mason et al. mas79 (1979); Seward sew90 (1990)). G 65.3+5.7 has also been observed by ROSAT both during the All-Sky survey (Aschenbach asc94 (1994)) and in pointed mode (Lu & Aschenbach lu02 (2002)). The analysis of these data showed that the X–ray emission is clumpy and that the forward shock is traveling with a velocity of $``$400 km s-1 about 26000 yrs after the supernova explosion, while the explosion energy and the interstellar medium density are $``$2$`\times `$1050 erg and $``$0.02 cm-3, respectively. These parameters were derived assuming a distance to the remnant of 1 kpc. Gorham et al. (gor96 (1996)) have performed a pulsar search of several known SNRs at 430 and 1420 MHz. Their observations were sensitive to pulsars with periods greater than 1 ms, and flux densities as low as 0.2 mJy but no pulsar was found to be associated with G 65.3+5.7 down to these limits.

In order to study in detail this extended and nearby supernova remnant we performed deep optical imaging observations in the low and medium ionization lines of $`[`$O ii$`]`$3727 Å and $`[`$O iii$`]`$5007 Å. We have also obtained flux calibrated images of a specific field of G 65.3+5.7 in H$`\alpha +[`$N ii$`]`$, $`[`$S ii$`]`$, $`[`$O ii$`]`$ and $`[`$O iii$`]`$. The imaging observations are supplemented by deep optical long–slit spectra at 8 different locations of G 65.3+5.7. Information about the observations and the data reduction is given in Sect. 2. In Sect. 3, 4 and 5 we present the results of our imaging observations and the long–slit spectra. Finally, in Sect. 6 we discuss the physical properties of the remnant.

## 2 Observations

### 2.1 Optical images

The observations presented in this paper were performed with the 0.3 m Schmidt Cassegrain telescope at Skinakas Observatory, Crete, Greece. The supernova remnant G~65.3+5.7 was observed on July 22, 23 and 24, 2001. A 1024 $`\times `$ 1024 SITe CCD was used which had a pixel size of 24 $`\mu `$m resulting in a 89′ $`\times `$ 89′ field of view and an image scale of 5″ per pixel. Several pointings were performed in order to cover the extended area of G 65.3+5.7. Each of the nine (9) different fields was observed for 2400 s. All fields were projected to a common origin on the sky and were subsequently combined to create the final mosaics in $`[`$O ii$`]`$ and $`[`$O iii$`]`$. The astrometric calculations utilized the HST Guide star catalog (Lasker et al. las99 (1999)). In addition, we obtained flux calibrated images of a specific field of G 65.3+5.7 in H$`\alpha +[`$N ii$`]`$, $`[`$S ii$`]`$, $`[`$O ii$`]`$ and $`[`$O iii$`]`$. The exposure time of these images is 1800 s with the exception of the $`[`$O ii$`]`$ image where we obtained two images, 1800 s each. This observation was performed in July 13, 1999 with the same telescope and hardware setup as the latest observations. The filter characteristics can be found in Mavromatakis et al. (mav00 (2000)). All coordinates quoted in this work refer to epoch 2000.

We employed standard IRAF and MIDAS routines for the reduction of the data. Individual frames were bias subtracted and flat-field corrected using well exposed twilight flat-fields. The spectrophotometric standard stars HR5501, HR8634 and HR7950 were used for absolute flux calibration in the 1999 observation of G 65.3+5.7.

### 2.2 Optical spectra

Long–slit spectra were obtained on July 15, 1999, and July 24–28, 2001 using the 1.3 m Ritchey–Cretien telescope at Skinakas Observatory. The data were taken with a 1300 line mm-1 grating and a 800 $`\times `$ 2000 Site CCD covering the range of 4750 Å – 6815 Å. The slit had a width of 7$`\stackrel{}{.}`$7 and a length of 7′.9 and, in all cases, was oriented in the south-north direction. The coordinates of the slit centers along with the number of spectra and their individual exposure times are given in Table 1. The spectrophotometric standard stars HR5501, HR7596, HR7950, HR9087, and HR718 were observed in order to calibrate the spectra.

## 3 Imaging of G 65.3+5.7

### 3.1 The $`[`$O ii$`]`$ and $`[`$O iii$`]`$ line mosaics

The wide field covered by G 65.3+5.7 contains a number of extended sources like bright nebulae (e.g. LBN 064.90+06.74; Lynds lyn65 (1965)), dark nebulae (e.g. LDN 817, LDN 816; Lynds lyn65 (1965)) as well as a few planetary nebulae (e.g. PN G064.7+05.0, PN G065.4+03.1). In Fig. 1 we show the mosaic in the low ionization line of $`[`$O ii$`]`$ along with letters designating the areas which are described in the text in detail. The $`[`$O ii$`]`$ image appears generally less filamentary than the $`[`$O iii$`]`$ image (Fig. 2). In addition, the former image contains more diffuse emission than the latter. We start the detailed morphological study of G 65.3+5.7 in these filters from area A (Fig. 1). In this area diffuse emission is seen in the $`[`$O ii$`]`$ image, while strong filamentary emission is detected in the medium ionization line image of $`[`$O iii$`]`$. The filament is $``$16′ long and its typical projected width is $``$1′. Immediately to the south–west of area A we find area B which contains the bright binary $ϕ$ Cyg and has been studied by Sitnik et al. (sit83 (1983)). They observed a narrow field around $`\varphi `$ Cyg ($``$6′.5 $`\times `$ 6′.5) and found a stratification in the emission from the $`[`$O iii$`]`$, $`[`$N ii$`]`$ and $`[`$S ii$`]`$ lines with the $`[`$O iii$`]`$ emission lying closer to the outer edge of the remnant. Here, the $`[`$O ii$`]`$ emission appears more filamentary than the $`[`$O iii$`]`$ emission. The $``$37′ long filament seen in $`[`$O ii$`]`$ is broader in the south and becomes narrower as it extends further to the north, while the $`[`$O iii$`]`$ image shows roughly a constant projected thickness. The projected thickness close to the northern tip of the $`[`$O ii$`]`$ filament is $``$32″, while at the same location the thickness of the $`[`$O iii$`]`$ filament is $``$65″. We also observe here an increasing offset in the outer boundaries of the $`[`$O ii$`]`$ and $`[`$O iii$`]`$ images with the $`[`$O iii$`]`$ emission leading and the $`[`$O ii$`]`$ emission following. Further to the west we come across area C for which we have performed flux calibrated imaging observations. This area will be discussed in the next section separately.

Significant differences between the low and medium ionization images are found in area D. There are two locations where the $`[`$O ii$`]`$ emission is quite intense and filamentary, while at the same location the $`[`$O iii$`]`$ emission is diffuse and weak. In this area there is also a location where the inverse is observed, i.e. filamentary emission in $`[`$O iii$`]`$ and diffuse, very weak emission in $`[`$O ii$`]`$. Area E contains a complex network of filaments and diffuse emission in $`[`$O ii$`]`$, while the $`[`$O iii$`]`$ image provides a sharper view of the area. LBN 148 (Lynds lyn65 (1965)) is also found in this area, and shows up as a 25′ long filament in $`[`$O ii$`]`$, while it appears much shorter in the $`[`$O iii$`]`$ emission line image ($``$18′). Emission in both images is detected further to the west of LBN 148 and seems to extend outside our field. We note here that our coverage to the west marginally reaches $`\alpha `$$``$19h24m, while the ROSAT data show that emission extends up to $`\alpha `$$``$19h23m (Lu & Aschenbach lu02 (2002)).

In area F we find a long filamentary structure both in $`[`$O ii$`]`$ and $`[`$O iii$`]`$ (Figs. 1 and 2) extending from the west edge of our field up to the central areas of the remnant. There is a wealth of short scale structures to the north of this long filament which are only seen in the low ionization image of $`[`$O ii$`]`$. However, further to the north in the $`[`$O iii$`]`$ image faint, filamentary and diffuse structures are discovered, around $`\alpha `$$``$19h27m and $`\delta `$$``$33°00′, which seem to extend outside the north boundary of our field and are projected outside the main body of the soft X–ray emission (Lu & Aschenbach lu02 (2002)). Based on the extent of the remnant as seen in the radio and X–ray observations it is likely that this emission is not part of G 65.3+5.7. Substantial differences between the two mosaics are also seen in area G. Although the nature of the $`[`$O iii$`]`$ structure observed is filamentary all along the $``$42′ of its length, the $`[`$O ii$`]`$ emission appears quite faint and completely diffuse. Another area with great differences among the low and medium ionization images is area H. The structures are much better defined in $`[`$O iii$`]`$ than in $`[`$O ii$`]`$. The $`[`$O ii$`]`$ image of the area seems to consist of several small scale filamentary structures and some of them do not have $`[`$O iii$`]`$ counterparts. Finally, in area I we find that the $`[`$O iii$`]`$ emission is better structured compared to the $`[`$O ii$`]`$ emission line image which, in addition, looks more clumpy. Note that around $`\alpha `$ $``$ 19h41m and $`\delta `$ $``$ 32°20′ faint $`[`$O iii$`]`$ emission is detected in an arc–like shape. Its typical brightness is $``$1/3 of the bulk of the emission seen a few arcminutes to the south.

## 4 Area C in detail

Long observing runs would be required in order to map this remnant in the major optical emission lines. The only area observed in more detail is area C (Fig. 3) which was observed through the lines of H$`\alpha +[`$N ii$`]`$, $`[`$S ii$`]`$, $`[`$O ii$`]`$ and $`[`$O iii$`]`$. These observations are flux calibrated and show that the emission, in all filters, is strong in absolute terms. However, there are striking differences between the low ionization lines and $`[`$O iii$`]`$ also in this area. The major difference being the double filament seen in the center of this field (Fig. 3). It is best defined in the medium ionization line of $`[`$O iii$`]`$, while diffuse emission is present in the $`[`$O ii$`]`$ image but is hardly detected in the images of H$`\alpha +[`$N ii$`]`$ and $`[`$S ii$`]`$. The two filaments are separated by $``$1′, while their typical projected thickness is less than $``$35″. This angular thickness corresponds to 0.17 pc, at a distance of 1 kpc, which is much smaller than the total length of this structure of $``$19 pc. We note here that only a small part ($``$16′) of this filamentary structure is present in the low ionization images strongly suggesting the absence of complete recombination zones. This seems to be a common characteristic of the various areas of G 65.3+5.7. In order to explore the properties of the double filament, we formed the ratio of the $`[`$O ii$`]`$ to the $`[`$O iii$`]`$ image. This new image could help to discriminate between projection effects and actual changes in the shock parameters. We find that the measured values of the ratio along each filament change by $``$20%, while the average values of the ratios at the two filaments agree to within the 1.5$`\sigma `$ range. The relatively small difference in the ratio between the two filaments might point to simple projection effects but images of higher spatial resolution would be required to check this issue in more detail.

Finally, we also formed the ratio of the $`[`$O iii$`]`$ emission line image to the H$`\alpha +[`$N ii$`]`$ image which would allow to establish areas with complete and incomplete recombinations zones (Raymond et al. ray88 (1988)). Indeed, this image clearly demonstrates that the double filament as well as the strong emission seen in the south–west are characterized by incomplete structures, while the emission from the south–east regions originates from complete shock structures.

## 5 The long–slit spectra from G 65.3+5.7

We have obtained long–slit spectra at eight different locations of G 65.3+5.7. The absolute H$`\alpha `$ flux covers a wide range of values from 0.5 to 22 $`\times `$ $`10^{17}`$ erg s-1 cm-2 arcsec-2 (Table 2). All spectra suggest that the emission detected during these spectroscopic observations originates from shock heated gas ($`[`$S ii$`]`$/H$`\alpha `$ $``$0.7-1.8; Table 2). This conclusion also applies to the single spectrum taken from LBN 148 showing that this object is related to the remnant and is not an H ii region. The ratio of the sulfur lines emitted at 6716 Å and 6731 Å allow us to estimate the electron density (e.g. Osterbrock ost89 (1989)) at the position observed. The electron densities were calculated with the nebular package within the IRAF software (Shaw and Dufour sha95 (1995)). All ratios tend to the high end of the allowable range of values suggesting electron densities below $``$150 cm-3 and even below $``$30 cm-3. Consequently, we cannot directly measure the preshock cloud density but can only place upper limits. The H$`\alpha `$/H$`\beta `$ ratio is used to estimate the logarithmic extinction c towards a source of line emission assuming an intrinsic ratio of 3. A signal to noise weighted average of the logarithmic extinction c is 4.3 for the different locations around G 65.3+5.7. Higher values are found among the spectra but are also accompanied by larger errors. Fesen et al. (fes85 (1985)) obtained spectra at two locations and found H$`\alpha `$/H$`\beta `$ ratios around 3.3 ($`\pm `$ 35%) in agreement with our measurements within the 1$`\sigma `$ range.

The most interesting result involves the very large $`[`$O iii$`]`$/H$`\beta `$ ratios observed in most of the spectra. Models of complete shock structures predict values below $``$6 for this ratio (Cox & Raymond cox85 (1985), Hartigan et al. har87 (1987)). This limit is easily exceeded in case of shocks with incomplete recombination zones (Raymond et al. ray88 (1988)). Our measured values for the $`[`$O iii$`]`$/H$`\beta `$ ratio range from 6 to 110 (Table 2) indicating that shocks with complete and incomplete recombination zones are present. Fesen et al. (fes85 (1985)) measured this ratio to be $``$2 (their position 1) and $``$40 (their position 2). Our area 1 lies close to their position 2 ($``$90″apart) and there the measured ratios of $`[`$O iii$`]`$/H$`\beta `$ are 46 and 86 along the slit (Table 2). These spectra were taken at different positions and thus, cannot be directly compared but they do point to the presence of incomplete recombination zones in this part of the remnant. We note here that the signal to noise ratios quoted in Table 2 do not include calibration errors which are less than $``$10%.

As is evident from Table 2 only a small number of the extracted spectra points to shocks with complete recombination zones. The spectra from positions 1, 2, 3a, 4, 6, 7, and 8 strongly support the presence of incomplete recombination zones ($`[`$O iii$`]`$/H$`\beta `$ $``$10–100). Examination of the theoretical estimates of the $`[`$O iii$`]`$/H$`\beta `$ and $`[`$O i$`]`$/H$`\beta `$ ratios as a function of column density (Raymond et al. ray88 (1988)) shows that shock speeds higher than $``$100 km s-1, most likely around 120 km s-1, and column densities in the range of 1017 – 1018.5 cm-2 are required to match the observations. Especially, at position 1 a shock velocity of 120 km s-1 and a column density of $``$1017.5 cm-2 is estimated, which is well below that of $``$1018.7 cm-2 required for a complete recombination zone.

## 6 Discussion

The extended remnant G 65.3+5.7 was observed for the first time with a highly efficient CCD camera in the low and medium ionization lines of $`[`$O ii$`]`$3727 Å and $`[`$O iii$`]`$5007 Å, respectively. The $`[`$O ii$`]`$ interference filter being in the far blue part of the optical wavelength band provides sharp images since the emission from field stars is significantly reduced, while any morphological differences between the two images can not be attributed to abundance variations. In addition, we performed the first flux calibrated imaging observations of a specific area in the south of G 65.3+5.7 (area C, §4) in major optical emission lines. These observations as well as the raster observations reveal in sufficient detail what is already known for parts of the remnant, i.e. strong differences between the lower and higher ionization lines (Rosado ros81 (1981), Fesen et al. fes83 (1983), Sitnik et al. sit83 (1983)).

The object under study is an excellent example of the class of remnants displaying prominent $`[`$O iii$`]`$ emission and is this line that provides the sharpest view of the system (e.g. Fesen et al. fes97 (1997), Mavromatakis et al. mav02 (2002)). In the following we will assume a distance of 1 kpc to the remnant and an angular radius of 1°.8 implied by the X–ray (Lu & Aschenbach lu02 (2002)) and our optical data. The thinnest filamentary structures seen in the $`[`$O ii$`]`$ images are characterized by typical projected widths, actually the full width at half maximum (fwhm), of $``$20″ or 0.1 pc, while the shortest filaments are $``$0.5 pc long. The longest ones are of the order of several pc. We also find gaps in certain filaments which otherwise look continuous. The shortest gaps are typically around 0.2 pc suggesting this length as the typical length of small scale inhomogeneities in this cloudy environment.

A shock velocity of $``$120 km s-1 and a column density of $``$1017.5 –1018 cm-2 are estimated for area 1 from the optical spectra and the modeling results of Raymond et al. (ray88 (1988)). The measured fwhm of 0.14 pc of the filament in this area and the estimated column density allow us to estimate an average filament density of $``$0.8–2.5 cm-3. The actual preshock cloud density must be lower than this by, at least, a factor of 4. An IUE spectrum taken a few arcseconds off our position shows emission from C iii\], N iii\], O iii\], but weak C iv and maybe Si iv, and O iv (Raymond J. 2002, private communication). Thus, it is not possible to fully explore the physical conditions in this area as Raymond et al. (ray88 (1988)) have done for a bright filament in the Cygnus Loop. Using the same approach to other areas of the remnant we find average densities in the range of $``$0.3–4 cm-3. Specifically, for area 5 a shock velocity of 90 km s-1 and column densities in the range of 1018 – 1018.5 cm-2 are estimated. Extraction of the intensity profiles from the $`[`$O ii$`]`$ and $`[`$O iii$`]`$ images along the declination axis shows that the positions of maximum intensity differ by $``$0.07 pc with the $`[`$O iii$`]`$ profile leading, while the fwhm of the fitted Gaussian profiles are 0.13 and 0.19 pc, respectively. This behaviour of the lower and higher ionization lines is also verified by the long slit spectra which are characterized by a spatial scale of 1/350 pc compared to that of 1/40 pc of the wide field images. Assuming a typical length of 0.19 pc and a column density of $``$1018.3 cm-2, we estimate the average density behind the radiative shock as $``$3.5 cm-3. These simplified calculations show that in several areas of the remnant the particle density in the recombination region is of the order of a few atoms per cm3.

In areas where complete recombinations zones are suggested by the observational data we can utilize the relations

$$\mathrm{n}_{[\mathrm{S}\text{ii}]}45\mathrm{n}_\mathrm{c}\mathrm{V}_\mathrm{s}^2,$$
(1)
and

$$\mathrm{E}_{51}=2\times 10^5\beta ^1\mathrm{n}_\mathrm{c}\mathrm{V}_\mathrm{s}^2\mathrm{r}_\mathrm{s}^3,$$
(2)
given by Fesen & Kirshner (fes80 (1980)) and McKee & Cowie (mck75 (1975)), respectively, to limit the preshock cloud density and explosion energy. The factor $`\beta `$ is of the order of 1–2, $`\mathrm{n}_{[\mathrm{S}\text{ii}]}`$ is the electron density derived from the sulfur line ratio, nc is the preshock cloud density, $`\mathrm{E}_{51}`$ is the explosion energy in units of 1051 erg, Vs is the shock velocity into the clouds in units of 100 km s-1, and rs is the radius of the remnant in pc. Assuming an upper limit of $``$70 cm-3 on the electron density (Table 2; areas 3b, 4, 5) we find that $`\mathrm{n}_\mathrm{c}\mathrm{V}_\mathrm{s}^2`$ $`<`$ 1.6 and consequently, $`\mathrm{E}_{51}`$ $`<`$ 1 for rs $``$ 32 and $`\beta `$ 1. Although this limit is rather broad, it is consistent with the explosion energy of 0.2 estimated by Lu & Aschenbach (lu02 (2002)) from the Sedov–Taylor solution of the X–ray data. However, in order to derive Eqs. (1) and (2) it is assumed that the magnetic field component perpendicular to the shock velocity vector is negligible (see also Bocchino et al. boc00 (2000)). If this is not the case then the preshock cloud density will be expected to be higher than what Eq. (1) implies. Comparison of the shock velocities estimated from the optical data with the shock velocity of the primary blast wave suggests density contrasts $``$10–20 between the interstellar medium and the interstellar clouds. The preshock medium density of 0.02 cm-3 given by Lu & Aschenbach (lu02 (2002)) implies typical preshock cloud densities $``$0.2–0.4 cm-3, which are consistent with the upper limits obtained earlier.

The observed (H$`\alpha `$/H$`\beta `$)obs ratio can be used to estimate the logarithmic extinction c towards a source of line emission. The relation we use is c = 1./0.331 * log((H$`\alpha `$/H$`\beta `$)obs/3.) and is based on the interstellar extinction curve of Kaler (kal76 (1976)) as implemented in the nebular package (Shaw and Dufour (sha95 (1995)) within the IRAF software. The minimum extinction of 0.24 ($`\pm `$0.04) is measured in area 5, while the maximum, statistically significant measurement of 0.8 ($`\pm `$0.2) is found in area 1a. The absorption at area 1 and 2 is higher than in other parts of the remnant. This may be related to the presence of excess emission seen in the IRAS 60 $`\mu `$m maps (Wheelock et al. whe94 (1994)) to the east of 19h40m provided that the dust emission lies to the foreground of G 65.3+5.7. The observed variations of the logarithmic extinction may indicate the presence of intrinsic absorption although these variations are at the 3–4 $`\sigma `$ level. Note that variations in the column density, determined from X–ray measurements, are reported by Lu & Aschenbach (lu02 (2002)). The probable non–uniform distribution of the column density and the logarithmic extinction would suggest a complex, inhomogeneous environment into which the primary blast wave propagates. In order to estimate the column density implied by the optical spectra we use the relation NH = 5.4($`\pm `$0.1) $`\times `$ $`10^{21}`$ E(B-V) cm-2 mag-1 (Predehl and Schmitt et al. pre95 (1995)). Assuming a logarithmic extinction of 0.36($`\pm `$0.03), the average from areas 3, 4, and 5, which is equivalent to an E(B-V) of 0.24 (E(B–V) = 0.664 c ; Kaler kal76 (1976), Aller all84 (1984)), a column density of 13($`\pm 1`$) $`\times `$1020 cm-2 is calculated. This column density exceeds that of 3–5$`\times `$1020 cm-2 measured in the ROSAT pointed data (Lu & Aschenbach lu02 (2002)) by factors in the range of 2.6 - 4.3. It is not clear at the moment if this implies the presence of excess intrinsic absorption affecting mainly the optical wavelengths or is due to a statistical effect since these factors are just compatible with the scattering around the straight line fit of Predehl and Schmitt (pre95 (1995)).

The supernova remnant G 65.3+5.7 was observed in the soft X–ray band by ROSAT through a large number of pointed observations as well as during the All–Sky survey (Lu & Aschenbach lu02 (2002)). In Fig. 4 we show the $`[`$O iii$`]`$ emission together with the contours of the soft X–ray emission detected by ROSAT. The outer $`[`$O iii$`]`$ filaments are rather well correlated with the X–ray emission mainly in the south, and north, north–east areas. The southern $`[`$O iii$`]`$ filaments probably lie close to the primary shock wave, while X–ray emission is detected ahead of the northern $`[`$O iii$`]`$ filaments. The clumpiness of the X–ray emission in the central areas of the remnant does not allow a reliable comparison with the optical emission. The temperature or column density maps do not provide any clear evidence for a correlation or anticorrelation with the optical data. There is also no well defined X–ray emission outside the shell of the remnant around $`\alpha `$ $``$ 19h28m, $`\delta `$ $``$ 33°00′ suggesting that the detected optical emission may not be related to it.

## 7 Conclusions

The supernova remnant G 65.3+5.7 was observed in the lower and medium ionization lines of $`[`$O ii$`]`$3727 Å and $`[`$O iii$`]`$5007 Å. The $`[`$O iii$`]`$ emission line mosaic is more filamentary than the $`[`$O ii$`]`$ image. Several new structures are detected, while some of them are found outside the X–ray shell indicating that those may not be related to G 65.3+5.7. The morphology of the $`[`$O ii$`]`$ and $`[`$O iii$`]`$ images suggest the presence of small scale inhomogeneities which are of the order of $``$0.2 pc. Larger scale inhomogeneities are also present. The long slit spectra clearly demonstrate the very strong $`[`$O iii$`]`$ emission relative to H$`\beta `$ suggesting the presence of incomplete recombination zones. The shock velocities are found in the range of 90–140 km s-1, while the column densities behind the radiative shock cover the range of $``$ 10$`{}_{}{}^{17.0}10^{18.5}`$ cm-2. The extinction variations seen in the optical spectra may indicate the presence of internal absorption.

###### Acknowledgements.

The authors would like to thank J. Raymond and the referee B. Aschenbach for their useful comments and suggestions. Skinakas Observatory is a collaborative project of the University of Crete, the Foundation for Research and Technology-Hellas and the Max-Planck-Institut für Extraterrestrische Physik. This research has made use of data obtained through the High Energy Astrophysics Science Archive Research Center Online Service, provided by the NASA/Goddard Space Flight Center.


text:
# Warped products and Reissner-Nordstrom metric

## Abstract.

We study a multiply warped products manifold associated with the Reissner-Nordstrom metric to investigate the physical properties inside the black hole event horizons. It is shown that, different from the uncharged Schwarzschild metric, the Ricci curvature components inside the Reissner-Nordstrom black hole horizons are nonvanishing, while the Einstein scalar curvature vanishes even in the interior of the charged metric. Introducing a perfect fluid inside the Reissor-Nordstrom black hole, it is also shown that the charge plays effective roles of decreasing the mass-energy density and the pressure of the fluid inside the black hole.

###### Key words and phrases:

multiply warped products, Reissner-Nordstrom metric

I. Introduction

Since the pioneering work in 19761, thermal Hawking effects on a curved manifold2 have been studied as an Unruh effect in a higher flat dimensional spacetime. Following the global embedding Minkowski space approach3-6, several authors recently have shown that this approach could yield a unified derivation of temperature for various curved manifolds in (2+1) dimensions7-13 and in (3+1) dimensions7,14-16. However all these higher dimensional embedding solutions have been constructed outside the event horizons of the metrics.

On the other hand, the concept of a warped products manifold was introduced by Bishop and O’Neill17, where it served to provide a class of complete Riemannian manifolds with everywhere negative curvature18. The connection with general relativity was first made by Beem, Ehrich, and Powell, who pointed out that several of the well-known exact solutions to Einstein field equations are pseudo-Riemannian warped products19. Beem and Ehrich further explored the extent to which certain causal and completeness properties20 of a spacetime maybe determined by the presence of a warped products structure. After first developing the general theory of warped products to spaces, O’Neill then applied the theory to discuss, in turn, the special cases of Robertson-Walker and Schwarzshild spacetime. The role of warped products in the study of exact solutions to Einstein’s equations is now firmly established, and it appears that these structures are generating interest in other areas of geometry. Very recently, Choi has represented the interior Schwarzschild spacetime as a multiply warped products spacetime with warping functions21 to yield the Ricci curvature in terms of $`f_1`$, $`f_2`$ for the multiply warped products of the form $`M=R^1\times _{f_1}R^1\times _{f_2}S^2`$.

In this paper we will further analyze the multiply warped products manifold associated with the Schwarzschild metric by including the charge degrees of freedom so that we can investigate the physical properties inside the event horizons of the Reissner-Nordstrom black hole.

In section II, we will briefly recapitulate the multiply warped products manifold and the corresponding Ricci curvature components. In section III, we will apply the multiply warped products manifold scheme to the Reissner-Nordstrom metric to explicitly obtain the Ricci and Einstein curvatures inside the event horizons of the metric and to study the corresponding charge effects. Introducing a perfect fluid inside the Reissor-Nordstrom black hole, we will discuss the charge effects on the mass-energy density and the pressure of the fluid inside the black hole.

II. Multiply warped products and Ricci curvature

In this section, we briefly review a multiply warped products manifold to investigate Ricci curvature inside the black hole horizons.

A Lorentzian manifold $`(M,g)`$ is a connected smooth manifold of $`d`$-dimension $`(d2)`$ with a countable basis together with a Lorentzian metric $`g`$ of signature $`(,+,+,\mathrm{},+)`$. Let $`(F_i,g_i)`$ be Riemannian manifolds, and let $`(B,g_B)`$ be either a spacetime, or let $`B`$ be $`R^1`$ with $`g_B=dt^2`$. Let $`f_i>0`$, $`i=1,\mathrm{},n`$ be smooth functions on $`B`$. A multiply warped products spacetime with base $`(B,g_B)`$, fibers $`(F_i,g_i)`$ $`i=1,\mathrm{},n`$ and warping functions $`f_i>0`$ is the product manifold $`(B\times F_1\times \mathrm{}\times F_n,g)`$ endowed with the Lorentzian metric:

$`(2.1)`$
$$g=\pi _B^{}g_B+\underset{i=1}{\overset{n}{}}(f_i\pi _B)^2\pi _i^{}g_idt^2+\underset{i=1}{\overset{n}{}}f_i^2g_i$$

where $`\pi _B`$, $`\pi _i`$ $`i=1,\mathrm{},n`$ are the natural projections of $`B\times F_1\times \mathrm{}\times F_n`$ onto $`B`$ and $`F_1`$,…,$`F_n`$, respectively.

Thus, warped product spaces are extended to richer class of spaces involving multiply products. Multiply warped products spaces were studied by Flores and Sánchez22. The conditions of spacelike boundaries in the multiply warped products spacetimes were studied by Harris23. The Kasner metric was studied as a cosmological model by Schücking and Heckmann24. Choi has investigated the curvature of a multiply warped product with $`C^0`$-warping functions21.

From a physical point of view, these spacetimes are interesting, first, because they include classical examples of spacetimes: when $`n=1`$ they are generalized Robertson-Walker spacetimes, standard models of cosmology; when $`n=2`$ the intermediate zone of Reissner-Norsdström spacetime and interior of Schwarzschild spacetime appear as particular cases.

III. Reissner-Nordstrom black hole as a multiply warped product manifold

In this section, to investigate a multiply warped product manifold for the Reissner-Nordstrom interior solution, we start with the four-metric inside the horizon

$`(3.1)`$
$$ds^2=N^2dt^2N^2dr^2+r^2d\mathrm{\Omega }^2$$

where the lapse function for the interior solution is given by

$`(3.2)`$
$$N^2=1+\frac{2m}{r}\frac{Q^2}{r^2},$$

with a mass $`m`$ and a charge $`Q`$, and $`d\mathrm{\Omega }^2=d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2`$. Note that, for the nonextremal case, there exist two event horizons $`r_\pm (Q)`$ satisfying the equations $`0=1+2m/r_\pm Q^2/r_\pm ^2`$ such that

$`(3.3)`$
$$r_\pm =m\pm (m^2Q^2)^{1/2}.$$

Furthermore the lapse function can be rewritten in terms of these outer and inner horizons as follows

$`(3.4)`$
$$N^2=\frac{(r_+r)(rr_{})}{r^2}$$

which, for the interior solution, is well defined in the region $`r_{}<r<r_+`$.

Proposition 3.1 Let $`M`$ be a manifold with the Reissor-Nordstrom metric solution $`ds^2=N^2dt^2N^2dr^2+r^2d\mathrm{\Omega }^2`$ (where $`d\mathrm{\Omega }^2=d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2)`$, then we have $`M`$ as a multiply warped manifolds with warping functions $`f_1,f_2`$

$`(3.5)`$
$$ds^2=d\mu ^2+f_1^2(\mu )d\nu ^2+f_2^2(\mu )d\mathrm{\Omega }^2$$

where

$`(3.6)`$
$$\mu =2m\mathrm{cos}^1\left(\frac{r_+r}{r_+r_{}}\right)(r_+r)^{1/2}(rr_{})^{1/2}=F(r),$$

$`(3.7)`$
$$f_1(\mu )=\left(1+\frac{2m}{F^1(\mu )}\frac{Q^2}{F^2(\mu )}\right)^{1/2},$$

$`(3.8)`$
$$f_2(\mu )=F^1(\mu ).$$

Proof: Define a new coordinate $`\mu `$ by

$`(3.9)`$
$$d\mu ^2=N^2dr^2,$$

which can be integrated to yield

$`(3.10)`$
$$\mu =_r_{}^r\frac{dxx}{(r_+x)^{1/2}(xr_{})^{1/2}}+\mu (r_{}).$$

We establish $`\mu (r_{})=0`$ to yield the analytic solution (3.6). Note that in (3.6) $`dr/d\mu >0`$ implies $`F^1(\mu )`$ is well-defined function. Exploiting the above new coordinate (3.6) and redefinition $`\nu =t`$, we rewrite the metric (3.1) as a multiply warped product $`M`$ as in the form of (2.1) to yield (3.5) with $`f_1`$ and $`f_2`$ in (3.7) and (3.8). ∎

Note that in the vanishing $`Q`$ limit the above solution is reduced to that of the uncharged Schwarzschild case21. Moreover, we have the following boundary conditions for $`F(r)`$

$`(3.11)`$
$$\mathrm{lim}_{rr_+}F(r)=m\pi ,\mathrm{lim}_{rr_{}}F(r)=0.$$

After some algebra, we obtain the following nonvanishing Ricci curvature components

$$\begin{array}{cc}& R_{\mu \mu }=\frac{f_1^{^{\prime \prime }}}{f_1}\frac{2f_2^{^{\prime \prime }}}{f_2},\hfill \\ & R_{\nu \nu }=\frac{2f_1f_1^{^{}}f_2^{^{}}}{f_2}+f_1f_1^{^{\prime \prime }},\hfill \\ & R_{\theta \theta }=\frac{f_1^{^{}}f_2f_2^{^{}}}{f_1}+f_2f_2^{^{\prime \prime }}+f_2^{}_{}{}^{}2+1,\hfill \\ & R_{\varphi \varphi }=\left(\frac{f_1^{^{}}f_2f_2^{^{}}}{f_1}+f_2f_2^{^{\prime \prime }}+f_2^{}_{}{}^{}2+1\right)\mathrm{sin}^2\theta ,(3.12)\hfill \end{array}$$

which hold also in the Schwarzschild metric case21. Note that, as shown in (3.6)-(3.8), the argument $`\mu `$ of $`f_1`$ and $`f_2`$ in (3.12) and $`f_1`$ itself are however different from those of the Schwarzschild black hole, since $`\mu `$ is described in terms of $`r_\pm `$ possessing the charge and $`f_1`$ has the additional charge term.

Proposition 3.2 Let $`M`$ be a multiply warped product manifold $`M=R^1\times _{f_1}R^1\times _{f_2}S^2`$ for the Reissor-Nordstrom metric solution $`ds^2=d\mu ^2+f_1^2(\mu )d\nu ^2+f_2^2(\mu )d\mathrm{\Omega }^2`$ (where $`d\mathrm{\Omega }^2=d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2)`$ with warping functions $`f_1,f_2`$, then we have Ricci curvature components as follows,

$$\begin{array}{cc}& (1)R_{\mu \mu }=\frac{Q^2}{f_2^4},(2)R_{\nu \nu }=\frac{Q^2f_1^2}{f_2^4},\hfill \\ & (3)R_{\theta \theta }=\frac{Q^2}{f_2^2},(4)R_{\varphi \varphi }=\frac{Q^2}{f_2^2}\mathrm{sin}^2\theta .\hfill \end{array}(3.13)$$

Proof: Using the explicit expressions for $`f_1`$ and $`f_2`$ in (3.7) and (3.8), we can obtain identities for $`f_1`$, $`f_1^{^{}}`$ and $`f_1^{^{\prime \prime }}`$ in terms of $`f_1`$, $`f_2`$ and their derivatives

$$\begin{array}{cc}& f_1=f_2^{^{}},\hfill \\ & f_1^{^{}}=\frac{m}{f_2^2}+\frac{Q^2}{f_2^3},\hfill \\ & f_1^{^{\prime \prime }}=\frac{2f_1f_1^{^{}}}{f_2}\frac{Q^2f_1}{f_2^4}.(3.14)\hfill \end{array}$$

Substituting (3.14) into (3.12), we evaluate $`R_{\mu \mu }`$ as follows,

$$R_{\mu \mu }=2\left(\frac{f_1^{^{}}}{f_2}+\frac{Q^2}{2f_2^4}\right)\frac{2f_2^{^{\prime \prime }}}{f_2}=\frac{Q^2}{f_2^4}.$$

Similarly, we establish the other components $`R_{\nu \nu }`$, $`R_{\theta \theta }`$ and $`R_{\varphi \varphi }`$. ∎

Here one notes that, differently from the Schwarzschild case where we have the flat Ricci curvature components as shown with $`Q=0`$, we have the nonvanishing Ricci tensor components for $`r_{}<r<r_+`$. Moreover, it is amusing to see that, exploiting the metric (3.5), the Einstein scalar curvature is given as follows

$`(3.15)`$
$$R=0,$$

even in the interior of the charged Reissner-Nordstrom black hole horizons. However, due to the charge of the black hole, the vacuum solution in the uncharged Schwarzschild case does not hold any more as shown in (3.13).

Next, in order to investigate the charge effects inside the Reissor-Nordstrom black hole, we assume a perfect fluid which is a continuous distribution of matter with stress-energy tensor $`T_{ab}`$ of the form,

$`(3.16)`$
$$T_{ab}=\rho u_au_b+P(g_{ab}+u_au_b),$$

where $`u^a`$ is a unit timelike 4-velocity vector field of the fluid, $`\rho `$ is the mass-energy density and $`P`$ is the pressure of the fluid as measured in its rest frame.

Proposition 3.3 In the interior ($`r_{}<r<r_+`$) of charged Reissor-Nordstrom metric solution, we obtain the following equations.

$$\begin{array}{cc}& (1)R_{\mu \mu }8\pi T_{\mu \mu }=\frac{Q^2}{f_2^4}8\pi Pf_1^2=0,\hfill \\ & (2)R_{\nu \nu }8\pi T_{\nu \nu }=\frac{Q^2f_1^2}{f_2^4}+8\pi \rho =0,\hfill \\ & (3)R_{\theta \theta }8\pi T_{\theta \theta }=\frac{Q^2}{f_2^2}8\pi Pf_2^2=0,\hfill \\ & (4)R_{\varphi \varphi }8\pi T_{\varphi \varphi }=\frac{Q^2}{f_2^2}\mathrm{sin}^2\theta 8\pi Pf_2^2\mathrm{sin}^2\theta =0.\hfill \end{array}(3.17)$$

Proof: Substituting into the Einstein equations of motion

$`(3.18)`$
$$R_{ab}\frac{1}{2}Rg_{ab}=8\pi T_{ab},$$

the Ricci components given by (3.13), the Einstein scalar curvature (3.15) and the stress-energy tensor (3.16) given in terms of the warped products, we obtain the above results.∎

Here one notes that in the case of the vanishing $`Q`$ limit which is the Schwarzschild case, one has contributions from the mass-energy density and the pressure of the fluid in the interior of the black hole. On the other hand, in the charged Reissor-Nordstrom case, the existence of nonvanishing charge $`Q`$ plays effective roles of decreasing the mass-energy density and the pressure of the fluid inside the black hole, since the terms associated with the charge are positive definite in the Einstein equations of motion. From a physical viewpoint, the above result is quite consistent with the phenomenology that the charge generates a repulsive force inside the charged black hole to dilute the density in the perfect fluid itself.

IV. Conclusions

We have studied a multiply warped product manifold associated with the Reissor-Nordstrom metric to evaluate the Ricci curvature components inside the charged black hole horizons. Differently from the uncharged Schwarzschild metric where both the Ricci and Einstein curvatures vanish inside the horizon, the Ricci curvature components inside the Reissner-Nordstrom black hole horizons are nonvanishing, even though the Einstein scalar curvature vanishes in the interior of the charged metric. Introducing a perfect fluid inside the Reissor-Nordstrom black hole, we have also explicitly evaluated the Einstein equations of motion inside the horizons to investigate the charge effects that the charge decreases the mass-energy density and the pressure of the fluid inside the black hole.

ACKNOWLEDGMENTS

One of us (JC) would like to thank Prof. John K. Beem for helpful discussions and kind concerns. STH, JC and YJP would like to acknowledge financial support in part from Korea Science and Engineering Foundation (R01-2000-00015), (R01-2001-00003), and Ministry of Education, BK21 Project (D-1099), respectively.

1 W.G. Unruh, Phys. Rev. D14, 870 (1976); P.C.W. Davies, J. Phys. A8,

609 (1975). 2 S.W. Hawking, Comm. Math. Phys. 42, 199 (1975). 3 E. Kasner, Am. J. Math. 43, 130 (1921); C. Fronsdal, Phys. Rev. 116, 778

(1959). 4 J. Rosen, Rev. Mod. Phys. 37, 204 (1965). 5 H. F. Goenner, General Relativity and Gravitation, edited by A. Held (Plenum,

New York, 1980) 441. 6 H. Narnhofer, I. Peter and W. Thirring, Int. J. Mod. Phys. B10, 1507 (1996). 7 S. Deser and O. Levin, Class. Quant. Grav. 14, L163 (1997); Class. Quant.

Grav. 15, L85 (1998); Phys. Rev. D59, 0640004 (1999). 8 M. Banados, C. Teitelboim and J. Zanelli, Phys. Rev. Let. 69, 1849 (1992). 9 M. Banados, M. Henneaux, C. Teitelboim and J. Zanelli, Phys. Rev. D48,

1506 (1993). 10 S. Carlip, Class. Quant. Grav. 12, 2853 (1995). 11 D. Cangemi, M. Leblanc and R.B. Mann, Phys. Rev. D48, 3606 (1993). 12 S.T. Hong, Y.W. Kim and Y.J. Park, Phys. Rev. D62, 024024 (2000). 13 S.T. Hong, W.T. Kim, Y.W. Kim and Y.J. Park, Phys. Rev. D62, 064021

(2000). 14 L. Andrianopoli, M. Derix, G.W. Gibbons, C. Herdeiro, A. Santambrogio and

A. V. Proeyen, Class. Quant. Grav. 17, 1875 (2000). 15 S. W. Hawking and H. S. Reall, Phys. Rev. D61, 024014 (1999). 16 Y.W. Kim, Y.J. Park and K.S. Soh, Phys. Rev. D62, 104020 (2000). 17 R.L. Bishop and B. O’Neill, Am. Math. Soc. 145, 1 (1969). 18 R.L. Bishop and B. O’Neill, Trans. A.M.S. 145, 1 (1969) 19 J.K. Beem, P. E. Ehrlich and K. Easley, Global Lorentzian Geometry (Marcel

Dekker Pure and Applied Mathematics, New York, 1996); T.G. Powell, Ph.D

Thesis, University of Missouri-Columbia (1982); K. Easley, Ph.D Thesis

University of Missouri-Columbia (1980). 20 J.K. Beem and P. E. Ehrlich, Math. Proc. Camb. Phil. Soc. 85, 161 (1979). 21 J. Choi, J. Math. Phys. 41, 8163 (2000). 22 J.L. Flores and M. Sánchez, math.DG/9909075; math.DG/0106174. 23 S.G. Harris, Class. Quant. Grav. 17, 551 (2000). 24 C.W. Misner, K. S. Thorne and J. A. Wheeler, Gravitation (W. H. Freeman

and Company, New York, 1970).