text:
# Subthreshold antiproton production in proton-carbon reactions

## 1 Introduction

Subthreshold particle production in a nuclear reaction is understood as production below the energy threshold of the considered process in a free nucleon-nucleon (NN) collision. It is thus a nuclear phenomenon which may be explained by rather different assumptions on the properties of nuclear matter and on the interaction dynamics. Hereby it is an open question to what extent subthreshold particle production is governed by properties of the nuclear ground state wavefunction and to what extent by the dynamical properties of nuclear matter, not reflected in the ground state description. The problem is far from a final solution at present and evidently requires a systematic study of high-momentum transfer processes, among them subthreshold particle production.

At COSY a research programme is in progress which is devoted to the systematic investigation of subthreshold particle production. Subthreshold production of $`\mathrm{K}^+`$ mesons  was one of the foundations for the design and building of the ANKE spectrometer , which allows for the measurement of $`\mathrm{K}^+`$ meson production cross sections under the condition of an extremely high background of other particles. Recently, first results  on $`\mathrm{K}^+`$ meson production have been published. Moreover, it was proposed  to extend the research programme to the subthreshold production of $`\mathrm{K}^{}`$ mesons. This made a considerable extension of the ANKE spectrometer necessary, because additional equipment for the registration of negative particles had to be installed.

Following this line, in this paper the question is examined if the investigation of subthreshold $`\overline{\mathrm{p}}`$ production is manageable at COSY energies. The highest energy of $`2.8\mathrm{GeV}`$ reached at COSY up to now is far below the threshold energy of $`5.63\mathrm{GeV}`$ for $`\overline{\mathrm{p}}`$ production in NN interactions. This makes the measurement a real challenge due to the low cross section and due to the high background of other particles. From the physical point of view, however, deep subthreshold production far below the threshold becomes especially interesting, because the collective aspects of the phenomenon can be expected to become more and more dominating. Indeed, the created mass in case of $`\overline{\mathrm{p}}`$ production ($`1.876\mathrm{GeV}`$) is essentially higher compared to $`\mathrm{K}^+`$ production ($`0.67\mathrm{GeV}`$), and the energy is $`2.83\mathrm{GeV}`$ below the NN threshold if a $`\overline{\mathrm{p}}`$ production experiment would be carried out at $`2.8\mathrm{GeV}`$. Such a difference is significantly larger than the value of $`0.58\mathrm{GeV}`$ achieved so far at ANKE in the $`\mathrm{K}^+`$ production experiments at an energy of $`1\mathrm{GeV}`$ .

First measurements of subthreshold $`\overline{\mathrm{p}}`$ production in proton-nucleus collisions at higher energies had been carried out  a few decades ago. Then investigations in nucleus-nucleus collisions  and more recent studies of proton-nucleus reactions  followed. Most of the descriptions proposed so far are based on transport calculations , thermodynamical considerations  or multi-particle interactions . All these approaches consider the $`\overline{\mathrm{p}}`$ spectra without any relation to measurements of the other reaction channels. At KEK , however, the spectra of $`\pi ^\pm `$ and $`\mathrm{K}^\pm `$ mesons were measured together with those of the antiprotons. The present approach is distinguished by considering simultaneously all these reaction channels. This allows for a more comprehensive determination of the model parameters and makes a prediction for COSY more profound.

It should be stressed that the investigation of subthreshold production in connection with the other reaction channels is also of great importance for the understanding of this process. Fragmentation of the target residue is such a distinguished channel having a large cross section. In a participant-spectator picture there seems to be rather a weak connection between the interaction of the participants leading to particle production and the excitation of the residual nucleus during this process. At subthreshold energies, however, the competition between these two energy consuming processes may heavily influence the cross sections observed. In  this question was considered for proton-induced subthreshold production of $`\mathrm{K}^+`$ mesons. There, the interplay of subthreshold production and the fragmentation of the residual nucleus was investigated in the Rossendorf collision (ROC) model, which allows the treatment of hadronic and nuclear reactions in a unified way. Here, we extend these considerations to the case of subthreshold $`\overline{\mathrm{p}}`$ production. Both the fragmentation of the residual nucleus and the interaction of the nucleons participating in the scattering process are treated on the basis of analogous assumptions. Even more important, the phase-space of the complete final state including the reaction products of the projectile-participant interaction as well as the fragments of the decay of the spectator system is exactly calculated. To the best of our knowledge this feature seems to be unique to the ROC model.

In the ROC model the nuclear residue becomes excited during the reaction due to the distortion of the nuclear structure by the separation of the participants from the spectators and due to the passage of the reaction products through the spectator system. In this global way final-state interactions are taken into account without making special assumptions concerning re-absorption, re-scattering, self-energies, potentials etc. for the various particle types. The discussion of these properties is rather controversial. If one considers e.g. re-absorption of $`\overline{\mathrm{p}}`$ then the treatment reaches from no explicit consideration  over a global factor 0.1 for p+Cu  to taking the mean free path into account . For Si+Si the survival rates of antiprotons are according to  in the region of several %, while in  a reduction of about 50% for p+Au and Si+Au from an analysis of $`\overline{\mathrm{p}}`$ production at higher energies  has been reported. The authors  conclude from an analysis of the target-mass dependence of their data that re-absorption of $`\overline{\mathrm{p}}`$ is much smaller than expected from the corresponding $`\overline{\mathrm{p}}`$N cross section. Thus, it is far from being clear to what degree antiprotons are absorbed in the nuclear medium.

The plan of the paper is as follows. In section 2 the main ingredients of the ROC model are explained, which is used for the calculations to be presented. Section 3 contains a comparison of theoretical and experimental results for particle production and fragmentation with special emphasis on subthreshold $`\overline{\mathrm{p}}`$ production. A prediction of the antiproton cross section at COSY energies and a comparison with the cross sections for $`\pi ^{}`$ and $`\mathrm{K}^{}`$ production is made in section 4. A summary is given in section 5.

## 2 The model

The ROC model is implemented as a Monte Carlo generator which samples complete events for hadronic as well as nuclear reactions. It makes no detailed assumptions on the intra-nuclear development of the interaction process, but calculates instead the statistical weights of the possible final states. The dynamics of the reaction is taken into account in form of empirical functions which modify the population of the final states. This is in contrary to transport models, where the interaction of the projectile with target nucleons and the subsequent interactions of particles originating from primary collisions with further target nucleons are modelled. The ROC approach does not need any parameterisation of elementary cross sections. Its applicability is not restricted at higher energies by the growing number of unknown elementary cross sections. All dynamic information is gathered in a few parameters which are either constant or change smoothly with energy and/or target mass. The energy necessary for subthreshold production stems from interactions of the projectile with few-nucleon groups, also called clusters in the following, and from the Fermi motion of these clusters.

The model was successfully tested for pp interactions up to ISR energies in , while nuclear reactions were considered in the papers . In the following the basic ideas of the ROC model relevant for this paper are summarised.

The cross section of the interaction of an incoming hadron with a target nucleus $`(A,Z)`$ consisting of $`AZ`$ neutrons and $`Z`$ protons is calculated as an incoherent sum over contributions from a varying number $`a`$ of nucleons (thereof $`z`$ protons) participating in the interaction

$$\mathrm{d}\sigma (s)=\underset{a=1}{\overset{A}{}}\underset{z=max(0,aZ)}{\overset{min(a,Z)}{}}\sigma _{az}\frac{\mathrm{d}W(s;𝜶_{az})}{_{𝜶_{az}}dW(s;𝜶_{az})}.$$
(1)
Here, $`s=P^2`$ denotes the square of the centre-of-mass energy of the projectile-target system with $`P=(E,\stackrel{}{P})`$ being the total four-momentum, $`\sigma _{az}`$ stands for the partial cross section of the interaction of the projectile with a cluster $`(a,z)`$, and the quantities $`\mathrm{d}W(s;𝜶_{az})`$ describe the relative probabilities of the various final channels $`𝜶_{az}`$.

The partial cross sections $`\sigma _{az}`$ account for sequential collisions between projectile and $`a`$ nucleons of the target. Such a sequence is treated in the ROC model as a collision between the incident particle and a cluster consisting of $`a`$ nucleons capable of sharing their energy in close analogy with the virtual clusters of the cooperative model . The $`\sigma _{az}`$ are calculated using a modified version of the Monte Carlo code  which is based on a probabilistic interpretation of the Glauber theory . We use the profile function

$$\mathrm{\Gamma }_A(b)=\left[1\underset{i=1}{\overset{A}{}}(1p_i)\right]\underset{i=1}{\overset{A}{}}\rho _A(\stackrel{}{r}_i)\mathrm{d}^3r_i$$

of the considered nucleus, which depends on the nucleon density $`\rho _A(\stackrel{}{r}_i)`$ and the probability

$$p_i=\mathrm{exp}(d_i^2\pi /\sigma _{NN})$$

for an interaction of the projectile and the $`i`$th target nucleon with $`d_i`$ being the distance between the interacting particles. The nucleon density

$$\rho _A(\stackrel{}{r})(1+\eta [1.5(f^2e^2)/f^2+e^2r^2/f^4])\mathrm{exp}(r^2/f^2)$$
(2)
of light nuclei $`A<20`$ can be derived from a standard shell model wavefunction with $`\eta =(A4)/6`$ , $`f^2=e^2(11/A)`$ and $`f=1.55\mathrm{fm}`$. Then the NN cross section $`\sigma _{NN}`$ is adapted such that the integral of the profile function over the impact parameter $`b`$ reproduces the total inelastic $`pA`$ cross sections

$$\sigma _{pA}^{\mathrm{in}}=\mathrm{d}^2b\mathrm{\Gamma }_A(b),$$

which are approximately constant in the energy region under consideration. The same calculation yields also the partial cross sections $`\sigma _{az}`$ we are interested in (for further details see ).

The relative probabilities of the various final channels $`𝜶_{az}`$ in (1)

$$\mathrm{d}W(s;𝜶_{az})\mathrm{d}L_n(s;𝜶_{az})\rho _A(\stackrel{}{P}_R)T^2(𝜶_{az})$$
(3)
are given by the Lorentz-invariant phase-space factor $`\mathrm{d}L_n(s;𝜶_{az})`$ multiplied by the square of the empirical reaction matrix element $`T^2(𝜶_{az})`$ responsible for the collision dynamics. The Fermi motion is implemented via the momentum distribution of the residual nucleus $`\rho _A(\stackrel{}{P}_R)`$, which is made a function of the number of participants $`a`$. It is taken as a Gaussian having a width of

$$\sigma _a=\sqrt{a(Aa)/5/(A1)}p_F$$
(4)
in accordance with the independent particle model  with $`p_F`$ being the Fermi-limit of the nucleus considered. No special high-momentum component  is used in this paper.

The Lorentz-invariant phase-space is defined as the integral over the four-momenta of the $`n`$ primarily produced final particles with energy-momentum conservation taken into account

$$\mathrm{d}L_n(s;𝜶_{az})=\underset{i=1}{\overset{n}{}}\mathrm{d}^4p_i\delta (p_i^2m_i^2)\delta ^4(P\underset{i=1}{\overset{n}{}}p_i).$$
(5)
Here, the four-momentum of the $`i`$-th particle is denoted by $`p_i=(e_i,\stackrel{}{p}_i)`$ with $`p_i^2=m_i^2`$.

For numerical calculations the $`\delta `$-function in equation (5) has to be removed by introducing a new set of $`3n4`$ variables to replace the $`3n`$ three-momentum components. It is reasonable to choose a set of variables, which reflects the underlying physical picture of the interaction process. Using recursion  equation (3) can be rewritten in the form (see Appendix)

$$\mathrm{d}W(s;𝜶_{az})\frac{\mathrm{d}^3P_R}{2E_R}\rho _A(\stackrel{}{P}_R)\mathrm{d}W_R(M_R^2)\mathrm{d}W_P(M_P^2).$$
(6)
Here, the four-momenta of the nuclear residue $`P_R=(E_R,\stackrel{}{P}_R)`$ and of the participants

$$P_P=PP_R=(E_P,\stackrel{}{P}_P),$$
(7)
yield the corresponding invariant masses according to $`M_R^2=P_R^2`$ and $`M_P^2=P_P^2`$. The integral over the Fermi motion $`\mathrm{d}^3P_R\rho _A(\stackrel{}{P}_R)/2E_R`$ separates the phase-space of the $`n_R`$ nuclear fragments

$$\mathrm{d}W_R(M_R^2)=\mathrm{d}M_R^2\mathrm{d}L_{n_R}(M_R^2)T_R^2(𝜶_{az}^R)$$
(8)
from the phase-space

$$\mathrm{d}W_P(M_P^2)=\mathrm{d}L_{n_P}(M_P^2)T_P^2(𝜶_{az}^P)$$
(9)
of the $`n_P=nn_R`$ final particles arising from the participant system. In (8) and (9) the matrix element

$$T^2(𝜶_{az})=T_R^2(𝜶_{az}^R)T_P^2(𝜶_{az}^P)$$
(10)
is split into two factors describing residue fragmentation into channel $`𝜶_{az}^R`$ and participant interaction resulting in channel $`𝜶_{az}^P`$, accordingly. There is, however, a strong kinematic link between participants and spectators, since invariant mass $`M_P`$ of the participant system, invariant mass $`M_R`$ of the residue and the relative kinetic energy $`\sqrt{s}M_PM_R`$ of these two particle groups are connected by energy-momentum conservation. For particle production to proceed the invariant mass of the participants must exceed the corresponding threshold value $`M_P^{th}`$ which in turn depends on the number of participants. The heavier the effective target is the more energy is available for particle production. Another way to reach the threshold value goes via the Fermi motion, since the participant mass is a function of the momentum vector of the nuclear residue. The excitation energy of the target residue comes usually into play via the spectral function (see e.g. ) derived from electron scattering data. Unique to the ROC model, however, is the treatment of the spectator system in close analogy with the participant subsystem. The ROC model calculates the complete final state of both the participants and the spectators. Thus, the huge amount of final channels of the spectator fragmentation influences directly the final state of the participant system and vice versa.

The term $`\mathrm{d}W_P(M_P^2)`$ in (6), (9) describes the interaction of the incoming proton with the group of participating target nucleons. Such a projectile-cluster reaction is treated in complete analogy to a hadronic reaction. In a first step intermediate particle groups called fireballs (FBs) are produced, which decay into so-called primary particles. The primary particles define the channels for which the weights (9) are calculated. Among them are resonances, which decay subsequently into stable hadrons. In case of the presence of several nucleons among the decay products of a FB they may form nuclear fragments with equal probability for all possible channels.

The dynamical input of the reaction is implemented by the empirical transition matrix element

$$T_P^2(𝜶_{az}^P)=T_\mathrm{i}^2T_{\mathrm{qs}}^2T_{\mathrm{ex}}^2T_\mathrm{t}^2T_\mathrm{l}^2T_{\mathrm{st}}^2,$$
(11)
which describes the interaction process $`T_\mathrm{i}^2`$ resulting in the production of a varying number $`N`$ of FBs ($`N2`$), the production of hadrons $`T_{\mathrm{qs}}^2`$ via the creation of quark-antiquark ($`q\overline{q}`$) pairs, the invariant-mass distribution of the FBs $`T_{\mathrm{e}x}^2`$, the transverse $`T_\mathrm{t}^2`$ and longitudinal $`T_\mathrm{l}^2`$ momentum distribution of the FBs, and, finally, some factors $`T_{\mathrm{st}}^2`$ necessary for the calculation of the statistical weights. The interaction is assumed to proceed via colour exchange leading to the removal of valence quarks or of gluons from the interacting hadrons. Additional up, down and strange quark pairs are created in the ratio

$$u:d:s=1:1:\lambda _s$$
(12)
with $`\lambda _s=0.15`$. They form the varying number of FBs, which subsequently decay into the final hadrons. The transverse momenta of the FBs are restricted by an exponential cut-off (longitudinal phase-space) according to

$$T_\mathrm{t}^2=\underset{I=1}{\overset{N}{}}\mathrm{exp}(\gamma P_{\mathrm{t},\mathrm{I}})$$
(13)
with the mean $`\overline{P}_\mathrm{t}=2/\gamma `$. Two leading FBs, the remnants of the incoming hadron and of the cluster, get in the mean larger longitudinal momenta than the central FBs by weighting the events with

$$T_\mathrm{l}^2=(X_1X_2)^\beta .$$
(14)
Here, the light-cone variables $`X_1=(E_1+P_{\mathrm{z},1})/(e_b+p_{\mathrm{z},b})`$ and $`X_2=(E_2P_{\mathrm{z},2})/(E_PP_{\mathrm{z},P})`$ are used with the four-momentum of the projectile given by $`p_b=(e_b,\stackrel{}{p}_b)`$ and that of the participants by (7). Each FB is characterised by two parameters, a temperature $`\mathrm{\Theta }_{FB}`$ and a volume $`V_{FB}`$. The temperature determines the relative kinetic energy of the particles the FB decays into via

$$T_{\mathrm{e}x}^2(\mathrm{\Theta }_{FB})=\underset{I=1}{\overset{N}{}}(M_I/\mathrm{\Theta }_{FB})K_1(M_I/\mathrm{\Theta }_{FB}),$$
(15)
while the volume defines the interaction region and influences mainly the particle multiplicity via the statistical factor $`T_{\mathrm{st}}^2V_{FB}^{n_P1}`$. In (15) $`K_1`$ stands for the modified Bessel function. Final hadrons are built-up by random recombination of the available quarks during the decay of the FBs. This procedure ensures automatically the conservation of all internal quantum numbers. Resonances decay later on until the final state consisting of stable particles is reached. For a more detailed discussion of the hadronic matrix element the reader is referred to . In this paper we use the same set of parameters as in  for the description of the interaction of the projectile with a single target nucleon. Most of the other terms in (1), however, contain clusters consisting of several target nucleons. The basic parameters of the FBs emerging from such an interaction are fixed by scaling the volume with the number of cluster nucleons $`a`$ and the temperature parameter with $`a^{1/3}`$ according to

$$V_{FB}=V_{FB}^0a\text{and}\mathrm{\Theta }_{FB}=\mathrm{\Theta }_{FB}^{max}a^{1/3}.$$
(16)
It remains to consider the target residue, the structure of which is strongly disturbed by the projectile-participant reaction and subsequent final-state interactions. This leads to the excitation and decay of the spectator system, which is characterised by the two parameters temperature $`\mathrm{\Theta }_R`$ and volume $`V_R`$ in the same way as the FBs emerging from the projectile-participant interaction. The part of the matrix element responsible for the residue fragmentation

$$T_R^2(𝜶_{az}^R)=T_{\mathrm{ex}}^2(\mathrm{\Theta }_R)T_{\mathrm{st}}^2(𝜶_{az}^R)$$

is identical with the corresponding factors (11) applied to the hadronic FBs. In order to restrict the excitation energy transferred to the residue we use the asymptotic approximation of (15) for large mass $`M_R`$ and small temperature $`\mathrm{\Theta }_R`$

$$T_{\mathrm{e}x}^2(\mathrm{\Theta }_R)=\sqrt{M_R/\mathrm{\Theta }_R}\mathrm{exp}(M_R/\mathrm{\Theta }_R).$$
(17)
An impact parameter dependence is assumed for the temperature parameter. This seems to be reasonable, because a peripheral collision with only few participating nucleons should excite the nuclear residue much less than a central collision with many participants. As a first guess we use

$$\mathrm{\Theta }_R=\mathrm{\Theta }_R^{max}[1\mathrm{exp}(a/\overline{a}A^{1/3})]$$
(18)
with $`\overline{a}`$ as parameter, here fixed to $`\overline{a}=0.5`$, which determines how fast the maximal temperature $`\mathrm{\Theta }_R^{max}`$ is reached with increasing number of participants $`a`$.

All factors still necessary for a correct calculation of the relative weights of the various channels are collected in complete analogy with the corresponding factors for the FBs in the term

$$T_{\mathrm{st}}^2(𝜶_{az}^R)=g(𝜶_{az}^R)\left(\frac{V_R}{\left(2\pi \right)^3}\right)^{n_R1}\underset{i=1}{\overset{n_R}{}}(2\sigma _i+1)\mathrm{\hspace{0.17em}2}m_i.$$
(19)
It contains the spin degeneracy factors $`(2\sigma _i+1)`$, the volume $`V_R`$ in which the particles are produced with $`V_R=4\pi (Aa)R_R^3/3`$ determined by the radius parameter $`R_R`$. The quantity $`g(𝜶_{az}^R)`$ is the degeneracy factor for groups of identical particles in the final state of the residue decay and prevents multiple counting of identical states.

Temperature and volume determine as in the case of hadronic reactions the number of final particles and their relative energy. The main difference consists in the value of the temperature parameter, $`\mathrm{\Theta }_{FB}^{max}300\mathrm{MeV}`$ for hadronic and $`\mathrm{\Theta }_R^{max}10\mathrm{MeV}`$ for nuclear systems. New particles can be produced in hadronic reactions, while for nuclear systems, due to the much lower temperature, the nucleons of the initial state are recombined into various fragments without producing new hadrons. The volume parameter defines the distance between the fragments where the strong interaction ceases to work. This volume goes into the calculation of the weight of the final state. Due to the long range of the Coulomb repulsion Coulomb energy is still stored in the system at those distances. This energy is calculated in the Wigner-Seiz-approximation (see ). The random distribution of the Coulomb energy among the charged fragments yields the final momenta of the fragments at large distances.

In the calculation we assume that the fragments emerging from the residue decay are stable against particle decay. The possibility of hot fragments cooling down by subsequent particle emission is also implemented in the model (see ), but it turned out that a readjustment of the volume and/or the temperature parameter yields quite similar results in both cases. So we decided in favour of the easier approach. Moreover, for a light nucleus as <sup>12</sup>C this question is of minor importance.

In concluding this section it should be mentioned, that it is the calculation of the relative probabilities $`\mathrm{d}W(s;𝜶_{az})`$, equation (3), which makes the difference between the cooperative  and the ROC model, although the notion of clusters is similar in both approaches. In the cooperative model the relative probabilities of the various channels are calculated in the pure phase-space limit and the excitation and decay of the target residue is completely neglected. In contrast, the ROC model is far from being a pure statistical approach, because it contains specific information on the properties of nuclei and the interaction of the projectile with few-nucleon groups (clusters) in the nucleus. The model includes strong assumptions on the amplitude of the projectile-cluster interaction which lead to the restriction of excitation energy and momentum transfer to clusters compared to pure phase-space. A further nontrivial assumption concerns the various possibilities for the dissipation of the transferred energy inside the hit clusters either as relative kinetic energy of the decay products or as newly produced particles. In extreme cases the whole energy can be completely accumulated for the production of new particles. The other extreme is the emission of very fast secondaries with energies outside the kinematical limits of a projectile-(single)nucleon interaction without producing new particles. Also the supposition that clusters behave themselves like hadrons with respect to quark statistics is of great importance for the relative weights of the various final reaction channels. All these assumptions form a new approach to the consideration of local excitations and may, therefore, give essential and definite information on the properties of nuclear matter (excitation probability, distribution of deposited energy and transferred momenta, decay modes of few-nucleon fireballs etc).

## 3 Comparison with data

In order to fix the parameters describing the Fermi motion and the decay of the residue we consider the energy spectra of fragments from the reaction of 2.1 GeV protons with carbon measured at $`90^{}`$ (see figure 1). For the volume a radius $`R_R=1.7\mathrm{fm}`$ is taken. This corresponds to about $`1/3`$ of normal nuclear density and is within the region of break-up densities used for example by . In figure 1 the heaviest fragment, $`C`$, arises mainly from the quasi-free interaction of the incoming proton with a single target neutron and may gain energy only from the Fermi motion. Thus, the spectrum is highly sensitive to the width $`\sigma _a`$, equation (4), of the Gaussian for the Fermi motion. On the other hand, the lightest fragments are decay products of the residue and their kinetic energy in the laboratory is a superposition of Fermi motion and relative kinetic energy of the residue fragments. A satisfactory description of the fragment spectra in figure 1 is achieved using $`\mathrm{\Theta }_R^{max}=12\mathrm{MeV}`$, equation (18), and $`p_F=320\mathrm{MeV}/\mathrm{c}`$, equation (4). No attempt was made to improve the description of the spectra by implementing an additional high-momentum component, although especially the high-momentum parts of the spectra of the heavier fragments are underestimated. The spectator-participant picture is an idealisation, and not only the nuclear wavefunction but also secondary interactions between spectators and participants might be responsible for the observed deviations. All the following calculations are carried out with the above fixed values for $`\mathrm{\Theta }_R^{max}`$ and $`p_F`$. It should be stressed in this connection that temperature $`\mathrm{\Theta }_R`$ and volume parameter $`R_R`$ of the nuclear residue are highly correlated, a fact already observed in  for the FB parameters in hadronic reactions.

In figure 2 the momentum spectra of $`\pi ^\pm `$ and $`\mathrm{K}^\pm `$ mesons as well as of $`\overline{\mathrm{p}}`$ are compared to the KEK data . The overall agreement is quite satisfactory in view of the large region of projectile energies and the variety of ejectile species calculated with one fixed parameter set. Particle yields are influenced by the suppression factor $`\lambda _s=0.15`$, equation (12), of strange quarks and by the algorithm for creating the final hadrons from the quarks produced in the first stage of the interaction process. Hadrons are built up in each FB independently according to the rules of quark statistics  by randomly selecting sequences of $`q`$’s and $`\overline{q}`$’s. A $`q\overline{q}`$ gives a meson, while baryons or antibaryons are formed from $`qqq`$ or $`\overline{q}\overline{q}\overline{q}`$. From a given sequence of quarks the different hadrons are formed according to the tables of the particle data group . There is no parameter which directly determines the ratio between meson and baryon production as e.g. in the PYTHIA-LUND model . Only an indirect influence via the temperature and the volume parameter is possible, which change the relative weights of the FBs in dependence on their invariant mass and final particle multiplicity, respectively.

It is therefore quite remarkable that the general trend of the data is well reproduced by using the parameters determined in  from the consideration of hadronic interactions. The moderate increase with energy of the cross sections for the light mesons, the steep increase of that for the antiprotons as well as the shift of the maximum in the $`\overline{\mathrm{p}}`$ spectra towards higher momenta is well reproduced by the calculations.

At the highest energies the ROC calculations obviously tend to underestimate the data. Therefore, we compare in figure 3 as a kind of cross-check a similar data set from  with ROC-model calculations. Here, the spectra are measured at a similar angle but in a much wider momentum region compared to the KEK data. In this case the tendency of the ROC results goes in the opposite direction and overestimates the cross-sections for the heavier ejectiles $`\mathrm{K}^\pm `$ and $`\overline{\mathrm{p}}`$.

The calculated cross sections are sensitive to the assumed temperature of the residual nucleus (see figure 5). It should be stressed that the temperature derived from the fragment spectra yields also a reasonable reproduction of the $`\overline{\mathrm{p}}`$ cross sections. This point is of special importance since it demonstrates that the link between the realms of fragmentation and of hadron production made by the ROC approach is obviously correct.

The authors  interpret their $`\overline{\mathrm{p}}`$ data by using the “first-chance NN collision model” from  where the internal nucleon momenta were extracted from backward proton production  as a superposition of two Gaussian distributions. In this way the momentum spectra and the incident-energy dependence could be successfully reproduced by adapting one normalisation parameter. For the similar case of subthreshold $`\mathrm{K}^+`$ production it has been, however, argued  that the contribution from the high-momentum component should be negligible.

In the ROC model the “first-chance NN collision” corresponds to the term in (1) describing the interaction of the projectile with a single target nucleon. In table 1 the relative contributions of the terms with a given number of participants to the total $`\overline{\mathrm{p}}`$ production cross section are summarised. It can be seen that the interaction with a single nucleon yields the main contribution only at the highest energy. The lower the energy the more contribute the terms with several target nucleons. At the lowest energy, $`2.8\mathrm{GeV}`$, considered here the main contributions to $`\overline{\mathrm{p}}`$ production arise from interactions with three and four target nucleons, while quasi-free collisions are negligible. At the highest energy, $`12\mathrm{GeV}`$, there are still considerable contributions from proton-cluster reactions, although quasi-free collisions became the dominating process. Thus, the interpretation of the data by the ROC model is quite contradictory to the assumptions of a “first-chance NN collision model”. From the viewpoint of the ROC model the key quantity for understanding subthreshold particle production is the number of participating nucleons. A direct experimental determination of this number for subthreshold $`\overline{\mathrm{p}}`$ production is highly desirable as discussed in  for the case of $`\mathrm{K}^{}`$ production.

A certain similarity can be, however, found between the ROC approach and the multi-nucleon mechanism of , where “the incoming proton and the interacting nucleons in the target act as sources of pions that merge to produce a nucleon-antiproton pair”. In this approach the number of participating nucleons obviously plays a similar role as in the ROC model.

Another interesting aspect of subthreshold particle production is the observation that the formation of light nuclei in the final channel alongside the produced particle(s) has to be taken into account for achieving a good reproduction of the measured cross sections. The formation of light nuclei leads to the freezing of relative kinetic energy of the nucleons, which is then available for particle production. We are speaking here about light fragments emerging from clusters, not about the decay of the residual nucleus. This effect has been discussed in  for pion, in  for photon, in  for $`\mathrm{K}^+`$ production in nucleus-nucleus reactions on the basis of the cooperative model, and in , within the ROC model, for subthreshold $`\mathrm{K}^{}`$ production in proton-nucleus interactions. In order to demonstrate the importance of this effect for subthreshold $`\overline{\mathrm{p}}`$ production we have in table 1 in addition to the “normal” ROC results $`\sigma _{\overline{\mathrm{p}}}`$ also the cross sections $`\sigma _{\overline{\mathrm{p}}}^{\mathrm{nf}}`$ with no freezing listed. At the lowest energy considered the difference amounts to about one order of magnitude. This is relatively small compared to the three orders of magnitude reported by  for the reaction $`{}_{}{}^{12}\mathrm{C}(85\mathrm{M}\mathrm{e}\mathrm{V}/A)+{}_{}{}^{12}\mathrm{C}\pi ^++X`$. The energy of 2.8 GeV is still far above the total production threshold of about 2.2 GeV for the coherent production at the whole target nucleus. Contrary to the energy considered in  also the number of competing channels is still rather high. This diminishes the influence of the neglected channels. The slight increase of $`\sigma _{\overline{\mathrm{p}}}^{\mathrm{nf}}`$ compared to $`\sigma _{\overline{\mathrm{p}}}`$ at the higher energies is due to the smoother energy dependence of the relative probability $`\mathrm{d}W(s;𝜶_{az})`$ of the considered channel and the smaller denominator in equation (1) as a result of neglecting a large number of final channels.

In  calculations have been carried out in the framework of transport theory, where questions concerning baryon self-energies, elementary $`\overline{\mathrm{p}}`$ production amplitudes, $`\overline{\mathrm{p}}`$ potentials and $`\overline{\mathrm{p}}`$ re-absorption are of importance for the reproduction of the experimental data. In the ROC model there is no equivalent of these intra-nuclear properties implemented so far. Nuclear properties enter the ROC calculations only via the scaling (16) of the temperature and the volume of the FBs and via the excitation of the target residue caused by secondary participant-spectator interactions. Nevertheless, the spectra of all particle species are comparably well reproduced. This result should be further explored by considering heavier target nuclei, where the problem of secondary interactions becomes more important. On the other hand, the experimental investigation of $`\overline{\mathrm{p}}`$ production on the lightest nuclei as $`{}_{}{}^{3}He`$ and $`{}_{}{}^{4}He`$ is of special interest, since the model predicts in the energy range between 3 and 5 GeV just the 3- and 4-nucleon groups as the main source of $`\overline{\mathrm{p}}`$ production. Secondary interactions of the produced particles and the influence of the excitation of the target residue are practically absent in this case. This makes the interpretation of such data more clear than in the case of heavier target nuclei.

## 4 Antiproton cross section at COSY energies

Figure 4 shows a ROC-model estimate of the antiproton spectrum for a bombarding energy of 2.8 GeV averaged over a polar angle $`\mathrm{\Theta }<10^{}`$. In order to emphasise the experimental difficulties of a measurement at such a deep subthreshold energy, the spectra of the other negative stable particles are depicted, too. In the region of the maximum in the spectra around 500 MeV/c about ten orders of magnitude more $`\pi ^{}`$ mesons than antiprotons must be expected what makes a possible measurement in addition to the small $`\overline{\mathrm{p}}`$ cross section a difficult task.

In figure 5 the energy dependence of the total production cross section for $`\overline{\mathrm{p}}`$ in $`\mathrm{p}+{}_{}{}^{12}\mathrm{C}`$ interactions is plotted. In the energy region considered the cross section changes by more than seven orders of magnitude. The importance of the temperature of the residual nucleus at low energies is demonstrated by a calculation with the temperature arbitrarily decreased by a factor of two. This increases the cross section by about half an order of magnitude at the lowest energy. With increasing energy the influence of this parameter diminishes and becomes negligible above $`10\mathrm{GeV}`$ projectile energy. The Fermi motion is of similar importance for subthreshold $`\overline{\mathrm{p}}`$ production. A decrease of $`p_F`$, equation (4), by about 30% makes the $`\overline{\mathrm{p}}`$ cross section by nearly one order of magnitude smaller. Again the influence becomes negligible at high energies. These results underline once more how important it is to use a model which allows the verification of these aspects by independent measurements as has been done in this paper. It should be stressed that the selected parameter range does not reflect the uncertainty of the considered parameters which can be rather precise determined from the fragmentation data  (see figure 1). Instead, it is the aim to demonstrate the sensitivity of the results to the inclusion of the nuclear residue into the considerations. All results in the subthreshold region from approaches neglecting the excitation of the nuclear residue should be, therefore, taken with caution.

But also the present prediction contains of course considerable uncertainties. From our experience with existing experimental data, in particular those shown in figures 2 and 3, we know that deviations between data and calculations are usually less than a factor of two. Since any drastic change of the model feasibility is not seen if the energy decreases from 3.5 GeV to 2.8 GeV we may expect the same level of uncertainty at this lowest energy. Certainly, for the present case of a prediction rather far below the region where experimental data are available it cannot be excluded that some processes or factors not included in the ROC approach may increase their significance with the energy reduction. A possible example could be the coherent production at the whole target nucleus . To be on the sure side we, therefore, assume at least the order of magnitude of the predicted cross section to be correct.

At the ANKE spectrometer an effective detection of particles produced at angles $`10^{}`$ in the momentum range $`0.21.0\mathrm{GeV}/\mathrm{c}`$ is possible. Assuming a luminosity of $`10^{33}\mathrm{cm}^2\mathrm{s}^1`$ available with a carbon strip target and a cross section of $`110^9\mathrm{mb}`$ expected at $`2.8\mathrm{GeV}`$ one gets a counting rate at the level of five hundreds events per week of beam-time. This estimate shows the feasibility of studying $`\overline{\mathrm{p}}`$ production even at such a low subthreshold energy as considered here, if the experimental equipment allows for an efficient suppression of the negative mesons.

## 5 Summary

Subthreshold particle production is a collective phenomenon which is far from being completely understood. Data on subthreshold particle production can be reproduced within the ROC model by considering the interaction of the projectile with few-nucleon systems in complete analogy to the interaction with a single nucleon, also with respect to high-momentum transfer processes. It is the simultaneous consideration of the data for all particle types measured, which distinguishes the present approach from the previous attempts to interpret solely the subthreshold $`\overline{\mathrm{p}}`$ spectra from KEK. The reproduction of the data requires no strong distortion of antiprotons in nuclear matter, which might be expected due to the fact that the free $`\overline{\mathrm{p}}N`$ cross section is larger than the corresponding ones for the other particle species considered. This finding needs further confirmation by measuring $`\overline{\mathrm{p}}`$ spectra at lower momenta than the region covered by the KEK data where the $`\overline{\mathrm{p}}`$N cross section increases further.

In this paper the feasibility of measuring subthreshold $`\overline{\mathrm{p}}`$ production at COSY-ANKE is demonstrated. It is argued that new insight into the mechanism of subthreshold production could be gained from an experiment nearly $`3\mathrm{GeV}`$ below the NN threshold, where the ROC-model calculations predict increasing contributions to the cross section from few-nucleon groups. The reliability of the predictions is verified by considering together with $`\overline{\mathrm{p}}`$ production also the results for other particles ($`\pi ^\pm `$, $`\mathrm{K}^\pm `$) which could be well described in a wide energy region.

One of the authors (HM) would like to thank W Enghardt for the promotion of this study.

## Appendix

The $`n`$ particles of the phase-space (5) of the final state of the reaction are divided into two groups, the decay products of the residual nucleus and the particles emerging from the projectile-cluster interaction. In a first step the decay products are separated by introducing their invariant mass $`M_R`$ and four-momentum $`P_R`$ via the following identities:

$$1=dM_R^2\delta (P_R^2M_R^2)$$
(20)
$$1=\mathrm{d}^4P_R\delta ^4(P_R\underset{i=1}{\overset{n_R}{}}p_i)$$
(21)
with the definition

$$P_R=\underset{i=1}{\overset{n_R}{}}p_i.$$
(22)
Inserting (20) $`\mathrm{}`$ (22) into (5) the phase-space factor becomes

$`\mathrm{d}L_n(s)`$ $`=`$ $`\mathrm{d}M_R^2\mathrm{d}^4P_R\delta (P_R^2M_R^2)`$ (23)
$`{\displaystyle \underset{i=1}{\overset{n_R}{}}}\mathrm{d}^4p_i\delta (p_i^2m_i^2)\delta ^4(P_R{\displaystyle \underset{i=1}{\overset{n_R}{}}}p_i)`$
$`{\displaystyle \underset{i=1}{\overset{n_P}{}}}\mathrm{d}^4p_i\delta (p_i^2m_i^2)\delta ^4(P_P{\displaystyle \underset{i=1}{\overset{n_P}{}}}p_i).`$
Taking into account that for on-shell particles the identity

$$\mathrm{d}^4p_i\delta (p_i^2m_i^2)=\frac{\mathrm{d}^3p_i}{2e_i}$$
(24)
is valid and considering the definition (5) of the n-particle phase-space equation (23) can be written as

$$\mathrm{d}L_n(s)=\mathrm{d}M_R^2\frac{\mathrm{d}^3P_R}{2E_R}\mathrm{d}L_{n_R}(M_R^2)\mathrm{d}L_{n_P}(P_P^2),$$
(25)
with $`n_P=nn_R`$ being the number of particles emerging from the participant interaction. The invariant mass of this participant group is given by $`M_P=\sqrt{P_P^2}`$ with the four-momentum $`P_P=PP_R`$. The separation into two groups makes sense, because the distribution of the momentum $`\stackrel{}{P}_R`$ can be derived from the internal momentum distribution $`\rho (\stackrel{}{P}_R^2)`$ of the target nucleus.

Inserting the matrix element (10) and the momentum distribution $`\rho (\stackrel{}{P}_R^2)`$ into equation (25) and taking the definitions (8) and (9) into account we arrive at the expression (6) for the relative probabilities of the final channels $`\mathrm{d}W(s;𝜶_{az})`$.

## References


format: text

text:
# Contents

## 1 Introduction

The recent measurements of both neutrino and antineutrino production of dimuon final states (charm signal) by the CCFR and NuTeV collaborations  provide the first promising direct experimental constraints on the strange and anti-strange quark distributions of the nucleon, $`s(x)`$ and $`\overline{s}(x)`$. In addition to the intrinsic interest in nucleon structure , the strange asymmetry ($`s\overline{s}`$) has important implications on the precision measurement of the Weinberg angle in deep inelastic scattering of neutrinos. We report here the first global QCD analysis that includes the new dimuon data, using the methods developed by the CTEQ collaboration, specifically to explore the strange and anti-strange parton parameter space.aaaA preliminary version of this study was reported at the Lepton Photon 2003 International Symposium (LP2003), Fermilab, August 2003. Cf. P. Gambino  and R. Thorne , published in the Proceedings of LP2003.

In previous global analyses, information on $`s`$ and $`\overline{s}`$ has resided only in inclusive cross sections for neutral and charged current DIS. The reliability of the extraction of the quite small $`s`$ and $`\overline{s}`$ componentsbbbThe strangeness content of the nucleon, as measured by the momentum fraction carried by $`s`$ or $`\overline{s}`$, is of order 3% at $`Q=1.5`$ GeV. (from differences of large cross sections measured in different experiments) was always in considerable doubt. For this reason, most global fits adopted the assumption $`s(x)=\overline{s}(x)=\kappa (\overline{u}+\overline{d})/2`$; this approximation was inferred from the earlier combined neutrino and antineutrino dimuon experiments  which extract a consistent value of

$$\kappa \frac{𝑑xx\left[s(x,Q^2)+\overline{s}(x,Q^2)\right]}{𝑑xx\left[\overline{u}(x,Q^2)+\overline{d}(x,Q^2)\right]}0.4$$
(1)
at some low value of $`Q`$. The recent high-statistics dimuon measurements of  provide greater accuracy, and the high purity of separate neutrino and anti-neutrino events offers the first opportunity to study the difference $`s(x)\overline{s}(x)`$. Neutrino induced dimuon production, $`(\nu /\overline{\nu })N\mu ^+\mu ^{}X`$, proceeds primarily through the subprocesses $`W^+sc`$ and $`W^{}\overline{s}\overline{c}`$ respectively, and hence provides independent information on $`s`$ and $`\overline{s}`$.

We present the first global QCD analysis that includes this new dimuon data. The new results demonstrate, first of all, that the strangeness asymmetry, as measured by the momentum integral

$$[S^{}]_0^1x[s(x)\overline{s}(x)]𝑑x,$$
(2)
is indeed more sensitive to the dimuon data than to the other DIS data. We then use the recently developed Lagrange multiplier method of global analysis to explore the range of uncertainty of $`[S^{}]`$. In this first report, we concentrate on $`[S^{}]`$, the integrated strangeness asymmetry, which represents a new parton degree of freedom in the nucleon heretofore largely unexplored, and which has immediate impact on precision electroweak physics because of the NuTeV anomaly. A full exploration of the density functions $`s(x)`$ and $`\overline{s}(x)`$ will be presented later.

We begin by describing the general features of the strangeness sector of the nucleon structure in the QCD framework, and our general parametrization of that sector. A brief review and discussion of the pQCD calculations that are relevant to the interpretation of the dimuon data is then followed by the main results of the global analysis, with emphasis on concrete representative global fits relevant for probing the strangeness asymmetry, and on Lagrange multiplier results for the integrated momentum asymmetry (2). The paper concludes with a summary of the extensive studies performed beyond the examples given, comparisons to previous work on strangeness asymmetry, and conclusions.

## 2 General properties of $`s(x)\overline{s}(x)`$ and its first two moments

Before we discuss the concrete fitting procedure and the results in the following Sections, we feel it is instructive to formulate the qualitative expectations that are based on general QCD requirements (before any data).

### 2.1 Requirements

For each $`Q,`$ let us define the strangeness *number densities* $`s^\pm (x)`$ and their integrals $`[s^^\pm ]`$ by

$$[s^\pm ]_0^1s^\pm (x)𝑑x_0^1[s(x)\pm \overline{s}(x)]𝑑x,$$
(3)
and the *momentum densities* $`S^\pm (x)`$ and integrals $`[S^\pm ]`$ by

$$[S^\pm ]_0^1S^\pm (x)𝑑x_0^1x[s(x)\pm \overline{s}(x)]𝑑x.$$
(4)
In the QCD parton model, certain features of these quantities are necessary:

* The parton distributions $`s(x,Q^2)`$ and $`\overline{s}(x,Q^2)`$cccAs we had done already in Eqs. (3), (4) above, we will suppress the obvious $`Q^2`$ dependence from now on. \[or equivalently $`s^\pm (x)`$\], are parametrized at some low (but still perturbative) scale $`Q_0`$; the full $`Q`$-dependence is then determined by DGLAP evolution.dddFor definiteness, the CTEQ PDFs evolve from $`Q_0=m_c=1.3\mathrm{GeV}`$.
* The strangeness number sum rule for the nucleon requires

$$[s^{}]=0\text{(for all }Q\text{)}.$$
(5)
A necessary corollary is that the density function $`s^{}(x)`$ must be less singular than $`1/x`$ as $`x0`$ for all $`Q`$. As DGLAP evolution preserves Eq. (5) (which is a consequence of the conservation of the strange vector current $`J_s^\mu =\overline{\psi }_s\gamma ^\mu \psi _s`$), it suffices to impose it at $`Q_0`$.
* The momentum sum rule requires

$$[S^+]=1\mathrm{\Sigma }_0\text{(for all }Q\text{)},$$
(6)
where $`\mathrm{\Sigma }_0`$ represents the momentum fraction of all non-strange partons. Through this condition the global inclusive DIS and other data—which directly affect $`\mathrm{\Sigma }_0`$—indirectly also constrain $`[S^+]`$.
* In the limit $`x0`$ (high energy and fixed $`Q`$), Regge considerations and the Pomeranchuk theorem predict $`s^{}(x)/s^+(x)0`$.

### 2.2 Expectations

From the above general constraints, we draw the following conclusions:

* The number sum rule, Eq. (5), implies that a graph of $`s^{}(x)`$ must cross the $`x`$-axis at least once in the interval $`0<x<1`$; and the areas bounded by the curve above and below the $`x`$-axis must be equal.
* Assuming a simple scenario in which there is—as supported by theoretical models  based on Lambda-Kaon fluctuations—only one zero crossing, either $`s^{}(x)<0`$ in the low $`x`$ region and $`s^{}(x)>0`$ in the high $`x`$ region, or vice versa. The two possibilities imply $`[S^{}]>0`$ or $`[S^{}]<0`$, respectively, because the momentum integral suppresses the small $`x`$ region and enhances the large $`x`$ region.
* The low-$`x`$ behaviour of $`\left|s^{}(x)\right|x^\beta _{}`$ will then be correlated with the size of $`\left|[S^{}]\right|`$: The steeper the function $`s^{}(x)`$, the larger $`\left|S^{}(x)\right|`$ has to be at large $`x`$.

To illustrate point (ii), we preview a representative CTEQ fit (to be discussed in full detail in Sec. 4) in Fig. 1 and juxtapose it with results from previous literature. According to the parametrizations of $`s(x)`$ and $`\overline{s}(x)`$ used by the CCFR-NuTeV dimuon study , $`s^{}(x)`$ is negative in the $`x`$ range covered by the experiment ($`0.01<x<0.3`$).eeeSince the preliminary version of this work was reported, the CCFR-NuTeV collaboration has emphasized that their most recent analysis favors an integrated strangeness asymmetry that is consistent with zero (K. MacFarland and P. Spentzouris, communications at the LP03 Symposium and WIN03 Workshop). More definitive studies are needed to clarify the situation. See further discussions in Sec. 5. A previous detailed global analysis of inclusive data by Barone et al.  (BPZ), finds that $`s^{}(x)`$ is positive in the large $`x`$ region. These previous results are shown in Fig. 1 as the dot-dashed and dashed curves, respectively. In light of the theoretical constraints discussed above, both these results would hint at the first possibility mentioned above, i.e. $`[S^{}]>0`$. This is not the case, however, for the CCFR-NuTeV curve, as it does not cross the $`x`$-axis; it violates the strangeness number sum rule and cannot be taken to represent the physics outside a limited window in $`x`$. The BPZ curve does satisfy the sum rule; we note (Fig. 1) it has two zero-crossings. The solid curves (Class B) in the two plots of Fig. 1 provide a concrete example of $`s^{}(x)`$ and $`S^{}(x)`$ following the above requirements and expectations.

To illustrate point (iii), we preview various classes of solutions (again, to be discussed in detail below) in the upper plot of Fig. 2. Because the experimental constraints are weak or non-existent in the very small $`x`$ region, say $`x<0.01`$, the detailed behavior of $`s^{}(x)`$ is unconstrained in this region. However, this uncertainty at small $`x`$ is considerably reduced in $`S^{}(x)`$, as demonstrated by the curves of the lower plot. Thus the above observations concerning $`[S^{}]`$ are affected only mildly by the uncertainty of the very small $`x`$ behavior, unless that behavior is so extreme that the small $`x`$ region provides a significant contribution to the number sum rule. We will refrain from exploiting such a mathematical possibility as long as it does not seem to be motivated by any physics.fffA residual bias of the results on the assumed functional forms, or rather on the rejection of some forms as merely mathematical and not physical, is unavoidable in parton analysis. One can imagine the extreme (unphysical) scenario of spike-like structures escaping “detection” in the $`x0`$ region but, nevertheless, affecting sum rules.

## 3 General parametrization of the strangeness distributions

To explore the strangeness sector of the parton structure of the nucleon, we need a suitable parametrization of $`s(x)`$ and $`\overline{s}(x)`$ (or equivalently $`s^\pm (x)`$) at a fixed scale $`Q_0`$. This parametrization must satisfy the theoretical requirements specified above, and it should be as general as possible so that the allowed functional space can be fully explored. A general form is essential, so that our conclusions are not artifacts of the parametrization, but truly reflect the experimental and theoretical constraints. In the following we explain our choice of such a parametrization. Full details including explicit parameters are given in an Appendix.

It is more natural to parametrize the $`s^\pm (x,Q_0)`$ functions independently (rather than $`s`$ and $`\overline{s}`$) since they satisfy different QCD evolution equations: pure non-singlet for $`s^{}`$ and mixed singlet/non-singlet for $`s^+`$. We use the following parametrizations,

$`s^+(x,Q_0)`$ $`=`$ $`A_0x^{A_1}(1x)^{A_2}P_+(x;A_3,A_4,\mathrm{})`$ (7)
$`s^{}(x,Q_0)`$ $`=`$ $`s^+(x,Q_0)\mathrm{tanh}[ax^b(1x)^cP_{}(x;x_0,d,e,\mathrm{})]`$ (8)
where $`P_+(x;A_3,\mathrm{})`$ is a positive definite, smooth function in the interval $`(0,1)`$, depending on additional parameters $`A_3,\mathrm{}`$ such as are used for $`u,d,g,\mathrm{}`$ in most CTEQ  and other global analyses ; and

$$P_{}(x)=\left(1\frac{x}{x_0}\right)\left(1+dx+ex^2+\mathrm{}\right)$$
(9)
where the crossing point $`x_0`$ is determined by the strangeness number sum rule $`[s^{}]=0`$, and the parameters $`d,e,\mathrm{}`$ are optional, depending on how much detail is accessible with the existing constraints. Important features of this parametrization are the following:

* The strangeness quantum number sum rule, $`[s^{}]=0`$, is satisfied by the choice of $`x_0`$. The parameter $`x_0`$ has a physical interpretation: it is the “crossing point” where $`s^{}(x)=0`$. (If $`d`$ and/or $`e`$ are not zero, there can be additional zeros of $`s^{}(x)`$; in practice—as explained in the previous section and below in Sec. 4.2—we restrict our analysis to solutions with a single crossing only.)
* The fact that the $`\mathrm{tanh}`$ function has absolute value less than 1 ensures positivity of $`s(x)`$ and $`\overline{s}(x)`$. The fact that $`\mathrm{tanh}`$ is a monotonic function guarantees that the function $`s^{}(x)`$ can be made as general as necessary by the choice of $`P_{}(x)`$.
* The small-$`x`$ behavior of $`s^{}(x)`$ must be such that the integral $`[s^{}]`$ converges (before the root $`x_0`$ is determined). Let $`\beta _{}A_1+b`$; then Eq. (8) implies

$$s^{}(x)x^\beta _{}\mathrm{as}x0.$$
(10)
The convergence of $`[s^{}]`$ is guaranteed if $`\beta _{}>1`$, i.e. the parameter $`b`$ is chosen in the range $`b>1A_1`$.

Because $`P_\pm (x)`$ can be made as general as necessary, the choice in Eqs. (7)–(9) is capable of exploring the full strangeness parameter space allowed by data in the PQCD framework.

Detailed formulas for all flavors used in the actual analysis are presented in the Appendix.

## 4 Global Analysis

We now describe the global QCD analysis, which includes all relevant experimental data and implements the theoretical ideas outlined above. This may be considered an extension of the on-going CTEQ program of global analysis. Several new elements (compared to the latest CTEQ6M  analysis) are present. On the experimental side, we have added the CDHSW inclusive $`F_2`$ and $`F_3`$ data sets , and the CCFR-NuTeV dimuon data sets . On the theoretical side, we have expanded the parameter space to include the strangeness sector as discussed in Sec. 3.

Compared to the global analyses of BPZ , which also allow $`s\overline{s}`$, the major difference experimentally is our inclusion of the dimuon data, which provide a direct handle on $`s`$ and $`\overline{s}`$; and, theoretically, the generality and naturalness of our parametrization of the strange distributions.gggRef.  parametrizes $`s(x)`$ and $`\overline{s}(x)`$ rather than $`s^\pm (x)`$. Since the results of  on strangeness asymmetry rely on small differences of inclusive DIS charged-current and neutral-current measurements, BPZ performed the analysis at the cross section level, applying uniform procedures to treat data from different experiments in the comparison to theory. Considering small differences between inclusive cross sections, the possible strange asymmetry is but one of many sources that could lead to such differences.

As the dimuon data more directly constrain the strange PDFs, this is an important new element to our fit. For all inclusive DIS processes we use the standard procedure of comparing theory with the published $`F_2`$ and $`F_3`$ structure function data. In our analysis, the fit to charged-current (neutrino) inclusive structure functions is dominated by the high statistics CCFR data. Although we have included the earlier inclusive CDHSW data (which play a prominent role in the analysis of ), they have no discernible influence on the results presented below.

To include the CCFR-NuTeV neutrino and antineutrino dimuon production data in a global QCD analysis is not a straightforward task. The experimental measurement is presented as a series of “forward differential cross sections” with kinematic cuts, whereas the theoretical quantities that are most directly related to the parton distribution analysis are the underlying (semi-) inclusive “charm quark production cross sections”. The gap between the two is bridged using a Monte Carlo program that incorporates kinematic cuts as well as fragmentation and decay models. In our analysis, we use a Pythia program provided by the CCFR-NuTeV collaboration to do this efficiency-correction.hhhWe thank Tim Bolton and Max Goncharov, in particular, for providing this program, as well as assistance in its use. Their help was vital for carrying out this project. This Monte Carlo calculation is done in the spirit and the framework of leading-order (LO) QCD. CTEQ5L parton distributions and Peterson fragmentation functions were used. The parameters of the model were tuned to reproduce, as closely as possible, the detailed differential dimuon cross sections published in .

### 4.1 $`d\sigma ^{\nu NcX}`$ in QCD

At LO in PQCD, the cross section formula for $`\nu NcX`$ is

$`\xi s^{}(\xi ,Q^2)_{\mathrm{eff}}`$ $``$ $`{\displaystyle \frac{1}{2}}{\displaystyle \frac{\pi (1+Q^2/M_W^2)^2}{G_F^2M_NE_\nu }}\left|V_{cs}\right|^2{\displaystyle \frac{d^2\sigma ^{\nu NcX}}{dxdy}}`$ (11)
$`=`$ $`(1{\displaystyle \frac{m_c^2}{2M_NE_\nu \xi }})\xi s^{}(\xi ,Q^2)+𝒪(\alpha _s),`$
with the CKM matrix element $`\left|V_{cs}\right|`$ and where the Barnett-Gottschalk parameter

$$\xi x\left(1+\frac{m_c^2}{Q^2}\right)$$
(12)
approaches Bjorken-$`x`$ as $`Q\mathrm{}`$ (relative to $`m_c=1.3\mathrm{GeV}`$). The quantity $`s^{}(\xi ,Q^2)_{\mathrm{eff}}`$ in Eq. (11) includes Cabibbo suppressed contributions in neutrino scattering via

$$s^{}s+\frac{\left|V_{cs}\right|^2}{\left|V_{cd}\right|^2}d$$
(13)
with obvious adjustments for the anti-neutrino caseiii$`q\overline{q}`$ for $`q=s,d`$..

The NLO corrections to $`\xi s^{}(\xi ,Q^2)_{\mathrm{eff}}`$, defined via the perturbative series

$$\xi s^{}(\xi ,Q^2)_{\mathrm{eff}}=\underset{i}{}\alpha _s^i\xi s^{}(\xi ,Q^2)_{\mathrm{eff}}^{(i)},$$
(14)
are generically of the form

$$\xi s^{}(\xi ,Q^2)_{\mathrm{eff}}^{(1)}\underset{f=g,s^{}}{}fH_f$$
(15)
with $``$ denoting a convolution integral over parton momentum. These were first calculated more than 20 years ago . Later calculations  corrected minor typos and employed the modern $`\overline{\mathrm{MS}}`$ renormalization scheme and the ACOT treatment  of amplitudes with massive quarks ($`m_{s,c}0`$). Very recently, the NLO charm production contributions to the full set of electroweak structure functions were calculated , including terms that are suppressed by $`m_\mu ^2/ME_\nu `$. In order to apply detector acceptance corrections to the data , differential NLO distributions were calculated in  and  that provide the charm hadron (D meson) kinematics in terms of the fragmentation variable $`z`$ and rapidity $`\eta `$. The $`d\sigma /dxdydzd\eta `$ code DISCO , written by two of the authors of the present article in collaboration with D. Mason of NuTeV, exists as an interface to the NuTeV MC event generator. Detailed results can be found in the articles listed above. It suffices to say: (i) the NLO calculations all agree; and (ii) for the fixed target kinematics under investigation, the NLO corrections to the LO results are modest—no bigger than $`20\%`$ (see Fig. 1 in ).

As mentioned before, the global fits performed in our study are extensions of the full NLO CTEQ6 analysis with the addition of constraints due to neutrino dimuon production. For the latter process, we have done extensive studies using either the LO formula, Eq. (11), or the NLO treatment of , Eq. (15). The results obtained in the two cases are quite similar. For definiteness, the main results presented in Sec. 4.3 are those obtained by using Eq. (11), since the acceptance corrections made to the data set are currently based on a LO model, as mentioned earlier. Since we have determined that the NLO corrections to the hard cross section are small (compared to experimental errors, for instance), and since we found the uncertainty range of the main result (on $`[S^{}]`$ ) is much broader than the difference between the central values obtained by using Eq. (11) with or without the corrections in Eq. (15) (cf. Sec. 5), this approximation does not affect the outcome of our analysis.jjjIt is certainly desirable to have the inclusive cross sections corrected for acceptance based on NLO models (such as ), that can be compared to (15) in a full NLO global analysis. This is under active development by a theory (CTEQ)-experiment (CCFR-NuTeV) collaboration.

### 4.2 Procedure

Our analysis is carried out in several stages. First we must find appropriate starting values for the fitting parameters. For this purpose, we implement the following steps.

* We rerun the CTEQ6M global fit with the added CDHSW inclusive neutrino scattering data, keeping all other conditions the same. This intermediate fit is extremely close to the CTEQ6M one, since the fit to inclusive DIS data is totally dominated by the high statistics neutral current experiments on the one hand, and the CCFR charged current experiment on the other.
* We then fix all of the “conventional” parton parameters to their values in this intermediate fit, and fit the complete set of data, including the new dimuon data, by varying only the parameters associated with the new degrees of freedom in $`s^{}`$. We obtain results consistent with expectations:

+ Most of the data sets used in the previous analysis are not affected at all by the variation in $`s^{}`$.
    + A few fully inclusive cross sections are slightly affected by the variation of $`s^{}`$, mainly:

    - $`F_3`$ which depends on $`u\overline{u}+d\overline{d}+s\overline{s}\mathrm{}`$.
        - The $`W^\pm `$ charge asymmetry which receives a contributions from $`gsW^{}c`$.

    These sensitivities to $`s^{}`$ are weak.
    + The CCFR-NuTeV dimuon data sets are the most constraining ones for fitting $`s^{}(x)`$.

We obtain good fits using either the 3-parameter ($`a,b,c`$) or the 4- or 5-parameter ($`a,b,c,d,e`$) versions of Eqs.(8,9). There are not enough constraints to choose among these fits. The higher-order polynomials allow oscillatory behavior of $`s^{}(x)`$ which the 3-parameter form does not. We consider the number of crossings a distinctive property of the physical asymmetry $`s^{}(x)`$ rather than a volatile function of the continuous fit parameter space. As explained in Section 2, one zero-crossing is unavoidably enforced by the sum rule in Eq. (5). We are not aware of any solid theoretical argument that would suggest a second crossing; nor do we find that the fits show any preference for more than one crossing. We therefore restrict the search for best fits in this section to one crossing. This choice also seems to be a stable feature of models  based on baryon-meson fluctuations.
* Using these candidate fits as a basis, we perform a second round of fitting allowing the parameters associated with $`s^+`$, Eq.(7), to vary in addition to the $`s^{}`$ variables. This improves the fit to all data sets slightly. We observe that the shape of $`s^+(x)`$ now deviates from the starting configuration in which $`s^+(x)`$ was set proportional to $`\overline{u}(x)+\overline{d}(x)`$. Defining, as in Eq. (1), the strangeness suppression parameter $`\kappa `$ as the ratio of the momentum fraction carried by the strange quarks, $`[S^+]`$, to that carried by $`\overline{u}+\overline{d}`$ at $`Q_0=1.3\mathrm{GeV}`$, we find that $`\kappa `$ may vary in the range $`0.3`$ – $`0.5`$: $`\chi ^2`$ has a shallow minimum around $`\kappa =0.4`$. This value agrees with previous analyses .

Because the experimental constraints are not sufficient to uniquely determine all the $`s^{}`$ and $`s^+`$ parameters, we categorize several classes of equally good solutions based on the behavior of $`s^{}(x)/s^+(x)`$ as $`x0`$ or $`x1`$.
* We finalize these classes of solutions by allowing all parton parameters to vary so that the non-strange parton distributions can adjust themselves to yield the best fit to all the experimental data sets. (As one would expect, these final adjustments are generally small.) The differences in the $`\chi ^2`$ values between the various categories of solutions are not significant; i.e. we find nearly degenerate minima with distinctively different $`s^{}(x)`$ solutions, classified according to their small $`x`$ behaviour. These solutions do not correspond to isolated local minima in $`\chi ^2`$ space; rather they are to be thought of as specific examples of a class of acceptable fits that lie along a nearly flat “valley” along which $`\chi ^2`$ changes very slowly.

### 4.3 Central Results

The quality of the fits to the global data sets other than the CCFR-NuTeV dimuon data remains unaltered from the previous CTEQ6M analysis, so we focus our discussion on the strangeness sector. Specifically, we examine closely the asymmetry functions $`s^{}(x),S^{}(x)`$ and the momentum integral $`[S^{}]`$. The asymmetry functions from three typical good fits, with different behaviors at small $`x`$ (labeled as classes A,B,C), were previewed in Fig. 2 as illustrations.

In the accompanying table, for each sample fit we list the small-$`x`$ exponent $`\beta _{}`$ \[$`s^{}(x)x^\beta _{}`$, cf. Eq. 10\], the integrated momentum fraction $`[S^{}]`$, and the relative $`\chi ^2`$ values, normalized to the $`\chi ^2`$ of solution B ($`\chi _B^2`$), which we use as the reference for comparison purposes. (Under column “B”, we give the absolute $`\chi ^2`$’s in parentheses.)kkkThe $`\chi ^2`$ values of the dimuon data sets, like those of some other data sets, do not carry rigorous statistical significance, because the correlated systematic errors are not available and, hence, cannot be included. In the global analysis context, the $`\chi ^2`$ value is nevertheless used as the only practical “figure of merit” for the fit. The relatively small value of the total $`\chi ^2`$ for the dimuon data sets, compared to the number of data points, underlines this fact. Under this circumstance, it is common practice to use the normalized $`\chi ^2`$ values to compare the quality of different fits. To gain some insight on the constraints on the strangeness sector due to the various types of experiments, we show separately the $`\chi ^2`$ values for the dimuon data sets, the inclusive data sets (I) that are expected to be somewhat sensitive to $`s^{}`$ (consisting of the CCFR and CDHSW $`F_3(x,Q)`$ and the CDF $`W`$-lepton asymmetry measurements), and the remaining data sets (II) that are only indirectly affected by $`s^{}`$ (the rest of the inclusive data sets).

|  | \# pts | B+ | A | B | C | B$``$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`\beta _{}`$ | - | $`0.78`$ | $`0.99`$ | $`0.78`$ | 0 | $`0.78`$ |
| $`[S^{}]\times 100`$ | - | 0.540 | 0.312 | 0.160 | 0.103 | $`0.177`$ |
| Dimuon | 174 | 1.30 | 1.02 | *1.00* (126) | 1.01 | 1.26 |
| Inclusive I | 194 | 0.98 | 0.97 | *1.00* (141) | 1.03 | 1.09 |
| Inclusive II | 2097 | 1.00 | 1.00 | *1.00* (2349) | 1.00 | 1.00 |

Table 1. The representative parton distribution sets, arranged in order by the value of $`[S^{}]`$.

Focusing on the three good fits {A, B, C} first, we note the following features:

* All three solutions {A, B, C} feature positive $`[S^{}]`$; and the more singular the behavior of $`s^{}(x)`$ as $`x0`$, the higher the value of $`[S^{}]`$. These are natural consequences of the strangeness sum rule (equal $`+/`$ areas under the curve of $`s^{}(x)`$) and the small $`x`$ suppression of the momentum integral, as discussed earlier in Sec. 2.
* Solution B is slightly favored over the other two. This, plus the fact that its small-$`x`$ behavior lies in the middle of the favored range, motivates its use as the reference fit.
* We chose these examples among fits with the simplest parametrizations: all cross the $`x`$ axis only once. With 4- or 5-parameters, which can allow more than one crossing point, many solutions can be found that entail oscillatory $`s^{}(x)`$. But since the $`\chi ^2`$ values are essentially the same as for the simple case, we deem it premature to dwell on complicated behaviors, which may be mere artifacts of the parametrization rather than reflections of physical constraints. Further studies described in Sec. 4.5 reinforce this point.

To show how these fits compare with data, we plot in Fig. 3 the ratio of data/theory for the reference fit B. The four graphs correspond to the CCFR and NuTeV neutrino and antineutrino data sets respectively. The data points are sorted in $`x`$-bins, and within each $`x`$-bin, by $`y`$ value. We see that the quality of the fit is good, within the experimental uncertainties. There are no significant systematic deviations. (The CCFR antineutrino data set may appear to be systematically higher than theory. However, upon closer inspection the difference is not significant. The data points that lie above theory consist mostly of points with large error bars, which tend to catch the attention of the eye; whereas the fit is actually dominated by points with small errors, which closely bracket the theory line on both sides.lllThis becomes apparent if the data points are re-plotted ordered by the size of the error bars. The value of $`\chi ^2/N`$ for this data set is less than 1, comparable to those for the other sets.)

The parameters for all the fits described in this section are given in detail in the Appendix. As is already obvious from Fig. 2, solutions with nearly degenerate $`\chi ^2`$ may correspond to parametrizations of $`s^{}(x)`$ with quite different parameter values, so that a simple linear error analysis cannot be applied. This reflects the fact that $`s^{}(x)`$ is not well determined as a detailed function of $`x`$, even when the dimuon data are included in the fit. On the other hand, reducing the parameter space to even fewer parameters than our minimal set would risk introducing artifacts of an inflexible parametrization. In the next section we will, therefore, apply the Lagrangian multiplier method to deduce the integrated momentum asymmetry $`[S^{}]`$ and its uncertainty.

### 4.4 Range of $`[S^{}]`$ by the Lagrange Multiplier Method

Beyond the best fits (A, B, C), we can study the range of $`[S^{}]`$ consistent with our global analysis in a quantitative way by applying the Lagrange Multiplier (LM) method developed in . By varying the Lagrange multiplier parameter, this method explores the entire strangeness parameter space in search of solutions with specified values of $`[S^{}]`$, i.e., constrained fits. The B- solution listed in Table 1 was obtained by forcing $`[S^{}]=0.0018`$ (a relatively large negative value, but not as large as the value $`0.0027`$ quoted by ). The B+ solution was generated by forcing $`[S^{}]`$ to go in the other (positive) direction until the increment of the overall $`\chi ^2`$ became comparable to that of B-; this results in $`[S^{}]=0.0054`$.

We see from the relevant entries in Table 1 that: (i) the $`\chi ^2`$ values of the dimuon data sets increase by about 30% in both B± fits; (ii) the “inclusive I” data sets disfavor the negative $`[S^{}]`$; and (iii) the “inclusive II” data sets are completely neutral. These results are shown graphically in Fig. 4a, where the square points represent the (relative) $`\chi ^2`$ values of the dimuon data sets, and the triangle points of the “inclusive I” data sets. (Not shown are those for the “inclusive II” data sets, which remain flat (at $`1.00`$).) The LM fits are chosen from a large number of fits spanning the entire strangeness parameter space. The pattern of dependence of the $`\chi ^2`$ values for the dimuon data sets on the value of $`[S^{}]`$ is nearly parabolic. This is clear evidence that the dimuon measurement is indeed sensitive to the strangeness asymmetry as expected. Further discussion of this observation, including the contrast to the sensitivity of other experiments, will be given in Sec. 4.5.

We see from Fig. 4a that, in this series of fits, the dimuon data sets favor a range of $`[S^{}]`$ centered around 0.0017, whereas the “inclusive I” data sets disfavor negative values of $`[S^{}]`$. Fig. 4b shows the dependence of the combined $`\chi ^2`$ of the two categories of data sets on $`[S^{}]`$. (The $`\chi ^2`$ of the remaining data sets used in the global analysis are totally insensitive to $`[S^{}]`$, cf. Table 1, hence is not included in this plot.) We would like to determine a “range of uncertainty” of $`[S^{}]`$ from these results. This is far from straightforward because of well-known problems shared by all error assessments in global analysis (mainly due to the unquantified systematic errors that show up as a lack of statistical compatibility among the input data sets).mmmThese difficulties, and practical methods to handle them, are discussed in detail in .

One naive method is to apply the $`\mathrm{\Delta }\chi ^2=1`$ criterion. From the parabola in Fig. 4b, which comes from 368 data points, this “estimation-of-parameters criterion” corresponds to an uncertainty of $`[S^{}]`$ of $`\pm 0.0005`$. (Cf. the lowest horizontal line in Fig. 4b.) It has been known, however, that the $`\mathrm{\Delta }\chi ^2=1`$ criterion is unrealistic in global analysis when combining data sets with diverse systematic errors from many different experiments ; in this circumstance, the overall $`\chi ^2`$ function provides a simple measure of relative goodness-of-fit in the minimization process, but it does not have the strict statistical significance of a pure parameter fitting problem as presented in textbook examples. This estimate of the uncertainty of $`[S^{}]`$ is far too small.

Another often-used method to evaluate ”goodness-of-fit” is to apply the cumulative distribution function $`P`$ for the $`\chi ^2`$ distribution. One considers unacceptable values of $`\chi ^2`$ greater than $`\chi _{68}^2`$ (or $`\chi _{90}^2`$) where $`P(\chi ^2<\chi _f^2)=f`$. For 386 data points, the 68% (90%) criterion corresponds to $`\chi ^2/\chi _{min}^2=1.033(1.1)`$ respectively. These two criteria are represented by the two upper horizontal lines in Fig. 4b. The uncertainty range of $`[S^{}]`$ for these two cases are $`\pm 0.002(0.003)`$ respectively.

The extensive studies on quantifying uncertainties in the global analysis context  suggest that for this case, a realistic range should be somewhere between the two extreme cases shown in Fig. 4b. Hence we adopt the uncertainty range $`0<[S^{}]<0.004`$ by this analysis, which corresponds to the middle horizontal line in the graph (or the $`\chi _{68}^2`$ criterion). Whereas a very small strangeness asymmetry, consistent with zero, is not ruled out by this criterion, large negative values of $`[S^{}]`$ (such as -0.0027, cited in ) are strongly disfavored; cf. also Table 1. Additional sources of uncertainty will be discussed in the following section.

### 4.5 Additional Sources of Uncertainty

We have performed three series of studies to further assess the rebustness of our main results. These will help us to determine a better estimate of the overall uncertainty of $`[S^{}]`$.

#### Pure Leading Order Fits

Since the experimental analyses of the CCFR-NuTeV dimuon data have been done in LO QCD , we have carried out a whole series of purely LO global analyses, following the same procedures as outline above, in order to provide a basis for comparison. The results can be summarized as follows.

* The overall $`\chi ^2`$ for the global fit increased by $`200`$ over the comparable fits described above; while the $`\chi ^2`$’s for the dimuon data sets actually decreased slightly. This is not surprising, since the current state of global analysis, with precision data from many experiments, requires the use of NLO QCD theory – in particular for the collider data with typically large perturbative corrections. On the other hand, the new dimuon data still have comparably large experimental errors and the NLO corrections are small \[$`𝒪(20\%)`$\], such that an LO fit is adequate for them.
* We explored the allowed range of strangeness asymmetry $`[S^{}]`$ in this LO study under different assumptions on the $`x0`$ and $`x1`$ behavior of the $`s^+(x)`$ and $`s^{}(x)`$ functions. First, we found that the dependence of $`\chi _{\mathrm{dimuon}}^2`$ on $`[S^{}]`$ is generally parabolic, rather similar to Fig. 4. The width of the distribution is comparable to Fig. 4. The central value for $`[S^{}]`$ is within the range $`(0,0.0015)`$; the exact value depends on the $`x0`$ and $`x1`$ behavior of $`s^\pm (x)`$ assumed.
* We also found that the $`\chi _{\mathrm{inclusive}\mathrm{I}}^2`$ vs. $`[S^{}]`$ curve, while generally flatter, does “flop around” enough for the cases studied so that no clear pattern can be discerned. The specific shape of this curve shown in Fig. 3 is not a common characteristic of these fits.

#### Charm Mass Dependence

The CCFR-NuTeV dimuon analysis treated the charm mass as one of the fit parameters. Their analyses favored a rather high value of $`m_c=1.6`$ GeV (compared to, e.g., the PDG estimate of $`1.0\mathrm{GeV}<m_c<1.4\mathrm{GeV}`$). The CTEQ global analyses are usually done with a fixed value of $`m_c=`$ 1.3 GeV. To see whether the comparison between our results is strongly influenced by the choice of the charm mass, we have performed several series of fits with $`m_c`$ varying from 1.3 GeV to 1.7 GeV. Again, the general features, as described above, stay the same. The central value of $`[S^{}]`$ does vary with the choice of $`m_c`$ within a given series of fits, but the pattern is not universal. The range over which the central value wanders is of the order $`0.0015`$, comparable to the width of the parabola in Fig. 4. Unlike the specific analysis of CCFR-NuTeV, the overall $`\chi ^2`$ for the global analysis does favor a lower value of $`m_c`$.

#### Dependence on Decay and Fragmentation Model

To estimate the dependence of our results on the model used to convert the measured dimuon cross sections to structure functions for charm production, we repeated our analyses using an alternative conversion table provided by the CCFR-NuTeV collaboration.nnnWe thank Kevin MacFarland for supplying this table. This alternative table is based on Buras-Gaemers PDFs used in CCFR-NuTeV analyses with Collins-Spiller fragmentation functions. It is similarly tuned to detailed features of the measured dimuon cross sections as that described in Sec. 4.2.oooHowever, since our CTEQ6-like PDFs are rather different from the CCFR Buras-Gaemers PDFs, it is not clear how good the approximation is to use this conversion table. That is, the self-consistency of the procedure is not assured. The results obtained from the alternative fits are, again, similar to those described earlier. The $`\chi _{\mathrm{dimuon}}^2`$ vs. $`[S^{}]`$ parabola generally has the same width as in Fig. 4. The central value for $`[S^{}]`$ is in the range $`0<[S^{}]<0.0015`$—in the lower half of the range quoted at the end of Sec. 4.4. The dependence of $`\chi _{\mathrm{inclusive}\mathrm{I}}^2`$ on $`[S^{}]`$ shows no definitive trend.

Taken together, the results of the additional studies described in these three paragraphs lead to several conclusions. (i) The general features described in Secs. 4.3 and 4.4 are robust; (ii) The central value of $`[S^{}]`$ wanders around within a range that is consistent with the width of the $`\chi _{\mathrm{dimuon}}^2`$ vs. $`[S^{}]`$ parabola; and (iii) These additional results do not significantly change the estimates of the previous section, except to shift the central estimated value of $`[S^{}]`$ to a slightly lower value, and to extend the range of uncertainty on the lower side somewhat. The envelope of these additional uncertainties provides an estimated range of uncertainty of the strangeness asymmetry ofpppThe uncertainties from the various sources need not be combined in quadrature, because they are not statistically independent sources, but rather systematic uncertainties of the theory.

$$0.001<[S^{}]<0.004.$$
(16)
This large range reflects both the limit of current experimental constraints and the considerable theoretical uncertainty, as explicitly discussed in the text. The theoretical uncertainties can be reduced in a refined NLO analysis; the results remain to be seen. The limitations on the experimental constraints will remain, until new experiments are done .

## 5 Comparisons to previous studies

A comprehensive global QCD analysis with emphasis on the strangeness sector has been carried out previously by BPZ .qqqAs mentioned in Sec. 4, BPZ work directly with DIS cross sections (instead of structure functions), with detailed attention to systematic errors and other sources of uncertainties. Without the dimuon data, which are directly sensitive to strangeness, the results of BPZ implicitly rely on small differences between large neutral- and charged-current inclusive cross sections from different experiments. The latest representative $`s^{}(x)`$ and $`S^{}(x)`$ functions extracted by BPZ are shown in Fig. 1, along with our reference fit B. The main feature of the BPZ curves is a positive bump at rather large $`x`$.rrr The BPZ curve displayed includes the more recent analysis “with CCFR (inclusive data).” The original “without CCFR” solution in  has a more pronounced large-x bump and a smaller negative region. An $`s^{}(x)`$ function of such magnitude at large $`x`$ is in disagreement with the CCFR-NuTeV analysis; in addition, it would make $`s(x)`$ and $`\overline{s}(x)`$ behave quite differently compared with the non-strange sea quarks and the gluon, i.e., $`(1x)^p`$, with $`p`$ in the range $`510`$. This feature has been attributed to the influence of the CDHSW data, particularly when re-analyzed at the cross section level along with the other DIS experiments. Their conclusion that data favor a positive value of the momentum integral $`[S^{}]`$ is in general agreement with our detailed study based on the LM method. However, the different shapes of $`s^{}(x)`$ seen in Fig. 1 clearly underline the difference in inputs: (i) our results are mainly dictated by the CCFR-NuTeV dimuon data (which are not present in the BPZ analysis); (ii) their results rely on a delicate analysis of DIS cross-section data (not matched in our structure function analysis); and (iii) the difference in flexibility of the parametrizations of the non-perturbative input functions can influence the results.

The CCFR and NuTeV collaborations performed separate and combined analyses of $`s`$ and $`\overline{s}`$ , based on their own dimuon and inclusive cross sections. To parameterize the $`s(x)`$ and $`\overline{s}(x)`$ distributions, they chose the model formula

$$\left(\begin{array}{c}s(x,Q)\\ \overline{s}(x,Q)\end{array}\right)=\frac{\overline{u}(x,Q)+\overline{d}(x,Q)}{2}\left(\begin{array}{c}\kappa (1x)^\alpha \\ \overline{\kappa }(1x)^{\overline{\alpha }}\end{array}\right)$$
(17)
for all $`\{x,Q\}`$, where $`\kappa ,\overline{\kappa },\alpha ,\overline{\alpha }`$ are fitting parameters.

Curves representing the general behavior of the model (17) at $`Q^2=10`$ GeV2, with parameter values ($`\kappa ,\overline{\kappa },\alpha ,\overline{\alpha }`$) taken from , have been shown in Fig. 1 for comparison with the other distributions. While this model might be acceptable for a limited range of $`x`$ and $`Q`$, it leads to problems when extrapolated to general $`\{x,Q\}`$ values: (i) the strangeness number sum rule $`[s^{}]=0`$ is violated (in fact, the integral $`[s^{}]`$ diverges unless $`\kappa =\overline{\kappa }`$), as is the momentum sum rule, Eq. 6; (ii) the QCD evolution equations are violated by the $`Q`$-dependent paramaterizations of the PDFs.sssThese problems are independent of the specific issues of the strange quark asymmetry. The first problem is evident in Fig. 1.tttRe-analysis of the CCFR-NuTeV data by the experimental group, taking into account these issues, are underway. Initial results from partial implementations of the above-mentioned theoretical constraints were reported by P. Spentzouris at *International Workshop on Weak Interactions and Neutrinos 2003* (WIN03), Lake Geneva, October, 2003.

It could be argued that, since the experiment only covers a limited range of $`x`$, the enforcement of the sum rules is not critical in extracting limited information on $`s(x,Q)`$ and $`\overline{s}(x,Q)`$. Given this,the CCFR curve in Fig. 1 implies negative $`s^{}(x)`$ over most of the experimental $`x`$ range (0.01 – 0.3).uuuThis is consistent with the fact that (without the constraint of sum rules)  quotes $`[S]=0.0027\pm 0.0013`$. The uncertainty is smaller at the lower $`x`$ end because of better statistics. It is this feature of the data that we invoked in the general discussion of Sec. 2.

To describe the behavior of $`s^{}(x)`$ over the full $`x`$ range, the strangeness number sum rule must be enforced. The sum rule also provides a powerful theoretical constraint on the data analysis, as demonstrated in Secs. 2 and 4. For functions to be candidate *universal* parton distributions of the PQCD formalism, they must satisfy both the sum rules and the QCD evolution equations.

## 6 Conclusion

We find a range of solutions in the strangeness sector that are consistent with all relevant world data used in the global analysis. The dimuon data are vital in constraining the strangeness asymmetry parameters. The constraints provided by other inclusive measurements, labeled as “inclusive I” in the text, are consistent with those provided by dimuon data, although much weaker. The allowed solutions generally prefer the momentum integral $`[S^{}]_0^1x[s(x)\overline{s}(x)]𝑑x`$ to be positive. This conclusion is quite robust, and it follows from the basic properties of PQCD and from qualitative features of the experimental data. However, the size of this strangeness momentum asymmetry is still quite uncertain; we can only estimate that $`\left[S^{}\right]`$ lies in the range from $`0.001`$ to $`+0.004`$. The Lagrange Multiplier method explicitly demonstrates that both the dimuon data and the “inclusive I” data sets strongly disfavor a large negative value of $`[S^{}]`$, although they may still be consistent with zero asymmetry.

The fact that $`\left[S^{}\right]`$ has a large uncertainty has significant implications for the precision measurement of the weak mixing angle, $`\mathrm{sin}^2\theta _W`$, from neutrino scattering. This issue is studied separately in Ref. .

This paper marks the first global QCD analysis incorporating direct experimental constraints on the strangeness sector. We have so far focused only on the strangeness asymmetry, which represents a new frontier in parton degrees of freedom. Much still needs to be done to improve the treatment of the dimuon data (to true NLO accuracy), and to fully explore all the degrees of freedom associated with $`s^+(x)`$ and $`s^{}(x)`$. As progress is made on these fronts, the uncertainty on $`[S^{}]`$ will no doubt decrease as well.

Note added: After this manuscript was completed an investigation of the 3-loop perturbative strangeness asymmetry was presented in .

Acknowledgment

We thank members of the CCFR and NuTeV collaboration, in particular T. Bolton and M. Goncharov, for discussions about the dimuon data, and assistance in their use, and Kevin MacFarland, Dave Mason, and Panagiotis Spentzouris for useful comments and suggestions. We also thank Benjamin Portheault for discussions of the BPZ work and for assistance in generating their results for comparison. F.O. acknowledges the hospitality of MSU and BNL where a portion of this work was performed. This research was supported by the National Science Foundation (grant No. 0100677), RIKEN, Brookhaven National Laboratory, and the U.S. Department of Energy (Contract No. DE-AC02-98CH10886, and No. DE-FG03-95ER40908), and by the Lightner-Sams Foundation.

## Appendix

In this appendix, we provide some detailed information on the parametrization of the non-perturbative parton distribution functions used in our global analysis work, as well as the values of the parameters for the representative fits.

The non-strange parton distributions at $`Q=Q_0=1.3`$ GeV are parametrized as in CTEQ6. We use

$$xf(x,Q_0)=A_0x^{A_1}(1x)^{A_2}e^{A_3x}(1+e^{A_4}x)^{A_5}$$
(18)
with independent parameters for parton flavor combinations $`u_vu\overline{u}`$, $`d_vd\overline{d}`$, $`g`$, and $`\overline{u}+\overline{d}`$. To distinguish the $`\overline{d}`$ and $`\overline{u}`$ distributions, we parametrize the *ratio* $`\overline{d}/\overline{u}`$, as a sum of two terms:

$$\overline{d}(x,Q_0)/\overline{u}(x,Q_0)=A_1x^{A_2}(1x)^{A_3}+(1+A_4x)(1x)^{A_5}.$$
(19)
The strangeness sector is parametrized according to the description of Sec.(3) with

$$xs^+(x,Q_0)=A_0x^{A_1}(1x)^{A_2}e^{A_3x}(1+e^{A_4}x)^{A_5}$$

where the $`A_{1,2,3,4,5}`$ coefficients are either equated to those of $`\overline{u}+\overline{d},`$ or allowed to vary independently, depending on the particular fit being performed; and

$$s^{}(x,Q_0)=s^+(x,Q_0)\mathrm{tanh}[A_0x^{A_1}(1x)^{A_2}\left(1\frac{x}{A_3}\right)\left(1+A_4x+A_5x^2\right)]$$

As an example, the “standard fit” (B-fit) described in Sec.4.3 has the following coefficients.

B-fit $`A_0`$ $`A_1`$ $`A_2`$ $`A_3`$ $`A_4`$ $`A_5`$    $`d_v`$ 1.55891 0.62792 5.08865 0.70688 -0.26338 3.00000  $`u_v`$ 1.73352 0.55260 2.90090 -2.63846 1.45774 1.85987  $`g`$ 31.24912 0.52014 2.38230 4.23010 2.33765 -3.00000  $`\overline{d}/\overline{u}`$  10.23811 5.19599 14.85860 17.00000 8.68941  $`\overline{u}+\overline{d}`$ 0.06758 -0.29681 7.71700 -0.82223 4.43481 0.66711  $`s+\overline{s}`$ 0.04966 0.03510 7.44149 -1.44570 5.13400 0.59659  $`s\overline{s}`$ -0.08438 0.18803 2.36708 0.04590 0.00000 0.00000

And the other sample fits have coefficients:

| A-fit | $`A_0`$ | $`A_1`$ | $`A_2`$ | $`A_3`$ | $`A_4`$ | $`A_5`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`d_v`$ | 1.55485 | 0.62740 | 5.09550 | 0.70010 | -0.25130 | 3.00000 |
| $`u_v`$ | 1.73345 | 0.55260 | 2.90090 | -2.63740 | 1.45810 | 1.85910 |
| $`g`$ | 31.05518 | 0.51840 | 2.38230 | 4.24730 | -2.33790 | 3.00000 |
| $`\overline{d}/\overline{u}`$ |  | 10.23920 | 5.19690 | 14.85860 | 17.00000 | 8.68530 |
| $`\overline{u}+\overline{d}`$ | 0.06677 | -0.29810 | 7.71700 | -0.80540 | 4.45630 | 0.66420 |
| $`s+\overline{s}`$ | 0.05133 | 0.06567 | 7.59880 | -1.43972 | 5.13400 | 0.61715 |
| $`s\overline{s}`$ | -0.02000 | -0.06562 | -0.48403 | 0.01987 | -0.58837 | 0.00000 |

| C-fit | $`A_0`$ | $`A_1`$ | $`A_2`$ | $`A_3`$ | $`A_4`$ | $`A_5`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`d_v`$ | 1.39587 | 0.60594 | 4.75379 | -0.75291 | 0.31387 | 3.00000 |
| $`u_v`$ | 1.72080 | 0.55260 | 2.90090 | -2.36331 | 1.60950 | 1.59812 |
| $`g`$ | 29.40547 | 0.50423 | 2.38230 | 4.41409 | 2.34266 | -3.00000 |
| $`\overline{d}/\overline{u}`$ |  | 10.23202 | 5.19618 | 14.85860 | 17.00000 | 8.62620 |
| $`\overline{u}+\overline{d}`$ | 0.06479 | -0.30038 | 7.71700 | -0.68604 | 4.61325 | 0.62988 |
| $`s+\overline{s}`$ | 0.04289 | 0.00809 | 7.71700 | -1.22365 | 5.13400 | 0.62988 |
| $`s\overline{s}`$ | -2.36926 | 0.99191 | 8.23499 | 0.07884 | 0.00000 | 0.00000 |

| B+-fit | $`A_0`$ | $`A_1`$ | $`A_2`$ | $`A_3`$ | $`A_4`$ | $`A_5`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`d_v`$ | 1.43611 | 0.61170 | 4.73270 | -0.67420 | 0.24920 | 3.00000 |
| $`u_v`$ | 1.71921 | 0.55260 | 2.90090 | -2.39220 | 1.60430 | 1.61490 |
| $`g`$ | 29.76781 | 0.50800 | 2.38230 | 4.35570 | -2.33860 | 3.00000 |
| $`\overline{d}/\overline{u}`$ |  | 10.19560 | 5.16810 | 14.85860 | 17.00000 | 8.69760 |
| $`\overline{u}+\overline{d}`$ | 0.06729 | -0.29650 | 7.71700 | -0.75670 | 4.52290 | 0.64380 |
| $`s+\overline{s}`$ | 0.03456 | 0.00210 | 8.23420 | -1.26970 | 5.13410 | 0.72501 |
| $`s\overline{s}`$ | -0.40480 | 0.22103 | 3.40190 | 0.04701 | 0.31550 | 0.00000 |

| B--fit | $`A_0`$ | $`A_1`$ | $`A_2`$ | $`A_3`$ | $`A_4`$ | $`A_5`$ |
| --- | --- | --- | --- | --- | --- | --- |
| $`d_v`$ | 1.43611 | 0.61170 | 4.73270 | -0.67420 | 0.24920 | 3.00000 |
| $`u_v`$ | 1.71921 | 0.55260 | 2.90090 | -2.39220 | 1.60430 | 1.61490 |
| $`g`$ | 29.76781 | 0.50800 | 2.38230 | 4.35570 | -2.33860 | 3.00000 |
| $`\overline{d}/\overline{u}`$ |  | 10.19560 | 5.16810 | 14.85860 | 17.00000 | 8.69760 |
| $`\overline{u}+\overline{d}`$ | 0.06717 | -0.29650 | 7.71700 | -0.75670 | 4.52290 | 0.64380 |
| $`s+\overline{s}`$ | 0.04356 | 0.00210 | 7.33918 | -1.26970 | 5.13400 | 0.60083 |
| $`s\overline{s}`$ | 0.01781 | 0.22103 | -15.02691 | 0.22693 | -1.23666 | 0.00000 |