Ranking and validation of spallation models for isotopic production cross sections of heavy residua

The production cross sections of isotopically identified residual nuclei of spallation reactions induced by 136Xe projectiles at 500AMeV on hydrogen target were analyzed in a two-step model. The first stage of the reaction was described by the INCL4.6 model of an intranuclear cascade of nucleon-nucleon and pion-nucleon collisions whereas the second stage was analyzed by means of four different models; ABLA07, GEM2, GEMINI++ and SMM. The quality of the data description was judged quantitatively using two statistical deviation factors; the H-factor and the M-factor. It was found that the present analysis leads to a different ranking of models as compared to that obtained from the qualitative inspection of the data reproduction. The disagreement was caused by sensitivity of the deviation factors to large statistical errors present in some of the data. A new deviation factor, the A factor, was proposed, that is not sensitive to the statistical errors of the cross sections. The quantitative ranking of models performed using the A-factor agreed well with the qualitative analysis of the data. It was concluded that using the deviation factors weighted by statistical errors may lead to erroneous conclusions in the case when the data cover a large range of values. The quality of data reproduction by the theoretical models is discussed. Some systematic deviations of the theoretical predictions from the experimental results are observed.


Introduction
It is generally accepted that for a validation of spallation models one has to use objective criteria besides a qualitative inspection. Such a task seems to be straightforward, however some problems may arise because of specific properties of the cross sections and the quantitative methods used in the analysis. The validation of models is usually performed by the application of deviation factors which measure the agreement between data and model cross sections. Then one must take into consideration the fact that data are always affected by statistical errors. Therefore they have to be treated as random variables. The same is true for deviation factors which are functions of the data. Moreover, most of the modern models describing spallation reactions use Monte Carlo methods for the simulation of interaction processes which causes that the model cross sections are also random variables. This means that the theoretical cross sections are also affected by their statistia e-mail: ufkamys@cyf-kr.edu.pl (corresponding author) cal errors and obviously the same is true for the deviation factors.
The perfect agreement of the experimental and theoretical cross sections which are both random variables cannot be represented by equality of measured σ exp i and calculated σ cal i cross sections but rather by equality of their expectation values which should be fulfilled for the set of all N studied cross sections, i = 1, 2, . . . , N. In practice one introduces various deviation factors (cf. [1]) which contain information on the average distance between the set of N data σ exp i , i = 1, 2, . . . , N and the corresponding set of N calculated cross sections σ cal i . For example the popular H-factor defined as is a root mean square distance between two points in the N -dimensional space evaluated in units of the experimental errors Δ exp i . Assuming that the experimental cross sections are normal random variables with standard deviation equal to Δ exp i and neglecting statistical spread of the calculated values it is possible to find exact probability density function of the H-factor in the case of perfect agreement (see eq. (1)) of measured and model-predicted cross sections [2]. For large N (asymptotically) this function can be replaced by a Gauss function with expectation value E(H) and variance var(H) equal to It is important to note that the individual components of the sum present in the definition of the deviation factor should be invariant in respect to exchange of the experimental σ exp i and theoretical σ cal i cross sections. This guarantees that no cancelation appears of the components in which the data are simultaneously partly overestimated and partly underestimated by the model. The H-factor fulfills this condition. This is also true for the M -factor (see eq. (5)) which we have used successfully in the previous validation and ranking of the spallation models [3] applied to a large set of double differential cross sections d 2 σ/dΩ dE for intermediate mass fragments from p+Ag collisions measured by Green et al. [4] for a proton beam energy of 480 MeV, For the perfect agreement of the experimental and theoretical cross sections (see eq. (1)) and under the same assumptions concerning the data as for the H-factor it is possible to show [2] that the expectation value and the variance of the M -factor are equal to For large N the M -factor behaves like the normal random variable with expectation value and variance given above.
In the previous analysis [3] the H-and M -factors were modified in respect to definitions (see eqs. (2) and (5)) to take into account that the model calculations were performed by means of theoretical models which use Monte Carlo methods to evaluate the cross sections. Therefore the experimental error Δ exp i was replaced in formulae (2) and (5) by where Δ cal i is the statistical error of the theoretical cross section calculated by Monte Carlo method.
In the present work an extension of the investigation of intermediate mass fragments emission from p+Ag collisions at E p = 480 MeV is discussed for the total isotopic production cross sections of residual nuclei measured by Giot et al. [5] in reactions induced by 136 Xe projectiles on hydrogen targets at 500 AMeV. Since the silver nucleus has very close atomic and mass numbers (Z = 47, A = 107) to the xenon nucleus (Z = 54, A = 136) and the collision energy per nucleon was almost the same in both experiments -480 AMeV for p+Ag [4] and 500 AMeV for Xe+p [5]-the results of the present analysis should be compatible with that for the silver data [3].
The model calculations of the isotopic cross sections for heavy residual nuclei from 136 Xe + p collisions at energy 500 AMeV are presented in the second section of the paper. The qualitative comparison of the model cross sections with the data is discussed there as well as the validation of the applied models based on the H-and Mdeviation factors (formulae (2)-(8)).
In the next section a new deviation factor; the A-factor is proposed and applied to the analysis of the studied data.
Results are summarized and discussed in the last section of this paper.

Analysis of isotopic cross sections σ(A|Z) for production of heavy residual nuclei
The isotopic cross sections σ(A|Z) for production of heavy (xenon-like) nuclei have been calculated using the INCL4.6 model [6] coupled to four different models describing the second stage of the reaction: ABLA07 [7], GEM2 [8,9], GEMINI++ [10,11] and SMM [12][13][14][15]. A detailed description of the physical assumptions underlying these models may be found in appropriate references given above. Here we only recall that INCL4.6 describes the first stage of the reaction as a sequence of nucleonnucleon and nucleon-pion collisions accompanied by the emission of nucleons, pions and complex nuclei with mass number A smaller than 9 [6]. The emission of ejectiles built of several nucleons is calculated assuming the surface coalescence of the nucleons with the nucleon which escapes from the nucleus after the cascade. Each of the models of the second stage assumes that the intranuclear cascade leaves the residual, excited nucleus in a thermodynamical equilibrium. The GEM2 model describes the de-excitation of this nucleus by sequential evaporation of particles, GEMINI++ treats this process as a sequential emission of particles passing through a barrier similar to the fission barrier whereas the ABLA07 and SMM models allow for simultaneous multifragmentation of the excited nucleus in competition with the evaporation of particles. All the models take also into account a possibility of fission for very heavy nuclei but due to the low fissility of the current system in the present case such a contribution may be neglected. The calculations were performed with default values of the parameters of all models. Due to this it was possible to judge the predictive power of the applied models.
The experimental data of Giot et al. [5] are compared in fig. 1 with results of the above calculations. As it was explained above, the present data may be treated as a supplement to previously studied double differential cross sections for IMF (intermediate mass fragments; isotopes of Li-Mg elements) production [3]. Thus it may be conjectured that the models which were found to reproduce well the IMF data should also work well for heavy remnants of the reaction.
The following qualitative conclusions may be derived from inspection of fig. 1: i) The bell-like shape of the mass dependence of the cross section σ(A|Z) is reproduced by all the models for all observed elements, i.e., from Nb (Z = 41) up to Ba (Z = 56). Also the magnitudes of the theoretical and experimental cross sections agree reasonably well, especially in the neighborhood of the maximum of the mass dependence. It should be pointed out that these values vary significantly -more than two orders of magnitude from Z = 41 to Z = 53.
ii) The cross sections for elements with small atomic number (Z = 41-45: Nb, Mo, Tc, Ru and Rh) are systematically underestimated by the model calculations.
iii) The cross sections for elements with larger atomic number (Z = 46-52) are well reproduced for the smallest and average values of the mass number but are systematically overestimated for the largest mass numbers. This is especially pronounced for elements with Z = 47-51 (Ag, Cd, In, Sn, Sb) and for Z = 56 (Ba).
iv) The difference between predictions of individual models is not large but some systematic effects are observed. For example, the GEM2 and SMM cross sections are systematically smaller than those evaluated in the frame of ABLA07 and GEMINI++ for elements with the smallest atomic numbers: Z = 41 (Nb) to Z = 47 (Ag).
In the second step of the analysis the H-and Mdeviation factors were calculated according to formulae (2)- (8). To validate the models their values were standardized by subtracting the expectation values of the deviation factors and dividing by their standard deviations, i.e., by the square root of the variances. For the perfect agreement of the experimental and calculated cross sections (see eq. (1)) these standardized deviation factors

The A-deviation factor
To avoid the problems discussed above it is necessary to apply a deviation factor being independent of the statistical errors and therefore of the absolute values of the cross sections. We propose to introduce a new deviation factor which we called the A-factor: This factor has the following appealing properties: i) Its value is independent of the absolute value of the cross sections.
ii) It is invariant in respect to exchange of the experimental and calculated cross sections thus no canceling of the components of the sum appears due to the overestimation of one part of the data by the models together with the underestimation of the other part.
iii) Its value belongs to the interval [0, 1] and allows for intuitive reasoning: it is equal to zero when all calculated cross sections are equal to corresponding experimental cross sections and increases with increase of the differences between the data and theoretical values. Furthermore, very small values of the A-factor correspond to a situation where σ cal i ≈ σ exp i therefore the A-factor value may be interpreted as a half of the average relative distance between the experimental and theoretical cross sections. For example A = 0.05 appears when the average relative distance of data and calculated values is equal to 10%. This property of the A-factor allows to apply it to the validation of models in spite of the lack of a detailed knowledge of its probability distribution function for the case of a perfect agreement of the data and model cross sections. Of course, the A-factor may be also used for ranking of models because its value increases monotonically with increase of the relative distance of the above cross sections.
The A-factor values calculated for all data and corresponding combinations of models are presented in fig. 3 (the upper panel) together with the non-standardized values of the H-factor (the middle panel) and the M -factor (the lower panel). As can be seen the Z dependence of the A-factor which has a minimum in neighbourhood of Z = 53-55 and increases almost monotonically with decreasing value of Z is completely different from those of the H-factor and M -factor. On the other hand it agrees with that derived from the qualitative judgement of the data description (cf. fig. 1).
The smallest value of the A-factor is equal to A ≈ 0.15 which corresponds to the average relative distance between the data and calculated cross sections of about 30%. Therefore none of the models can be validated according to the A-factor values. This conclusion agrees with that derived from application of the H-factor and the M -factor. The next step of the analysis consisted in establishing the ranks of the applied theoretical models separately for each element -the reaction product using the A-factor value as the quality criterion. The ranks of all models are collected in table 1 for individual elements and the sum of the ranks over all elements is also depicted. In the case of equally good reproduction of the data by two or more models the tied rank (average rank) was listed. The sum of ranks over all elements was used to determine the average rank of each model.
The smallest values of the A-factor appear at Z = 55 and they are similar for all models with SMM being slightly superior and GEM2 inferior in respect to other models. Similar situation is also present for smaller Z values up to Z = 49. The situation changes for still smaller atomic numbers of the reaction products. Then GEM-INI++ gives systematically the smallest values of the Afactor and SMM the largest whereas those obtained with ABLA07 and GEM2 vary from element to element. The sum of ranks over all elements is equal to: 28, 40, 43.5 and 48.5 for GEMINI++, SMM, ABLA07 and GEM2, respectively. Thus the A-factor values lead to the following ranks of the models, averaged over all elements: 1) GEMINI++, 2) SMM, 3) ABLA07 and 4) GEM2.
Different ranks were obtained using H-deviation factor as well as M -deviation factor. According to first of them the sum of ranks is equal to: 27.5, 40.5, 45.5 and 46.5 for SMM, GEMINI++, ABLA07 and GEM2, respectively. The ranking according to M -deviation factor leads to the sum of ranks equal to: 28, 30.5, 46 and 55.5 for SMM, GEMINI++, GEM2 and ABLA07, respectively. Thus both these deviation factors give an order of the two first average ranks which is opposite to the A-deviation factor: 1) SMM, 2) GEMINI++. The average ranks of ABLA07 and GEM2 based on H-deviation factor are the same as those based on A-factor, i.e. 3) ABLA07, 4) GEM2 but the application of the M -deviation factor leads to opposite ranks of ABLA07 and GEM2: 3) GEM2, 4) ABLA07. It is worth noticing that the sums of ranks of  ABLA07 and GEM2 for the H-deviation factor are so close one to the other (45.5 and 46.5) that both models may be treated as equally good. The same situation appears for the sums of ranks of SMM and GEMINI++ (28 and 30.5) for the M -deviation factor, however it is not the case for ABLA07 and GEM2 (46 and 55.5). Thus H-and M -deviation factors seem to lead to different conclusions. This fact together with qualitatively wrong behavior of H-and M -factors for elements with smallest atomic numbers indicates that they cannot be treated as reliable tools for establishing the ranking of models for an investigated set of data. The conclusions based on the A-factor values seem to be more reliable.

Summary and discussion
In the present work the large set of isotopically resolved production cross sections σ(A, Z) of heavy reaction products from Xe+p collisions at 500 AMeV energy [5] has been compared with model calculations based on the twostep model. First step of the reaction was described by intranuclear cascade model INCL4.6 [6] whereas the second step by four different models of de-excitation of the heavy equilibrated residuum of the first stage of the process: ABLA07 [7], GEM2 [8,9], GEMINI++ [10,11] and SMM [12][13][14][15]. Default parameter values of all the models were applied with the aim to study predictive power of the models. Validation of the models was undertaken comparing qualitatively and quantitatively the data with the theoretical cross sections. The qualitative comparison of experimental and theoretical isotopic distributions σ(A|Z) is presented in fig. 1 for all elements from Nb (Z = 41) to Ba (Z = 56). The main properties of experimental distributions, i.e. their shape and magnitude were reasonably well reproduced by all the models. The quantitative analysis was done using H-and M -deviation factors [3]. It turned out that values of the factors were too large to validate the models.
The ranking of the models was undertaken to select the best of them. However, the inspection of dependence of Hand M -factors on atomic number of the reaction products suggested that they are too sensitive to the statistics of the data covering several orders of magnitude to be treated as a reliable tool for this purpose. New deviation factor, the A-factor (eq. (9)), independent of the statistics of the data, has been proposed. Application of this factor in the analysis lead to the same conclusion as concerns validation of the models as the H-and M -factors, i.e. values of the A-factor were too large to accept the model predictions as accurate enough. On the other hand the ranking of the models based on A-factor values (GEMINI++, SMM, ABLA07 and GEM2) was different from that due to Hand M -factors (SMM, GEMINI++, ABLA07, GEM2) and (SMM, GEMINI++, GEM2, ABLA07), respectively.
It is interesting to note that the above conclusions derived from the A-factor values (lack of validation and ranking of the models) agree with those obtained from the analysis of more exclusive data, i.e. double differential cross sections d 2 σ/dΩ dE [3] of many intermediate mass fragments from p+Ag collisions at proton beam energy 480 MeV. Since the Ag nucleus has very close atomic and mass numbers as the Xe nucleus and the 500 AMeV Xe beam energy in Xe+p collisions is almost equivalent to proton beam energy of 480 MeV in p+Ag collisions the present data may be treated as a part of the larger data set consisted of differential and total production cross sections.
The theoretical cross sections for residues with Z = 41-45 are systematically smaller in the present analysis of Xe+p collisions than the experimental ones. The result is common for all four models of the decay of excited remnants of the first stage of the process. This may lead to the hypothesis that the excitation energy of these remnants is underestimated by the intranuclear cascade model of the first step of the process, i.e. light charged particles and intermediate mass fragments emitted from the cascade stage are too energetic. Such a hypothesis is in agreement with observation made for p+Ag collisions [3] that the spectra of intermediate mass fragments predicted by intranuclear cascade including coalescence contain high energy tails larger than observed in the experimental data. The analysis of coincidence data as, e.g., that performed in refs. [16] for Au+Au or [17] for 12 C + Au would be necessary to solve this problem unambiguously.