text:
# 1 Introduction

## 1 Introduction

In an earlier paper  we presented the phenomenological evidence for the depletion of high-momentum gluons as a projectile proton traverses a target nucleus in $`pA`$ collisions. Here we present a more complete discussion of the subject with more technical details and with extension to nucleus-nucleus collisions.

Our approach is unconventional in that we do not make the usual assumption that the parton distributions in a proton remain unaltered as the proton propagates through a nucleus, even if a hard subprocess occurs deep in the nucleus. That assumption amounts to factorization, a property that has been proven for $`pp`$ collision, but clearly must fail for $`pA`$ collision if $`A`$ is infinitely large. For realistic nuclear sizes our alternative assumption that the parton distribution can be modified led to $`J/\psi `$ suppression  that cannot be distinguished from the effects of the usual mechanisms . What is different now is the appearance of new data on $`J/\psi `$ suppression in $`pA`$ collisions at large $`x_F`$ . Since those data cannot be explained in terms of hadronic absorption of the produced $`c\overline{c}`$ state , we find more direct support for the idea of gluon depletion before the hard subprocess that produces the $`c\overline{c}`$ state .

The data of the Fermilab E$`866`$ experiment  on the $`J/\psi `$ suppression in $`pA`$ collisions at $`800`$ GeVc are expressed as the ratio

$`R(x_F,A)=\sigma _A\left(x_F\right)/A\sigma _N\left(x_F\right)=A^{\alpha \left(x_F\right)1},`$ (1)
where $`\sigma _{N,A}`$ is the cross section for $`J/\psi `$ production by a proton on a nucleon $`(N)`$ or on a nucleus $`(A)`$. An analytical formula for $`\alpha (x_F)`$ is given in

$`\alpha \left(x_F\right)=0.96\left(10.0519x_F0.338x_F^2\right),`$ (2)
for $`0.1<x_F<0.9`$. Equation (2) differs significantly from the one given in their preprint , on the basis of which the analysis in  was done. Here we perform a reanalysis based on the new parameterization in Eq. (2).

In our view there are three sources that can contribute to the $`x_F`$ dependence, and we write them in a product form

$`R(x_F,A)=G(x_F,A)N(x_F,A)H(x_F,A),`$ (3)
where $`H(x_F,A)`$ represents the hadronic absorption of the $`c\overline{c}`$ state, $`N(x_F,A)`$ nuclear shadowing, and $`G(x_F,A)`$ the gluon depletion effect, defined as

$`G(x_F\left(x_1\right),A)=g(x_1,A)/g(x_1,0),`$ (4)
$`x_1`$ being the momentum fraction of the gluon in the projectile proton, and $`g(x_1,A)`$ being the effective gluon distributions of that projectile in a nucleus A at the point of $`c\overline{c}`$ production, averaged over the penetration depth at which the hard process occurs.

By $`H(x_F,A)`$ we mean the absorption factor that operates between the production of the $`c\overline{c}`$ state and the detection of the final $`J/\psi `$ due to any mechanism, including the interaction with comovers. Since no good arguments have been advanced to show that $`H(x_F,A)`$ can have a significant dependence on $`x_F`$ , we shall assume in the following that $`H(x_F,A)`$ is independent of $`x_F`$. This is not a serious limitation in our formalism. If a reliable $`x_F`$ dependence is found at a later date, we can easily incorporate it in our analysis. For now we adopt the usual Gerschal-Hüfner form.

$`H(A)=\mathrm{exp}\left[\sigma \rho L(A)\right]`$ (5)
where $`\sigma `$ is the absorption cross section, $`\rho `$ the nuclear density, and $`L(A)`$ the mean path length in $`A`$ that a $`c\overline{c}`$ state propagates.

In the next section the nuclear shadowing factor $`N(x_F,A)`$ will be discussed. We shall find a simple formula that can represent the change in the gluon distribution in the target nucleus due to shadowing and anti-shadowing. Such a formula facilitates the analysis, and offers a simple parametrization that is convenient to use, independent of the particular problem that we apply it to here.

The determination of the depletion factor $`G(x_F(x_1),A)`$ is the main theme of this paper. In Sec. 3 we shall go beyond a review of the content of Ref.  , not only because the data has changed from those in the preprint , resulting in numerical differences, but also because we shall improve on  in some technical details and include some new material.

In Sec. $`4`$ we shall extend the result from the study of the $`pA`$ problem to $`AB`$ collisions. We shall show how the enhanced depletion at large $`x_F`$ does not affect the $`J/\psi `$ suppression at small $`x_F`$, which is where the existing data for $`AB`$ collisions were collected. We shall have nothing to add on the subject of enhanced suppression observed in $`Pb`$-$`Pb`$ collisions beyond what we have advanced in Ref. .

## 2 Nuclear Shadowing and Anti-shadowing

The nuclear shadowing and anti-shadowing problem has been studied phenomenologically by Eskola et. al. . Instead of focusing on the physics of the origin of the problem in QCD, they analyzed the deep inelastic scattering data of nuclear targets. On the basis of DGLAP evolution  they can determine the parton distribution at any $`Q^2>2.25`$ GeV2. The results are given in terms of numerical parametrizations (called EKS98 ) of the ratio

$`N_i^A(x,Q^2)=f_{i/A}(x,Q^2)/f_i(x,Q^2),`$ (6)
where $`f_i`$ is the parton distribution of flavor $`i`$ for the free proton and $`f_{i/A}`$ is that for a proton in a nucleus $`A`$.

For the purpose of $`J/\psi `$ production we are interested in Eq. (6) for only the gluons. We consider only the dominant subprocess $`g_1(x_1)+g_2(x_2)c+\overline{c}`$, where $`x_1`$ is the momentum fraction of the projectile gluon whose distribution is $`g_1(x_1)`$, and $`x_2`$ is that of the target gluon whose distribution is $`g_2(x_2)`$. For the $`c\overline{c}`$ state produced with momentum fraction $`x_F=x_1x_2`$, the usual kinematical relations are

$`x_{1,2}=[(4\tau +x_F)^{1/2}\pm x_F]/2,x_1x_2=\tau M_{J/\psi }^2/s,`$ (7)
where it is assumed that the $`c\overline{c}`$ state that turns to $`J/\psi `$ is produced near threshold. Thus the virtuality of the subprocess $`g+gc+\overline{c}`$ is given by the value of $`M_{J/\psi }^2`$, or slightly higher. We shall take $`Q^2=10`$ GeV2, a value that is chosen in EKS98  to give explicit values of $`N_{i}^{}{}_{}{}^{A}(x_1,Q^2)`$. For simplicity we shall label the gluon distribution at $`Q^2=10`$ GeV2 by $`N(x_2,A)`$.

The numerical output of EKS98 for $`N(x_2,A)`$ is shown by the points in Fig. 1 for $`A=50`$, $`100`$, and $`200`$. We exhibit only the values for $`x_2`$ in the range $`0.01x_20.12`$, since that is the range relevant for the production of $`J/\psi `$ at $`800`$ GeV/c for $`0<x_F<0.8`$. For the purpose of convenience in using those values of $`N(x_2,A)`$ in analytic manipulation and computation, we propose a simple formula that contains the shadowing and anti-shadowing effects. Since the cross-over of the two effects occurs at $`x_2=0.02`$ where $`N(x_2,A)=1`$ for all $`A`$, it is sensible to use an auxiliary variable $`\xi `$, defined by

$`\xi =3.912+\mathrm{ln}x_2`$ (8)
which vanishes at $`x_2=0.02`$. Moreover, we notice that $`\mathrm{ln}N(x_2,A)`$ depends linearly on $`\mathrm{ln}A`$ to a good approximation, so it suggests a power-law behavior

$`N(x_2,A)=A^{\beta (x_2)}.`$ (9)
The exponent $`\beta (x_2)`$ can be determined by fitting the data for $`A=100`$. In terms of $`\xi `$ we find a good fit with the parametrization

$`\beta (\xi )=\xi (0.0284+0.0008\xi 0.0041\xi ^2),`$ (10)
the result of which is shown in Fig. 2, where the points are for $`\mathrm{ln}N(x_2,A)/\mathrm{ln}A`$ obtained from EKS98 at $`A=100`$. This convenient formula for $`\beta (\xi )`$ can then be used in conjunction with Eq. (8) to determine $`N(x_2,A)`$ for other values of $`A`$. The curves in Fig. 1 exhibit the good agreement between the data and our parameterization for $`A=50`$, $`100`$, and $`200`$.

In the following we shall simply use Eq. (9) as a summary of the effects of nuclear shadowing and anti-shadowing for problems in $`pA`$ collisions where the relevant range of $`x_2`$ is in the interval $`0.01x_20.12`$.

## 3 Evolution of Gluon Distribution in a Nucleus

We now may regard $`R(x_F,A)`$ and $`N(x_2,A)`$ as known phenomenologically. Thus, from Eq. (3) we may write

$`G(x_F,A)H(A)=A^{\alpha (x_F)\beta (x_2(x_F))1}.`$ (11)
Although the form of $`H(A)`$ is given by Eq. (5), we do not know the value of $`\sigma `$ in the present circumstance where we allow the possibility of gluon depletion. Hence, we treat $`H(A)`$ temporarily as unknown, along with the depletion factor $`G(x_F,A)`$. However, the $`x_F`$ dependence is completely known from the RHS of Eq. (11).

To proceed we need some theoretical input on the possible form of $`G(x_F,A)`$, or more directly the gluon distribution $`g(x_1,z)`$ in the projectile, where $`z`$ is the distance traversed in a nucleus. Note that we have refrained from referring to the projectile as the proton, since the possible modification of $`g(x_1,0)`$ for $`z>0`$ implies that the incident proton loses its usual identity, in particular, the nature of its partonic content, as the projectile, now identified only as a flux of partons, propagates in the nuclear medium. How $`g(x_1,z)`$ evolves in the nuclear medium is clearly a nonperturbative process that involves multiple scatterings of gluons and quarks at low virtualities. Nevertheless, for every incremental distance, $`dz`$, that a gluon travels the modification that $`g(x_1,z)`$ undergoes must be perturbative in that $`g(x_1,z+dz)g(x_1,z)`$ is small and is proportional to $`dz`$. It is therefore reasonable to adopt an evolution equation similar in spirit to that of DGLAP , but with the change in resolution scale $`d\mathrm{ln}Q^2`$ replaced by the change in penetration depth $`dz`$, so that we write

$`{\displaystyle \frac{d}{dz}}g(x,z)={\displaystyle _0^1}{\displaystyle \frac{dx^{}}{x^{}}}g(x^{},z)Q({\displaystyle \frac{x}{x^{}}}),`$ (12)
where the unknown kernel $`Q(x/x^{})`$ controls the gain and loss of the gluons in $`dz`$. $`Q(y)`$ cannot be determined by perturbative calculation, as the splitting functions in pQCD for $`Q^2`$ evolution. Equation (12) is similar to the nucleonic evolution equation proposed in , except that this is now at the parton level. In Eq. (12) the quark sector has been left out for simplicity. To be more complete one should include also the effects of the couplings of gluons with the quarks, a task that is deferred to the future. Thus what we can achieve now is the determination of an effective kernel $`Q(y)`$ that can account for the enhanced suppression of $`J/\psi `$ at large $`x_F`$.

To determine the $`z`$ dependence of $`g(x,z)`$ let us take the moments by defining

$`g_n(z)={\displaystyle _0^1}𝑑xx^{n2}g(x,z)`$ (13)
and

$`Q_n={\displaystyle _0^1}𝑑xx^{n2}Q(y).`$ (14)
Then by the convolution theorem, Eq. (12) becomes

$`dg_n(z)/dz=g_n(z)Q_n,`$ (15)
whose solution is

$`g_n(z)=g_n(0)e^{zQ_n}.`$ (16)
It is possible that Eq. (15) is valid only when $`z`$ is large enough, in which case Eq. (16) should be modified to read

$`g_n(z)=g_n(z_0)e^{(zz_0)Q_n}`$ (17)
for $`z`$ greater than some positive value of $`z_0`$.

To proceed, let us substitute Eq. (4) in (11) and define

$`J(x_1,A)g(x_1,0)A^{\alpha (x_F(x_1))\beta (x_2(x_1))1}`$ (18)
where the interrelationships among $`x_1`$,$`x_2`$ and $`x_F`$ are specified by Eq. (7). For $`g(x_1,0)`$ we use the canonical form

$`g(x_1,0)=g_0(1x_1)^5.`$ (19)
The final result is insensitive to its form and independent of $`g_0`$, which we shall set to be $`1`$. Thus we may regard $`J(x_1,A)`$ as known. Since we also have

$`J(x_1,A)=g(x_1,A)H(A),`$ (20)
its moments are, by virtue of Eqs. (5) and (16),

$`J_n(z)=g_n(0)\mathrm{exp}[z(Q_n\sigma \rho )].`$ (21)
Here and in the following we shall use $`z`$ (until Sec. 5) to denote the average penetration depth (i.e., $`z\overline{z}_A`$) in $`A`$ when a $`c\overline{c}`$ state is produced by $`gg`$ annihilation. It is then also the average length that the $`c\overline{c}`$ must travel in $`A`$ and be subject to hadronic absorption, i.e., $`z=L(A)`$.

We can determine $`J_n(z)`$ by taking the moments of $`J(x_1,A)`$, as expressed in Eq. (18). However, there is a problem in evaluating $`J_n(z)=_0^1𝑑x_1x_{1}^{}{}_{}{}^{n2}J(x_1,z)`$, since $`J(x_1,A)`$ is ill-defined at $`x_1=0`$. According to Eq. (7), $`x_2`$ diverges as $`x_10`$; thus limiting $`x_2`$ to $`1`$ implies that $`x_1`$ cannot be less than $`M_{J/\psi }^2/s`$, a small but nonvanishing value. Furthermore, $`x_F`$ becomes negative at small $`x_1`$, and we lose any knowledge about $`\alpha (x_F)`$ for $`x_F<0.1`$ . Also, $`\beta (x_2)`$ is not reliably known at large $`x_2`$, so the RHS of Eq. (18) cannot offer accurate determination of $`J(x_1,A)`$ as $`x_10`$. These defects can be suppressed by the factor $`x_{1}^{}{}_{}{}^{n2}`$ in the integrand, if we restrict $`n`$ to $`3`$. We shall therefore determine $`J_n`$ only for $`n3`$.

For convenience, let us define

$`K_n(z)\mathrm{ln}[J_n(z)/g_n(0)]=z(Q_n\sigma \rho ).`$ (22)
Using Eqs. (2), (10) and (19) in (18), we can calculate $`J_n(x)`$ and therefore $`K_n(z)`$ for $`n3`$. As mentioned earlier, Eq. (2) is different from a previous form of $`\alpha (x_F)`$ given in  and used in . The results are shown as discrete points in Fig. 3 for $`A=100`$ and $`200`$. The corresponding values of $`z`$ are halves of the average total path lengths of the nuclei, i.e., $`z=3R_A/4=0.9A^{1/3}`$ fm and are therefore $`z_1=4.177`$ and $`z_2=5.262`$ fm, respectively.

To extract the information contained in those points in Fig. 3 we need an analytical representation of $`K_n`$. We choose to fit $`K_n`$ by the following formula, different from the one used in ,

$`K_n={\displaystyle \underset{i=0}{\overset{3}{}}}k_in^i+k_4n^{1/2}.`$ (23)
The results of our fits are shown by the smooth curves in Fig. 3. The corresponding parameters are given in Table 1. Note that the fits allow us to extrapolate smoothly to $`n=2`$, where we could not calculate $`J_2`$.

Before we determine $`Q_n`$ from $`K_n`$ in Eq. (22), we need to verify the $`z`$ dependence. It should first be recognized that the experimental parameterization of the $`A`$ dependence, such as in Eq. (1), is not compatible with the theoretical expectation, such as in Eq. (5) and (21), except in a certain range of $`A`$. Since a power-law $`A^\gamma `$, expressed as $`e^{\gamma \mathrm{ln}A}`$, is approximately $`e^{\gamma ^{}z}`$, where $`zA^{1/3}`$, only when $`\mathrm{ln}A`$ is approximately $`A^{1/3}`$, the correspondence can only be for $`60<A<240`$. With that understanding, let us nevertheless calculate $`J(x_1,z)`$ for all $`z<6`$ using $`A=(z/0.9)^3`$ in Eq. (18), take the moments, and then determine $`K_n(z)`$ through the first half of Eq. (22). The result is shown as points in Fig. 4 for eleven values of $`z`$ between $`0.9`$ and $`5.9`$, corresponding to $`A`$ from 1 to 282, and for four representative values of $`n`$, viz., 3, 8, 13 and 20. The straight lines are linear fits of the last six points for each value of $`n`$. Evidently, the $`z`$ dependence of $`K_n(z)`$ is very nearly linear for $`3.4<z<5.6`$, which corresponds to $`54<A<240`$. Thus our theoretical formalism is consistent with the experimental data in the region where $`lnAA^{1/3}`$. At $`z=0.9`$, or $`A=1`$, all points converge to $`K_n=0`$, as they should. We cannot reliably apply our formalism to the collision problems where $`A<50`$. Fig. 4 also suggests that even when $`A`$ is large, say $`>100`$, the gluon evolution equation (12) may not be valid at small $`z`$, here used in the sense of penetration depth within the large nucleus, not the average depth. In the following we shall limit our consideration to only the linear portion of Fig. 4. In that region the second half of Eq. (22) should be treated as differentially correct, i.e.,

$`\mathrm{\Delta }K_n(z)/\mathrm{\Delta }z=Q_n\rho \sigma `$ (24)
where $`\mathrm{\Delta }K_n(z)=K_n(z+\mathrm{\Delta }z)K_n(z).`$

Since the values of $`k_i`$ in Table 1 are determined in the linear region, we can use them to obtain $`Q_n`$. If we write

$`Q_n={\displaystyle \underset{i=0}{\overset{3}{}}}q_in^i+q_4n^{1/2},`$ (25)
then we have from Eqs. (23) and (24)

$`q_0=\mathrm{\Delta }k_0/\mathrm{\Delta }z+\rho \sigma ,q_i=\mathrm{\Delta }k_i/\mathrm{\Delta }z(1i4).`$ (26)
For the two $`z`$ values in Table 1, we get (with $`\mathrm{\Delta }z=1.086`$ fm)

$`\begin{array}{ccc}q_0=0.0758+\rho \sigma ,\hfill & q_1=0.0268,\hfill & \\ q_2=7.13\times 10^4,\hfill & q_3=9.6\times 10^6,\hfill & q_4=0.0568\hfill \end{array}`$ (29)
in units of fm-1.

Since the absorption cross section $`\sigma `$ is unknown when gluon depletion is not negligible, $`q_0`$ is not fixed by Eq. (29). Whatever the dynamics of gluon depletion is, we require that the total gluon momentum does not increase with $`z`$. Since the gluon momentum is $`𝑑xg(x,z)=g_2(z)`$, that requirement implies in conjunction with Eq. (17) that $`Q_20`$. Choosing the upper bound $`Q_2=0`$ leads to the condition, on account of Eqs. (25) and (29),

$`q_0=(2q_1+4q_2+8q_3+\sqrt{2}q_4)=0.0295.`$ (30)
With all parameters in Eq. (25) now determined, the $`n`$ dependence of $`Q_n`$ can be exhibited, as shown in Fig. 5. Evidently, $`Q_n`$ is smoothly varying and the use of the polynomials in Eqs. (23) and (25), which are different from those used in , is justified.

Although it is more direct to proceed immediately to the use of Eq. (17) to the determination of $`g_n(z)`$ and therefore $`G(x_1,z)`$, which is our goal, there is some advantage in an attempt to find $`Q(y)`$ at this point, while we are on the subject of $`Q_n`$. The easiest way to do that is to put $`Q_n`$ in form

$`Q_n=c_0+{\displaystyle \underset{j}{}}c_j/(n+j1),`$ (31)
used in , so that it can imply directly

$`Q(y)=c_0\delta (1y)+{\displaystyle \underset{j}{}}c_jy^i,`$ (32)
To translate from Eq. (25) and (31), we fit the values of $`Q_n`$ at the 19 integer points, $`2n20`$, determined by (25), by use of the formula (31) with a suitable number of terms in the sum. It turns out that a good fit can be achieved with three terms: $`j=3,4`$ and $`5`$. The result is

$`c_0=0.1988,c_3=6.205,c_4=23.316,c_5=19.866.`$ (33)
We show in Fig. 6 both the discrete values of $`Q_n`$ at integral $`n`$ and the fitted curve using Eqs. (31) and (33). The corresponding function $`Q(y)`$, calculated using Eqs. (32) and (33), is shown in Fig. 7. It is evident that the $`c_0`$ and $`c_4`$ terms correspond to gluon depletion, while the $`c_3`$ and $`c_5`$ terms correspond to gluon regeneration.

To determine $`G(x_1,z)`$ as defined in Eq. (4), we have Eq. (17) that specifies the evolution in $`z`$ in the linear region (see Fig. 4) from $`z_0`$. At this point we have no formalism to extrapolate in the nonlinear region from $`z_0`$ down to $`0`$, and it is in reference to $`g(x_1,0)`$ that $`G(x_1,z)`$ is defined. However, in view of the fact that the hadronic absorption term $`H(A)`$ is known empirically to be an exponential, as in Eq. (5), i.e., $`\mathrm{ln}H(z)`$ being linear in $`z`$ for all $`z`$, one may regard the nonlinear portion of Fig. 4 to be due primarily to the mismatch between $`\mathrm{ln}A`$ and $`A^{1/3}`$ at low $`A`$. Then we adopt the approximation that Eq. (16) is adequate for relating $`g_n(z)`$ to $`g_n(0)`$. With $`Q_n`$ being known from Eqs. (25), (29) and (30), and $`g_n(0)=B(n1,6)`$ which is the beta function, we can calculate $`g_n(z)`$. The result is shown in Fig. 8 by the full (open) circles for $`A=100(200)`$, It is natural to fit the resultant $`g_n(z)`$ by a linear combination of beta functions in the form

$`g_n(z)={\displaystyle \underset{i=1}{\overset{3}{}}}a_i(z)B(n1,5+i),`$ (34)
The coefficients $`a_i(z)`$ are determined by fitting $`g_n(z)/g_n(0)`$ in order to reduce the range of variation. For the two values $`z_1`$ and $`z_2`$, corresponding to $`A=100`$ and $`200`$, we obtain

$`\begin{array}{ccc}a_1(z_1)=0.3526,\hfill & a_2(z_1)=1.44,\hfill & a_3(z_1)=0.78,\hfill \\ a_1(z_2)=0.2362,\hfill & a_2(z_2)=1.655,\hfill & a_3(z_2)=0.869.\hfill \end{array}`$ (37)
The curves in Fig. 8 are generated using Eqs. (34) and (37). Evidently, the fits are good.

The inverse transform of the moments in Eq. (34) is

$`g(x_1,z)={\displaystyle \underset{i=1}{\overset{3}{}}}a_i(z)(1x_1)^{4+i},`$ (38)
whose implication for

$`G(x_1,z)`$ $`=`$ $`g(x_1,z)/g(x_1,0)`$ (39)
$`=`$ $`a_1(z)+a_2(z)(1x_1)+a_3(z)(1x_1)^2`$
can readily be calculated using Eq. (37). The results for $`z_1`$ and $`z_2`$ are shown in Fig. 9. Clearly, there is significant depletion of gluons at large $`x_1`$, roughly 40% at $`x_10.8`$. There is a small amount of regeneration at small $`x_1`$. The cross-over is at around $`x_10.2`$. Although the enhancement at small $`x_1`$, is at the 2 to 3% level, in terms of the number of gluons in a small $`dx`$ interval it is not insignificant compared to the depletion at large $`x_1`$, because $`g(x_1,0)`$ is strongly damped at large $`x_1`$. That is how the condition $`Q_2=0`$ is satisfied.

## 4 Suppression of $`J/\psi `$ and $`\psi ^{}`$

Having determined $`G(x_1,A)`$ in the previous section, we can now return to the problem of charmonium suppression, including $`\psi ^{}`$. The Fermilab E866 experiment  gives data for both $`J/\psi `$ and $`\psi ^{}`$ in the form of $`\alpha (x_F)`$. In using Eq. (3) to calculate $`R(x_F,A)`$, and then $`\alpha (x_F)`$, we have Eqs. (39), (9) and (5) for $`G(x_1,A)`$, $`N(x_2,A)`$ and $`H(A)`$, respectively. For both $`J/\psi `$ and $`\psi ^{}`$ we use $`\sigma =6.5`$ mb in $`H(A)`$. The results on $`\alpha (x_F)`$ are shown as curves in Fig. 10 for $`J/\psi `$ and Fig. 11 for $`\psi ^{}`$. The two curves in each figure are for $`A=100`$ and $`200`$; they are sufficiently close to each other to be almost independent of $`A`$, thereby affirming our gluon evolution model for that range of $`A`$.

The most significant part of what we have learned from this work on the $`pA`$ data is that it is hard to reproduce the strong damping of the measured $`\alpha (x_F)`$ at large $`x_F`$ without a substantial amount of gluon depletion at large $`x_1`$, as seen in Fig. 9. Furthermore, if there is significant modification of the parton distribution as the penetration depth increases, then it is hard to justify the notion that a proton can traverse a large part of the target nucleus without any changes, or that it can be wounded upon the first collision with a nucleon and then remaining unchanged thereafter.

So far out attention has been given to the $`x_F>0`$ region only where data exist. However, the $`x_F<0`$ region has interesting physics also, and the data can be obtained either at RHIC, or at a fixed target experiment with proton being the target for a heavy-ion beam. At large negative $`x_F,x_1`$ would be small where Fig. 9 shows a very small enhancement, not depletion. On the other hand, $`x_2`$ would be large, where $`N(x_2,A)`$ would exhibit not only anti-shadowing, but also the EMC and Fermi motion effects. As an illustration of the $`x_F`$ dependence, we have performed the calculation for $`J/\psi `$ suppression in $`pW`$ collision at three energies and for all $`x_F`$, negative as well as positive. The absorption cross section is set at $`\sigma =8`$ mb, a value found in the next section. The values of $`N(x_2,A)`$ for the full $`x_2`$ range are obtained from EKS98 . The results are shown in Fig. 12. We see that the negative $`x_F`$ region shows very little dependence on energy and reflects mainly the property of $`N(x_2,A)`$.

Gluon depletion of the type discussed in this paper has no effect on the suppression of $`J/\psi `$ in the $`x_F<0`$ region. However, there is another type of gluon depletion, discussed in the second paper in Ref. , called nonlinear depletion, that can influence the survival probability in the $`x_F<0`$ region. When $`x_2`$ is large, the fast gluons in the rear part of the nucleus can catch up and interact with the slow gluons released by the $`p`$-$`N`$ collisions in the front part of the nucleus, leading to a depletion of gluon at large $`x_2`$. Although this involves the interaction between gluons in different nucleons in the nucleus, it is not nuclear shadowing, since the conventional nuclear shadowing at low $`x_2`$ does not require the invasion of an external proton to initiate the cascading interactions of gluons among the broken nucleons. We have no prediction on the nature of the effect of this type of depletion. If the experimental data on the production rate differ significantly from the curves in Fig. 12 in the $`x_F<0`$ region, then there will be strong motivation to consider this unusual type of gluon depletion.

## 5 Nucleus-nucleus Collisions

Having investigated the $`pA`$ collision problem above, and finding the necessity to consider the depletion of gluons in the projectile before the production of $`c\overline{c}`$ states in the $`x_F>0`$ region, it is natural to ask next what the implication would be for the $`J/\psi `$ suppression problem in $`AB`$ collisions. Even without any detailed calculations it is straightforward to infer that there will be enhanced suppression at large $`x_F`$. However, the currently available data on $`J/\psi `$ suppression in nuclear collisions are from CERN-SPS, and are limited to the central rapidity region. It is therefore our burden to show that the gluon depletion mechanism discussed in Sec. 3 is consistent with the existing nuclear collision data at small $`x_F`$.

In Sec. 3 we have used $`z`$ mainly as the average penetration depth in the target nucleus where $`c\overline{c}`$ is produced, although in the evolution equation (12) $`z`$ is the actual path length. Since in this section we need to average the production point over all impact parameters in the $`AB`$ collisions, we now restore $`z`$ to be actual path length of a gluon in a nucleus. Thus we allow it to vary from 1 to 12 fm. We use Eq. (16) to determine $`g_n(z)`$ at all odd integer points in the range $`1z12`$ (for the sake of simplicity in fitting), and then fit them at each such value of $`z`$ by the formula (34). The resulting values of $`a_i(z)`$ are shown by the points in Fig. 13, they are in turn fitted by quadratic equations of the form

$`a_i(z)=b_{i0}+b_{i1}z+b_{i2}z^2.`$ (40)
The result is shown by the curves in Fig. 13. The corresponding coefficients are

$`\begin{array}{ccc}b_{10}=0.974,\hfill & b_{11}=0.175,\hfill & b_{12}=0.0068,\hfill \\ b_{20}=0.0866,\hfill & b_{21}=0.402,\hfill & b_{22}=0.0202,\hfill \\ b_{30}=0.066,\hfill & b_{31}=0.226,\hfill & b_{32}=0.0142.\hfill \end{array}`$ (44)
Using Eqs. (40) and (44) in (38), we can evaluate $`g(x_1,z)`$ at all $`x_1`$ and $`z`$.

The cross section for the production of $`J/\psi `$ in $`AB`$ collisions can be calculated in the standard way. We shall just write it down as follows (see, e.g., Refs ):

$`\sigma _{J/\psi }`$ $`=`$ $`{\displaystyle d^2bd^2s𝑑z_A𝑑z_B\rho _A(\stackrel{}{s},z_A)\rho _B(\stackrel{}{s}\stackrel{}{b},z_B)}`$ (45)
$`{\displaystyle }{\displaystyle \frac{dx_1}{x_1}}{\displaystyle \frac{dx_2}{x_2}}g_A(x_1,L_Bz_B)g_B(x_2,L_Az_A)N(x_1,A)N(x_2,B)`$
$`e^{\sigma _a\left[\rho _A(L_A+z_A)+\rho _B(L_B+z_B)\right]}\widehat{\sigma }_{ggc\overline{c}}(x_1,x_2),`$
where we have included the $`z_{A,B}`$ dependences in the gluon distributions $`g_{A,B}(x_{1,2},L_{B,A}z_{B,A})`$ and the nuclear shadowing functions $`N(x_1,A)`$ and $`N(x_2,B)`$. $`L_{A,B}`$ are the path length through $`A(B)`$ at the distances $`\stackrel{}{s}`$ and $`\stackrel{}{s}\stackrel{}{b}`$, respectively, from the centers of the nuclei, i.e.

$`L_A=(R_A^2s^2)^{1/2},L_B=(R_B^2|\stackrel{}{s}\stackrel{}{b}|^2)^{1/2}.`$ (46)
The gluon distribution $`g_A(x_1,L_Bz_B)`$ is given by \[see Eq. (38)\]

$`g_A(x_1,L_Bz_B)={\displaystyle \underset{i=1}{\overset{3}{}}}a_i(L_Bz_B)(1x_1)^{4+i},`$ (47)
and similarly, for $`g_B(x_2,z_A)`$, the coefficients $`a_i(z)`$ in Eq. (38) are replaced by $`L_Az_A`$, which is the distance that a parton in $`B`$ travels in $`A`$ before the production of $`c\overline{c}`$ at $`z_A`$. The distance that a $`c\overline{c}`$ state travels in $`A`$ is $`L_A+z_A`$, and $`z_A`$ is integrated from $`L_A`$ to $`+L_A`$. For the energies at CERN-SPS, $`\sqrt{s}1820`$ GeV, and for $`x_F0`$, the hard cross section $`\widehat{\sigma }_{ggc\overline{c}}(x_1,x_2)`$ restricts the gluon momentum fractions to $`x_1x_2M_{J/\psi }/\sqrt{s}0.16`$.

The survival probability is

$`S_{J/\psi }^{AB}`$ $`=`$ $`\sigma _{J/\psi }^{AB}/\sigma _{J/\psi }^{AB(0)}`$ (48)
$`=`$ $`N_{AB}^1{\displaystyle d^2bd^2s_{L_A}^{L_A}𝑑z_A_{L_B}^{L_B}𝑑z_BW(\stackrel{}{b},\stackrel{}{s},z_A,z_B)}`$
where

$`W(\stackrel{}{b},\stackrel{}{s},z_A,z_B)`$ $`=`$ $`G_A(x_1,L_Bz_B)G_B(x_2,L_Az_A)N(x_1,A)N(x_2,B)`$ (49)
$`e^{\sigma _a\left[\rho _A(L_A+z_A)+\rho _B(L_B+z_B)\right]},`$
and $`N_{AB}`$ is the same integral in Eq. (48) but with $`W(\stackrel{}{b},\stackrel{}{s},z_A,z_B)`$ replaced by 1. $`G_A(x_1,L_Bz_B)`$ is as given in Eq. (39) except that $`z`$ is replaced by $`L_Bz_B`$. Because of both the gluon enhancement at $`x_{1,2}0.16`$ and the anti-shadowing, the absorption cross section $`\sigma _a`$ now has to be somewhat larger than before . An overall agreement with all the $`pA`$ and $`AB`$ collision data, except the $`Pb`$-$`Pb`$ case, can be achieved with the use of one value of $`\sigma _a=8`$ mb. The result is given in Table 2 for the various $`AB`$ cases. Fig. 14 shows

how those values compare with the experimental data  by the straightline segments that connect those successive points. It is evident that apart from the $`Pb`$-$`Pb`$ case the agreement with the data is satisfactory. Thus we can conclude that the gluon depletion mechanism used to treat the $`pA`$ problem leads to no disagreement with the $`AB`$ collisions—except that the $`Pb`$-$`Pb`$ case remains as an anomaly.

Finally, we compute the $`x_F`$ dependence of the $`J/\psi `$ suppression factor for just one nuclear-collision case as an example, which we take to be $`Pb`$-$`Pb`$. We consider two energies: $`E_{\mathrm{lab}}=160`$ GeV and $`\sqrt{s}=60`$ GeV for RHIC. Except for the kinematics in $`\widehat{\sigma }_{ggc\overline{c}}`$ that affects the values of $`x_1`$ and $`x_2`$, the cross section and survival factor can be calculated as before. The results are shown in Fig. 15, which exhibits a substantial degree of suppression at large $`x_F`$. Any data at large $`x_F`$ would put considerable constraint on the models that attempt to explain the anomalous suppression at small $`x_F`$.

## 6 Conclusion

In our attempt to understand the enhanced suppression of $`J/\psi `$ at large $`x_F`$, we have found that the depletion of gluons at large $`x_1`$ in the projectile is the most natural explanation for the effect. We have proposed an evolution equation for the gluon distribution as the gluons propagate in a nuclear medium. The depletion at high $`x_1`$ contributes to a mild growth of the gluon distribution at small $`x_1`$. However, that growth does not lead to any contradiction with the existing data on $`J/\psi `$ suppression at mid-rapidity. Indeed, we have gone further to show where to find informative clues on the dynamics of suppression at large (positive and negative) $`x_F`$ in both $`pA`$ and $`AB`$ collisions.

What we have done here is only a modest first step towards understanding parton evolution in nuclear matter. While concentrating on the gluons, we have ignored the influence of the quark sector, a subject to be investigated at a later point. The depletion of quarks at large $`x`$ reveals itself in the suppression of dilepton and leading meson production at large $`x_F`$, the experimental evidences for which exist, though in subtle ways. Because of the conservation of Fermion number, the degradation of the quark distribution at large $`x`$ cannot be substantial. Nevertheless, the influence on the gluon distribution at small $`x`$ may not be negligible.

An important implication of this work is that in $`pA`$ or $`AB`$ collisions the concept of a nucleon propagating through nuclear matter as an identifiable, fixed entity needs modification. The usual notion that in nuclear collisions the total transverse energy $`E_T`$ is proportional to the number of nucleon-nucleon collisions would seem to have difficulty in reconciling with the insistence that the nucleons remain unaltered, if each inelastic collision of the nucleons contributes a fraction of their energies to $`E_T`$. The wounded nucleon model  makes a crude approximation of what goes downstream as an average quantity that is different from the incident nucleon, but ignores the way it changes as it propagates. Our evolution equation indicates that the parton flux changes continuously and may emerge with a profile that cannot be identified with that of a free nucleon in any sensible comparison. The revelation made by this understanding will undoubtedly affect many aspects of high-energy nuclear collisions.

Acknowledgment

We are grateful to Kari Eskola for providing us with the Fortran codes for EKS98. This work was supported, in part, by the U.S.-Slovakia Science and Technology Program, the U. S. National Science Foundation under Grant No. INT-9319091 and by the U. S. Department of Energy under Grant No. DE-FG03-96ER40972.

## Figure Captions

Nuclear shadowing and anti-shadowing factor taken from EKS98  for $`Q^2=10`$ GeV2. The curves are fits by the simple formula in Eq. (9).

The $`\beta (\xi )`$ function used to fit EKS98 data points at $`A=100`$.

$`K_n`$ calculated from the moments at discrete $`n`$ and fitted by the formula in Eq. (23).

The $`z`$ dependences of $`K_n(z)`$. The lines are straight-line fits in the large $`z`$ region.

$`Q_n`$ as calculated from Eq. (25).

The points of $`Q_n`$ are determined from Eq. (25) or Fig. 5 at integer $`n`$ values; the curve is a fit using Eq. (31).

The kernal $`Q(y)`$ as calculated from Eq. (32).

$`g_n`$ for two values of $`A`$.

The ratio of gluon distributions, $`G(x_1,A)`$, for two values of $`A`$.

$`\alpha (x_F)`$ for $`J/\psi `$ production. The data points are from Ref. ; the curves represent our result.

$`\alpha (x_F)`$ for $`\psi ^{}`$ production. The data points are from Ref. ; the curves represent our result.

The survival probability for $`J/\psi `$ production in $`pW`$ collisions at different energies for the entire range of $`x_F`$.

The $`z`$ dependences of the coefficients $`a_i(z)`$.

The survival probability in $`AB`$ collisions. The data points are from Ref. ; the line is composed of straightline sections connecting the calculated points listed in Table 2.

The survival probability for $`J/\psi `$ in $`Pb`$-$`Pb`$ collision for all $`x_F`$ at two energies.


format: text

text:
# ESO Imaging Survey

## 1 Introduction

The successful application of photometric techniques to identify high-redshift galaxies ($`z>3`$) combined with large-aperture telescopes with multi-object spectrographs have provided the means to study large-scale structures at high-redshifts and to directly probe their evolution over a broad range of redshifts. Currently, the challenge is to construct well-defined statistical galaxy samples at faint magnitudes and use photometric redshift techniques to pre-select candidates for large spectroscopic surveys probing the high-redshift Universe. To achieve this goal, deep optical/infrared multi-passband surveys covering areas of the order of a few degrees are required.

Interest in this type of dataset is widespread in Europe given the imminent commissioning of the multi-object spectrograph VIMOS, part of the new generation of VLT instruments. Foreseeing this need the Working Group (WG) for public surveys at ESO recommended in early 1999 the undertaking of a deep, multi-passband optical survey, covering three separate regions of sky one square degree each, with some areas also covered in the infrared. These regions, referred to as Deep-1, 2 and 3, are located over a wide range of right ascension and each correspond to four adjacent pointings (a-d) of the Wide-Field Imager (WFI) mounted on the 2.2m MPG/ESO telescope at La Silla. The main goal of the survey is to produce suitable statistical samples to study large-scale structures in the redshift domains $`z<1`$ and $`z>3`$. Details of the position of the selected regions can be found at the URL “http://www.eso.org/science/eis/”. One of the selected regions includes the Chandra Deep Field South (CDF-S; Giacconi et al. 2001) chosen because of its low column-density of hydrogen and absence of bright stars, characteristics that make it well-suited for deep X-ray observations with the Chandra (AXAF) telescope. These observations have recently been completed for a total integration time of one million seconds producing the deepest, high-resolution X-ray image ever taken (Rosati et al. 2001).

In this paper the $`UBVRI`$ optical observations with the WFI of one of the selected fields (Deep-2c), covering the CDF-S, are presented. These data are used to make a first assessment of the procedures and overall performance of the optical part of the Deep Public Survey (DPS). The optical data presented here complement those in the near-infrared covering an area of $`0.11`$ square degrees (Vandame et al. 2001) encompassing the region covered by the Chandra observations. The new data reported here extend earlier optical/infrared observations of this field carried out by EIS at the NTT in 1998, covering a considerably smaller region. These older data are currently being re-analyzed to produce, together with the new data, a homogeneous dataset (Benoist et al. 2001a) which will supersede the earlier 1998 release (Rengelink et al. 1998). Section 2 reviews the overall observing strategy and describes the observations. Section 3 discusses the techniques employed in the data reduction and the photometric and astrometric calibration of the images. While other derived products are expected to be made available in the near future, only photometric parameters of the sources detected in each of the available passbands are presented for the moment in Section 4. Even though the goal of the present paper is not to interpret the data being released, Section 5 presents the results of comparisons between the present observations and those of other authors. This is done for the sole purpose of assessing the quality of the astrometry and photometry of the present data set and its adequacy to meet the main goals of the survey. Finally, a brief summary is presented in Section 6.

## 2 Observations

The optical observations of the CDF-S region ($`\alpha 3^h32^m`$, $`\delta 27^o48^{}`$) were carried out using the WFI camera at the Cassegrain focus of the MPG/ES0 2.2m telescope at the La Silla observatory. WFI is a focal reducer-type mosaic camera with $`4\times 2`$ CCD chips. Each chip is composed of a $`2048\times 4098`$ pixel array and covers a field of view of $`8.12\times 16.25`$ arcmin with a projected pixel size of 0.238 arcsec. The 8 CCDs are physically separated by gaps of width 23.8 and 14.3 arcsec along the right ascension and declination directions, respectively. The full field of view of the camera is thus $`34\times 33`$ arcmin, with a filling factor of 95.9%. Details on the WFI layout can be found, for instance, at the URL “http://www.eso.org/science/eis/eis\_proj/public/wfi\_chip.html” as well as at other locations on the ESO web pages.

Observations of the CDF-S have been performed using the WFI filters summarized in Table LABEL:tab:filters. The $`U`$band observations were originally conducted using the available narrow filter which was recently (October 2000) replaced by a broader filter to increase the efficiency of the survey. The transmission curves of these filters are plotted in Figure 1 together with their final throughputs taking into account the optics and the CCD quantum efficiency. It is worth mentioning that very little was known about the performance of the instrument at the time the survey was conceived, especially at the blue and red ends of the wavelength range considered.

The observations in each passband were split into observation blocks (OBs) consisting of a sequence of five exposures with the integration time for each individual exposure ranging from 5 to 15 minutes depending on the passband. Table LABEL:tab:strategy summarizes the planned observational strategy listing: in column (1) the passband; in column (2) the integration time of each individual exposure; in column (3) the number of exposures per OB; in column (4) the integration time per OB; in column (5) the number of OBs; and in column (6) the planned total integration time per band. According to the planned strategy a total of 210 WFI frames were required to complete a single pointing. In order to minimize the imprint of the inter-chip gaps the choice of pointings were done in the following way. First, the nominal position of each OB (in a given passband) was drawn from a pre-defined list of distinct positions, as given in Table LABEL:tab:mainpos. For instance, in the case of $`U`$band which requires 13 OBs, the dithering pattern used is illustrated in Figure 2. Next, each exposure in the OB sequence was dithered relative to the OB reference position using the offsets in right ascension ($`\mathrm{\Delta }\alpha \mathrm{cos}\delta `$) and in declination ($`\mathrm{\Delta }\delta `$) as listed in Table LABEL:tab:dither and illustrated in the inset of Figure 2. This dithering is essential for de-fringing as well as for efficiently removing cosmic rays and other CCD blemishes.

The optical observations span six observing runs over a period of one year (ESO observing periods 64 and 66), distributed over a total of 12 nights. The log of the observations is summarized in Table LABEL:tab:logs. The table lists: in column (1) the date of the observations; in column (2) the passband; in column (3) the integration time; and in columns (4) and (5) the range and average seeing during the night as measured by the seeing monitor. Not listed are the $`B`$-band observations taken during the commissioning phase of the WFI, in January 1999, and the $`U`$-band data of the COMBO survey (Wolf et al. 2000). Both data sets have been re-processed, for the sake of homogeneity, and used to produce the final combined images. Table LABEL:tab:obs summarizes the data available for the CDF-S used in this paper. The table lists: in column (1) the band; in column (2) the total integration time; in column (3) the total number of frames; and in column (4) the average seeing of the final co-added images.

To assess the effectiveness of the adopted dithering pattern Figure 3 illustrates the exposure time map associated with the final $`R`$-band image which is an example of one of the passbands with the shortest exposure and is thus the least homogeneous (see Section 5). In this particular case the flux ratio of the best sampled region to that at the center of the image is a factor $`2`$. This is due in part to the clipping ($``$ 200 pixels) of the CCD frames not taken into account in the original design of the dithering pattern. In the weight-map the effects of at least three satellite tracks can be seen. It is important to point out that some of the survey data have been obtained using dithering patterns different from the one shown in Figure 3 (e.g., the $`B`$-images taken during the commissioning of WFI). In fact, during the course of the survey the pattern has been changed at least once to improve the uniformity of the final image.

The present optical survey complements both recent infrared observations reported by Vandame et al. (2001) as well as optical/infrared observations conducted as part of the EIS-DEEP survey first released in 1998 and currently under revision (Rengelink et al. 1998). The coverage of this field by these different sets of data is shown schematically in Figure 4.

## 3 Data Reduction

### 3.1 Pipeline Processing

The WFI images were processed using the refurbished EIS pipeline (e.g., da Costa et al. 1999; da Costa 2001) which includes new operational modules to enable unsupervised pipeline reduction of WFI images as well as new software routines for carrying out the basic steps of image processing and calibration. Given the scope of the surveys being carried out by EIS, special attention has been given to the development of a survey software system to ensure the adequate and timely processing of large amounts of data for long periods of time in a uniform way. For instance, the complete DPS survey will consist of some 2800 WFI science frames representing over 1 Terabyte of data, to which should be added a comparable amount of data from the so-called Pre-Flames survey not discussed here. While a detailed account of the extent of these changes is beyond the scope of the present paper, some of the main features are briefly summarized below.

One of the main changes has been the complete removal of different IRAF packages (e.g., xccdred, drizzle) used earlier to process images from CCD mosaics and in the co-addition phase. Instead, WFI images are now split into individual chips at the very beginning of the process. The main advantages of adopting this procedure are that the pipeline uses a common library of routines, single and multi-chip images are treated in the same way and that the EIS pipeline is more easily adaptable to use PC-based Linux clusters for parallel-processing images which will eventually lead to a higher throughput.

A major effort has also been made to address specific problems encountered with the data files coming from WFI. These have constantly evolved which makes them difficult to handle from the point of view of a pipeline. One of the main challenges has been the lack of information in the FITS header regarding the reference pixel. This prevented the unsupervised reduction of the WFI data due to frequent failures in the astrometric calibration from run to run. As the EIS observations with WFI span a long period of time (over two years), during which several upgrades of the data acquisition were made, it was necessary to develop robust procedures to cope with these various changes automatically.

Using tools available from the multi-resolution visual model package (MVM) developed by Bijoui and collaborators, a special routine was developed to find a first-order astrometric solution for an image, searching a large ($`6\times `$) area surrounding the nominal position of the pointing. The method uses one of the science frames in the run and takes the available information in the image header regarding the pointing only as a first guess. This image is quickly reduced and decomposed using a wavelet transform. The lowest resolution component is then used to cross-correlate features against a low-resolution mock image representing the reference catalog over a much larger area than the original image. The low-resolution science image is then moved around in a spiral pattern starting from its nominal position as given in the image header. At each step a cross-correlation function between the features in the real and mock images is computed and the process resumed. At the end of the search the position with the largest amplitude of the cross-correlation is used to define a first approximation for the location of the reference pixel for that run. This procedure is aimed at correcting only for large shifts. A more refined astrometric calibration is performed at a later stage.

Another important generic routine was developed to recognize different type of exposures (calibration frames such as dome, sky flats and bias, photometric standard stars or science frames), reject bad frames and associate frames into pseudo-OBs to be reduced according to pre-defined prescriptions depending on their nature and passband. The association is based on the frame type and on the spatial separation and time interval between consecutive frames. Next, the removal of instrumental signatures such as bias subtraction, flat-fielding and de-fringing are carried out using identical procedures as those applied for infrared images as presented by Vandame et al. (2001). Some of the main features of the software worth mentioning are that the stacking of images is carried out in pixel space (see below), thereby allowing the optimal rejection of cosmic rays and bad pixels using sigma–clipping. Satellite tracks, quite frequent in wide-field images, are efficiently masked on the individual exposures using the Hough transform (Vandame 2001), thereby minimizing their impact on the final stacked image. The number of satellite tracks removed in a single passband ranges from 9 (e.g.$`R`$-band images) to 75 ($`U`$-band images) depending on the exposure time. In principle, other linear features such as diffraction spikes may also be treated in the same way. The net result of these new developments has been the production of remarkably clean final images, significantly superior to those of earlier releases. These results lead to the production of considerably more reliable source detections since the successful masking of features at the level of individual exposures dramatically minimizes the number of false detections in the final image. Finally, using the sky-flats a first evaluation of the relative gains of the different chips is conducted for each passband which it is later applied to the science frames. These relative gains vary from passband to passband and are, in some cases, especially in the UV, as large as 10%.

The stacking of the images is done in two distinct steps. First, the five images belonging to the same OB are warped to the first image of the OB, which is used as reference for the pixel based relative astrometry. All the images in the OB are then averaged. The second step involves warping all the co-added images corresponding to each OB to an absolute astrometric reference catalog (see next section). The resulting images are then co-added taking into account their respective weight-maps, noise and flux scales. For this particular data set the warping of the images has been performed using the nearest-neighbor approach. The advantage of this warping is that the noise of the final image is uncorrelated. However, the astrometric solution is less precise due to discreteness effects. A new version of the warping has been developed allowing for a suite of kernel functions.

In contrast to other instruments (e.g., FORS) the fringing observed in $`I`$-band with WFI has a peak-to-peak amplitude of $`10\%`$ relative to the background and varies on scales of hours. The de-fringing is performed in exactly the same way as for the infrared images, except that the number of frames in an OB is considerably smaller, making a good estimate harder. Currently the residual contribution after de-fringing is comparable to the image noise. However, as in the case of the infrared, it would be preferable to decrease the integration time of each exposure and increase the number of dithered images. Based on the current findings future $`I`$band observations will be done using a different observing strategy.

The implementations described above make the EIS pipeline a tool capable of handling data from different telescope/instrument setups with different implementations of the FITS header. Finally, it is worth pointing out that using the available software several new developments are planned for the future. These include seeing deconvolution, morphological classification and difference image photometry for variability studies.

### 3.2 Astrometric calibration

The astrometric calibration of the images was performed using as the reference catalog a pre-release version of the Guide Star Catalog-II (GSC-II). The GSC-II (McLean et al. 2001) is based on multi-passbands all-sky photographic plate surveys including the Palomar Observatory Sky Survey (POSS-I-II), the SERC and the ESO Red Survey. The astrometric calibration of GSC-II has been obtained relative to Hipparcos (Perryman et al. 1997) and Tycho catalogs (Høg et al. 1997) and the ACT (Urban, Corbin & Wycoff 1998). The astrometric calibration was done using the method developed by Djamdji et al. (1993) based on a multi-resolution decomposition of images using wavelet transforms (MVM). An implementation of these algorithms to stack images, referred to as MVM-astrometry, has been done for EIS as is described above. This package has proven to be efficient and robust for pipeline reductions.

### 3.3 Photometric Calibration

Magnitudes were calibrated to the Johnson-Cousins system using Landolt standards taken from Landolt (1992). While standard stars were observed on every night, to optimize the observing time most observations were used to monitor the zero-point of the night and to compute color terms. For the final calibration, three nights in October 2000, all photometric, were used to obtain measurements of standards over a broad range of airmasses. During these nights a set of dithered images in each passband was taken overlapping two adjacent pointings (Deep2-b and Deep-2c). These were then used to set the absolute flux scale for all passbands, except the $`I`$band due to problems with de-fringing. For this case the night of Nov 04 1999 was used to determine the $`I`$-band calibration.

Photometric solutions were obtained using the new EIS photometric pipeline which is fully integrated to the associated database. The tool provides a convenient environment to obtain, verify and store the results of the fits. In the interactive mode it is possible to fit the data after rejecting individual measurements, stars and chips. Alternatively, one can also easily obtain independent solutions for each chip. Unfortunately, the number of standards for the fields considered does not allow independent solutions for each of the eight chips. Therefore, the present calibration had to rely on solutions including all the available chips.

Magnitudes for Landolt stars were obtained using an aperture 6 arcsec in diameter, which proved to be adequate, by monitoring the growth curve of all the measured stars. The accuracy of the zero-points are estimated to be in the range $`\pm `$0.03 mag to $`\pm `$0.08 mag, with the $`U`$-bands having the largest uncertainties. At face value these results indicate that the gain correction applied to the different chips based on the flatfield exposures is adequate for most passbands with the possible exception of the $`U`$-band. A final determination of the relative gains will require a more careful monitoring of the amplitude of these corrections and chip-based solutions. A set of secondary stars for all the Landolt fields considered would also be extremely useful to always allow a chip-based calibration. Efforts in this direction are already underway by several groups.

The computed extinction coefficients compare well with those estimated from the mean extinction curve of the Chilean sites. Finally, the following color terms, relative to the Johnson-Cousins (JC) system, have been empirically determined:

$`\begin{array}{c}(U_{JC}U_{EIS}^{})=0.04\times (UB)_{JC},\hfill \\ (U_{JC}U_{EIS})=0.03\times (UB)_{JC},\hfill \\ (B_{JC}B_{EIS})=0.14\times (UB)_{JC},\hfill \\ (B_{JC}B_{EIS})=0.24\times (BV)_{JC},\hfill \\ (V_{JC}V_{EIS})=\text{}0.11\times (BV)_{JC},\hfill \\ (V_{JC}V_{EIS})=\text{}0.21\times (VR)_{JC},\hfill \\ (R_{JC}R_{EIS})=\text{}0.04\times (VR)_{JC},\hfill \\ (R_{JC}R_{EIS})=\text{}0.04\times (RI)_{JC},\hfill \\ (I_{JC}I_{EIS})=0.25\times (RI)_{JC},\hfill \end{array}`$

These equations can be inverted to a more useful form, expressing the color terms as a function of the observed colors yielding:

$`\begin{array}{c}(U_{JC}U_{EIS}^{})=0.04\times (U^{}B)_{EIS},\hfill \\ (U_{JC}U_{EIS})=0.03\times (UB)_{EIS},\hfill \\ (B_{JC}B_{EIS})=0.36\times (BV)_{EIS},\hfill \\ (V_{JC}V_{EIS})=\text{}0.08\times (BV)_{EIS},\hfill \\ (R_{JC}R_{EIS})=\text{}0.05\times (VR)_{EIS},\hfill \\ (I_{JC}I_{EIS})=0.35\times (RI)_{EIS},\hfill \end{array}`$

These values have a typical error of $`\pm 0.03`$ and are in good agreement with those computed using the response function of the filters considered, except for $`U^{}`$ for which the difference in the color term can be $`0.2`$ mag. This is because the UV-response of the system is at the time of writing still poorly determined.

The magnitudes were also corrected for galactic absorption, using $`E(BV)=0.0089`$ as derived from Schlegel, Finkbeiner & Davis (1998), yielding $`A_U=0.05`$ mag, $`A_B`$=0.04 mag, $`A_V`$=0.03 mag, $`A_R`$=0.02 mag and $`A_I`$=0.02 mag. All magnitudes have been converted to the $`AB`$ system, unless otherwise specified, using the following relations: $`U_{AB}^{}=U^{}+1.04`$; $`U_{AB}=U+0.80;B_{AB}=B0.11;V_{AB}=V;R_{AB}=R+0.19;`$ and $`I_{AB}=I+0.50`$.

### 3.4 Image Products

The images being released are fully astrometrically (conical equal-area projection; COE) and photometrically calibrated (normalized to 1 sec exposure time). The FITS files being provided include both the pixel and weight maps as FITS extensions. A stand-alone software tool to convert the COE images to best-fit TAN projection by modifying the image headers is available through the EIS web pages. The astrometric information is stored in the world coordinate system (WCS) keywords in the FITS headers. Similarly the photometric calibration is provided by the zero-point and its error available in the header keywords ZP and ZP\_ERR. The zero-point includes the normalized zero-point of the photometric solution and the atmospheric extinction correction. In the header one can also find a product identification number (P\_ID) which should always be used as reference. The headers also provide information on the number of stacked frames and the on-source total integration time. Additionally, the seeing obtained by measuring the FWHM of bright stars in the final stacked images is stored in the header keyword SEE\_IMA.

High-resolution images in the different passbands can be found at the URL “http://www.eso.org/science/eis/”. The image quality of the final images can be assessed by investigating the nature and the amplitude of the PSF distortions. These are measured by using the two components of the so-called polarization vector. These components are derived from the major- (A) and minor-axis (B) and the position angle ($`\theta `$) as measured by SExtractor for point-like sources using the following expressions:

$$e_1=\frac{AB}{A+B}\mathrm{cos}(2\theta )$$

$$e_2=\frac{AB}{A+B}\mathrm{sin}(2\theta )$$

The quantities $`e_1`$ and $`e_2`$ give the amplitude of the distortion. As an illustration, Figures 5 and 6 show the results obtained from the final $`R`$-band image, which consists 30 co-added WFI exposures. Figure 5 shows the distribution of the amplitude of the components, $`e_1`$ and $`e_2`$, for the individual point sources observed over the entire field and in Figure 6 the polarization vector field smoothed over a scale of 7 arcmin. The scale of the amplitude of the vectors is set by the length of the small thick line shown on the upper left corner of the figure, which corresponds to an amplitude of 0.0015 (0.15%). These results attest to the excellent quality of the images produced by WFI leading to a uniform and small PSF distortion throughout the frame with an $`rms`$ amplitude $`<2\%`$. Similar results are obtained for the other passbands.

## 4 Source Lists

Source extraction was performed using the latest version of SExtractor software (Bertin & Arnouts 1996; ver. 2.2.1). This version resolves a subtle problem, first identified in the prompt release of the Pilot Survey, which caused the creation of regions devoid of sources when applying SExtractor to WFI frames with weights (Bertin 2001). Detection was carried out separately using the co-added image of each passband and field. The main parameters in the detection are the smoothing kernel, taken to be a Gaussian with a FWHM equal to 0.8 of that of the PSF measured on the frame; the minimum number of connected pixels above the detection threshold, depending on the seeing taken to be between 8 and 12 pixels, and the SExtractor detection threshold, taken to be between 0.55 (for $`U`$) and 0.7 ($`R`$ and $`I`$) depending on the seeing. As an illustration, the tabulation of the first 40 entries in the $`R`$-band source catalog is presented in Table 7. All magnitudes are given in the $`AB`$ system. The table lists:

Column (1): the EIS identification name

Columns (2) and (3): right ascension and declination (J2000.0);

Columns (4)-(9): total, isophotal and aperture (3 arcsec diameter) magnitudes and respective errors. The first two magnitudes correspond to the mag\_auto and mag\_iso magnitudes measured by SExtractor. The magnitudes have been corrected for Galactic extinction taken from Schlegel et al. (1998). The errors are those estimated by SExtractor and include only the shot-noise of the measured source and background counts. Only objects detected with signal to noise $`S/N3`$ (based on the isophotal magnitude errors) are included.

Column (10): an estimate of the $`S/N`$ of the detection, from the inverse of the errors estimated for the isophotal magnitude;

Columns (11): the isophotal area $`A`$ of the object in square arcsec;

Column (12): the half-light radius $`r_h`$ in arcsec;

Column (13) and (14): minor to major-axis ratio and the position angle;

Column (15): the stellarity index computed by SExtractor;

Column (16): SExtractor flags (see Bertin 1998)

Column (17): EIS flags; these flags are used to identify objects in regions with very low weight or close to bright objects. Objects detected in regions with weight $`30\%`$ of the maximum weight (proportional to the total integration time) and not affected by bright stars have EIS flag=0. Objects with weight $`<30\%`$, located at the edges of the frame, have EIS flag=1. Finally EIS flag=2 indicates that the object is located in a region masked out due to the presence of a bright object.

Similar tables for the other passbands are available upon request in ASCII format. Catalogs for the infrared mosaic and optical/infrared color catalogs will be distributed as soon as they become available. Note that besides distributing catalogs in tabular form as done here, catalogs in FITS formats will become part of the distribution of final products. From the inspection of the $`R`$-band weight map (see Figure 3) one sees that in contrast to catalogs extracted from single chip detectors, describing source catalogs drawn from multi-passband dithered observations of CCD mosaics requires significantly more information than can reasonably be provided in a single ASCII catalog. While the latter is adequate if the derived objects list is to be used for optical cross-identification, it is clearly not sufficient if one is interested in drawing statistical samples with well-defined limiting magnitudes and color coverage. To cope with these deficiencies future releases will be made using FITS tables which will include: a FITS header with all the relevant information about the observations (pointing, filter, instrument), a FIELDS table which will include: the SExtractor parameters; the magnitude and stellarity index chosen for the star/galaxy classification; magnitude corrections (e.g., AB, extinction) appropriate for the filters used; a variety of limiting magnitudes (point source, turnover, completeness, spurious objects) including the brightest limiting magnitude to extract a homogeneous magnitude-limited sample; the effective area corresponding to this sample taking into account the area removed at the outer edges of the image based on the weight-map and the masks around bright objects both now being produced automatically by the pipeline according to user-specified parameters. At the same time no objects will be removed but will carry both their SExtractor and EIS flags, with the vertices of the corresponding masks also being included in the FITS catalogs in a specific fields table. This will also apply to color catalogs which, in addition to information similar to single passband catalogs, must also include others such as the definition of color context to clearly indicate the area of intersection of images in different passbands. In summary, the goal is to produce standard catalogs from which more complex analysis can be carried out with or without the images and associated weight maps.

The output of the new standardized procedures for single passband catalog production is shown in Figure 7 for the catalog extracted from the final $`R`$band image. The figure shows the outer bounding box defined based on the percentage of the peak exposure time one wants to consider and the masks around saturated or bright objects. Both the magnitude of the objects to be masked and the size of masks can be set by the user. In this example one sees that not all possibly spurious sources can be eliminated automatically. In particular, in the lower right side of the figure ($`RA=53.35`$, $`Dec=27.9`$ degrees) one sees the effect of a ghost image of a bright star which produces a ring-shaped excess of objects. While attempts are currently being made to identify ghosts and automatically eliminate these regions it is clear from the diversity of features that in practice, further by-hand clipping of the image is required depending on the specific application of the catalog. Tools for doing so are also implemented in the pipeline with the results being stored in the FITS catalog.

Each catalog produced by the pipeline goes through an automatic verification process which carries out several tests to check the reliability of the products by comparing them with other empirical data and model predictions. The results of some of these tests are described in the next section.

## 5 Discussion

### 5.1 Astrometry and Photometry

In order to externally verify the accuracy of the astrometric calibration the source lists extracted from the different passbands were compared with that available from the 2MASS survey (Cutri et al. 2000). The number of objects found in common ranges from 286 to 415 depending on the passband considered. Using “good” objects, considering the various flags, in both catalogs the relative positions of coincident pairs were computed and the distribution of the differences in the case of the $`R`$-band is shown in Figure 8. One finds an offset of about 0.12 arcsec in right ascension, a negligible shift in declination and an $`rms`$ in the differences of $`<0.23`$ arcsec, independent of the passband considered. This shows that the astrometric solutions for the different passbands are essentially identical and an overall accuracy of $`0.16`$ arcsec, probably limited to the internal accuracy of the reference catalog, is achieved. Note that the shift of the optical data relative to 2MASS is smaller than that obtained with the infrared data. This, however, does not imply any inconsistency in the relative astrometry of the infrared and optical EIS data but rather it reflects the problem described in Vandame et al. (2001) with the astrometry of the pre-release of GSC-II. This problem will be fixed with the official release of the GSC-II.

An independent check of the astrometry can also be carried out by comparing the present astrometry with that determined by Wolf et al. (2000; COMBO survey) who recently obtained deep exposures of the same field in different passbands using the same instrument. In particular, their $`R`$-band data reaches $`R=26.0`$ (Vega) at $`5\sigma `$ for a total integration time of 23700 sec. Within the overlapping area of both surveys, there are: $`36800`$ objects in common; $``$ 17500 detected only by COMBO; $``$ 4400 only by EIS; and $``$ 4200 have multiple associations, which are not considered below. Figure 9 shows the result of the comparison of the positions for about 8500 objects in common. These objects include only objects with stellarity index $`>`$0.9 in the COMBO $`R`$-band source list and satisfying $`R_{EIS}<25`$ mag and $`S/N(EIS)>`$4. From this comparison one finds a relative offset of 0.15 arcsec in right ascension and $``$0.24 arcsec in declination, but more importantly a remarkably small $`rms`$ of 0.12 arcsec in both directions, possibly suggesting an accuracy of $`<0.10`$ arcsec for each individual catalog. This is well within the requirements for slit/fiber positioning, a top requirement for the public survey. These results are insensitive to the exact way the sample is chosen.

The COMBO data also enables a comparison of the photometric calibration. Figure 10 shows the results of this comparison in the Vega system from which one finds a zero-point offset of $`0.01`$ with an $`rms`$ of about 0.01 mag for $`16<R_{EIS}<20`$ mag, while for objects brighter than $`R16`$ mag saturation sets in. For magnitudes fainter than $`R20`$ the scatter increases and beyond $`R22`$ one sees the Malmquist bias effect.

### 5.2 Number Counts

A simple statistic that can be used to verify the overall characteristics of the derived catalogs and the parameters chosen for the star/galaxy classification is to compare the number counts of stars and galaxies with model predictions and other empirical data.

Figure 11 show a comparison of the stellar counts as a function of the total magnitude (Vega system) with those obtained using the Galactic model described by Mendez & van Altena (1996), using the standard parameters described in their Table 1. Considering that no attempt has been made to choose model parameters to fit the observed counts the agreement between observations and the model predictions is very good. In the $`B`$-band, where the present survey is deepest, there seems to be an excess of stellar objects. This is possibly due to QSOs, which are estimated to have a surface density of 120 objects per square degree per 0.5 mag interval at $`B`$=24 (Hatziminaoglou et al. , 2001). This effect is not expected in $`V`$ and $`I`$ because the surface density of QSOs at the magnitude where the stellar counts drop is expected to be small, of order 30-40 per square degree per 0.5 mag interval. Currently, the comparison of the star counts is only possible for the three passbands considered. However, this problem is being addressed and soon a more generic galactic model will be available in the EIS pipeline to compare the predicted counts in all passbands using the transmission curves for the filters being used in a given survey setup (Girardi 2001).

Similarly, Figure 12 shows the comparison of galaxy counts as a function of the total magnitude, for each of the available passbands, with those obtained for other data sets. In this case galaxies were defined as objects with a stellarity index less than 0.95 (0.90 in the $`I`$band) or fainter than the classification limit of each passband. The sample includes only objects with $`S/N>3`$ and with the appropriate SExtractor and EIS flags and the normalization of the counts takes into account the effective area of each catalog. The $`U`$band counts of Guhathakurta et al. (1990) and $`I`$-band counts of Postman et al. (1998) have been corrected assuming an $`AB`$ correction of 0.8 and 0.45, respectively. As can be seen the observed counts are in good agreement with previously published results indicating that the catalogs being distributed are statistically consistent.

### 5.3 Survey Performance

The characteristics of the data obtained by the present survey are summarized in Table LABEL:tab:quality which lists: in column (1) the filter; in column (2) the seeing of the final co-added image; in columns (3) and (4) the $`5\sigma `$ and $`3\sigma `$ limiting $`AB`$ magnitude measured within an aperture $`2\times `$FWHM; and in column (5) the number of objects with $`S/N3`$, i.e., the number of entries in the catalogs being distributed. At the $`5\sigma `$ limiting magnitudes the fraction of spurious objects with $`S/N5`$ is estimated to be $`<2\%`$ for all bands. The fraction of spurious objects was estimated by creating catalogs from the survey images multiplied by $`1`$. The mock images were then used as input to the catalog production and from the comparison between the catalogs of the mock and real objects the fraction of false-positive detections was estimated.

The fact that the COMBO data reaches fainter magnitudes can also be used to obtain a first estimate of the completeness limit of the CDF-S data in the $`R`$band. Figure 13 shows the ratio as a function of the $`R_{COMBO}`$ of objects in both catalogs. From this plot one finds that the $`R`$band catalog is complete to $`R_{AB}24.7`$ and 50% complete at $`R_{AB}25.7`$. This estimate is probably conservative. A first indication is that there are a large number of objects detected only in the data presented here, most of them at the faint end $`R>25.2`$. Second, the number of spurious detections at $`R26`$ is still relatively small ($`<30\%`$). Clearly, a determination of the completeness of the catalog from the CDF-S data themselves is highly desirable and will be implemented in the near future.

The limiting magnitudes reached by the present observations are, in general, brighter than those originally proposed ($`U_{AB}=26.8;B_{AB}=26.0;V_{AB}=26.0;R_{AB}=26.3;I_{AB}=26.0`$) For $`BVR`$ an increase of the integration time by a factor of two is required to reach the originally stated $`5\sigma `$ limits, representing a 50% increase in the total integration time requested per pointing. Certainly more serious is the requirements in $`U`$ and $`I`$, even if more efficient filters are used. From the present data one finds that the new $`U^{}`$ filter has an efficiency at least 35% larger than the original $`U`$. A similar increase in efficiency in $`I`$ is possible using a filter similar to the one available for EMMI for the EIS-WIDE survey. However, even with the new filters the time required to reach the desire magnitudes in $`U`$ and $`I`$ involve a large increase (a factor $`>4`$) in their already long integration times. Such an increase would make the desired area coverage unattainable. On the other hand, the limiting magnitudes obtained by the present survey compare favorably with those reached by EIS-DEEP from which Lyman-break galaxy candidates drawn from those samples (da Costa et al. 1998; Rengelink et al. 1998) have been successfully confirmed spectroscopically at $`z>2.8`$ (Cristiani et al. 2000) strongly suggesting that the primary goal of the present survey is indeed possible. A more detailed discussion of this point will be presented elsewhere when attempts to combine the data obtained using the two different $`U`$ filters will be made (Arnouts et al. 2001, Benoist et al. 2001b).

Finally, as mentioned in Section 2, the CDF-S has already been observed both in optical and infrared passbands as part of the EIS-DEEP survey using the NTT (Rengelink et al. 1998). These old observations are being re-analyzed, using the new methods of the pipeline and compared with the present observations elsewhere (Benoist et al. 2001a). The infrared data from this survey is particularly important as it complements the infrared CDF-S survey recently completed (Vandame et al. 2001).

## 6 Summary

This paper presents the first optical results of the ongoing multicolor Deep Public Survey being conducted by the EIS program. The observations and the data for the first pointing completed as part of this program which corresponds to the Chandra CDF-S field are described. The set of products being publicly released includes fully calibrated pixel maps and source lists covering the CDF-S region. Given the current interest in this area of the sky both the optical data presented here, covering $`0.25`$ square degrees, and the infrared survey recently released, covering an area of $`0.11`$ square degrees, are being made available world-wide. Hopefully, this will, together with the Chandra X-ray data, contribute to making this region of the sky a natural target for multi-wavelength observations in the southern hemisphere and an ideal area for cosmological studies. Other derived products from the EIS program will be described in forthcoming papers in this series.

###### Acknowledgements.

We thank all of those directly or indirectly involved in the EIS effort. Our special thanks to K. Meisenheimer and C. Wolf for making their $`U`$-band images and $`R`$-band catalogs available to EIS, E. Bertin for his quick response to our request for help, A. Bijaoui for allowing us to use tools developed by him and collaborators over the years and past EIS team members for building the foundations of this program. We would also like to thank A. Renzini and past and present members of the Working Group for Public Surveys. The Guide Star Catalogue-II is produced by the Space Telescope Science Institute in collaboration with the Osservatorio Astronomico di Torino. Space Telescope Science Institute is operated by the Association of Universities for Research in Astronomy, for the National Aeronautics and Space Administration under contract NAS5-26555. Additional support is provided by the Association of Universities for Research in Astronomy, the Italian Council for Research in Astronomy, European Southern Observatory, Space Telescope European Coordinating Facility, the International GEMINI project and the European Space Agency Astrophysics Division. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation. This research has made use of the simbad database, operated at CDS, Strasbourg, France.


text:
# I INTRODUCTION

## I INTRODUCTION

A survey of the molecule-like structure is one of the most challenging subjects in light neutron-rich nuclei. In Be isotopes, recently, decay into fragments of He isotopes (4He, 6He, 8He) has been observed from the excited states of 10Be and 12Be, and the presence of a two-center configuration is suggested. From the theoretical side, these states are studied by various models by which the molecular-orbital nature of the weakly bound neutrons around the two $`\alpha `$-clusters has been revealed .

As for multi-cluster configurations in light $`\alpha `$-nuclei beyond the two-center systems, the existence of $`N\alpha `$ states has been predicted around the threshold energy in the so-called Ikeda diagram . For example, it has been suggested that the second $`0^+`$ state of 12C has a 3$`\alpha `$-like molecular configuration. However, according to many theoretical analyses, the state is not necessary to have a linear-chain of 3$`\alpha `$, but is described as a weak-coupling state with a triangular shape or the 8Be+$`\alpha `$ configuration. Furthermore, a rotational band of 16O with a very large moment of inertia has been observed around the 4$`\alpha `$ threshold energy region, and the 4$`\alpha `$ linear-chain have been discussed . However, the 4$`\alpha `$-chain has not been experimentally confirmed yet.

Recently, the discussions of the well-developed cluster structure are extended to the neutron-rich nuclei, and the role of valence neutrons which stabilize the linear-chain structure has been pointed out. For example, von Oertzen has extended his analyses for the molecular structure in Be isotopes  to C isotopes, and the linear-chain state consisting of 3$`\alpha `$ and valence neurons around it has been speculated. Even if the 3$`\alpha `$-system without valence neutrons (12C) does not have a linear-chain structure, the valence neutrons around it are expected to increase the binding energy and stabilize the linear-chain state.

In this paper, the molecular-orbital approach introduced for the Be isotopes is applied to a study of C isotopes (A=12, 14, 16). The appearance of prolonged structure with the $`\alpha `$-$`\alpha `$ core in 10Be and 12Be has been successfully described by the $`\alpha `$+$`\alpha `$+n+n+$``$ model, where the orbits for the valence neutrons are classified based on the molecular-orbit (MO) model. In the present study, MO is extended to C isotopes, we focus on the stabilization of the linear-chain structure with the 3$`\alpha `$-core to which valence neutrons are added.

## II Model

The total wave function of a microscopic $`\alpha `$+$`\alpha `$+$`\alpha `$+$`n`$+$`n\mathrm{}`$ model is fully antisymmetrized and expressed by a superposition of basis states centered to different relative distances between the $`\alpha `$-clusters ($`d`$) with various configurations of the valence neutrons ($`c1`$, $`c2`$ $``$) around the $`\alpha `$-clusters:

$$\mathrm{\Phi }_{MK}^J=\underset{d,c1,c2}{}P_{MK}^J𝒜[\varphi _1^{(\alpha )}\varphi _2^{(\alpha )}\varphi _3^{(\alpha )}(\varphi _1^{c1}\chi _1)(\varphi _2^{c2}\chi _2)\mathrm{}].$$
(1)
The projection to the eigen-states of angular momentum $`(P_{MK}^J)`$ is performed numerically. Each $`\alpha `$ cluster consisting of four nucleons is described by Gaussians ($`G_{\alpha i}`$) centered at $`R_{\alpha i}`$ and spin-isospin wave function$`(\chi )`$:

$$\varphi _i^{(\alpha )}=G_{R_{\alpha i}}^pG_{R_{\alpha i}}^pG_{R_{\alpha i}}^nG_{R_{\alpha i}}^n\chi _p\chi _p\chi _n\chi _n,i=1,2,3.$$
(2)
$$G_R=\left(\frac{2\nu }{\pi }\right)^{\frac{3}{4}}\mathrm{exp}[\nu (\stackrel{}{r}\stackrel{}{R})^2],\nu =1/2\beta ^2,$$
(3)
where, the oscillator parameter $`(\beta )`$ is equal to 1.46 fm. For the linear-chain state, the values of $`\{R_{\alpha i}\}`$ are $`d`$, $`0`$, and $`+d`$ on the $`z`$-axis. Each valence neutron ($`\varphi _i^{ci}\chi _i`$) around the $`\alpha `$-$`\alpha `$-$`\alpha `$ core is expressed by a linear combination of local Gaussians:

$$\varphi _i^{ci}\chi _i=\underset{j}{}g_jG_{R_j}\chi _i.$$
(4)
These valence-neutron orbits are classified according to the MO picture. The orbit of the valence neutron perpendicular to the $`z`$-axis of the 3$`\alpha `$ linear-chain core is called $`\pi `$-orbit, and one along the $`z`$-axis is called $`\sigma `$-orbit. The antisymmetrization imposes the forbidden space for the valence neutrons; the $`\pi `$-orbit must have at least one node perpendicular to the $`z`$-axis, and the $`\sigma `$-orbit must have at least three nodes since three $`\alpha `$-clusters along the $`z`$-axis already occupy the orbitals with $`n_z=0,1,2`$. The amplitude of the $`\pi `$\- and the $`\sigma `$-orbit are schematically drawn in Fig. 1.

In the present framework, each valence-neutron orbit is introduced to have a definite $`K^\pi `$ value at the zero limit of centers of local Gaussians $`(\{R^j\})`$ describing the spatial distribution of the orbit. The precise positions of $`\{R^j\}`$ are determined variationally before the angular-momentum projection. Since the values of $`\{R^j\}`$ are optimized to be finite, the orbits are not exactly the eigen-state of $`K^\pi `$, and are labeled as $`\overline{K}^\pi `$. For the $`\pi `$-orbit with $`\overline{K}^\pi =3/2^{}`$ ($`|3/2_\pi ^{}`$), the spatial part and the spin part of $`\overline{K}`$ are defined to be parallel ($`rY_{11}|n`$), for which the spin-orbit interaction acts attractively. At the same time, $`|3/2_\pi ^{}`$ is described as a linear combination of two orbits centered at right- and left-hand side of the system based on the MO picture.

$$|3/2_\pi ^{}=\frac{1}{\sqrt{N_\pi }}\{\frac{1}{\sqrt{2}}(p_x+ip_y)_{+a}+\frac{1}{\sqrt{2}}(p_x+ip_y)_a\}|n.$$
(5)
$$(p_x)_{\pm a}=G_{\pm a\stackrel{}{e}_z+b\stackrel{}{e}_x}G_{\pm a\stackrel{}{e}_zb\stackrel{}{e}_x},(p_y)_{\pm a}=G_{\pm a\stackrel{}{e}_z+b\stackrel{}{e}_y}G_{\pm a\stackrel{}{e}_zb\stackrel{}{e}_y}.$$
(6)
Here, $`(p_x+ip_y)_{\pm a}`$ denotes the $`p`$-orbit centered at $`\pm a`$ on the $`z`$-axis, and these variational parameters $`a`$ and $`b`$ are optimized by using the Cooling Method in antisymmetrized molecular dynamics (AMD)  for each basis state. Furthermore, the $`|1/2_\pi ^{}`$ orbit where the spin-orbit interaction acts repulsively can also be defined by changing the spin-direction of $`|3/2_\pi ^{}`$, where the spatial part and the spin part of $`\overline{K}`$ are anti-parallel ($`rY_{11}|n`$):

$$|1/2_\pi ^{}=\frac{1}{\sqrt{N_\pi }}\{\frac{1}{\sqrt{2}}(p_x+ip_y)_{+a}+\frac{1}{\sqrt{2}}(p_x+ip_y)_a\}|n.$$
(7)
The distribution of the $`\sigma `$-orbit is just along the 3$`\alpha `$ axis, then it is introduced to have three nodes. $`|1/2_\sigma ^{}`$ is represented as a linear combination of three orbits with $`\overline{K}^\pi =1/2^{}`$, whose centers are $`+a`$, $`0`$, and $`a`$ on the $`z`$-axis.

$$|1/2_\sigma ^{}=\frac{1}{\sqrt{N_\sigma }}\{(p_z)_{+a}(p_z)_0+(p_z)_a\}|n,$$
(8)
$$(p_z)_{\pm a}=G_{a\stackrel{}{e}_z+b\stackrel{}{e}_z}G_{a\stackrel{}{e}_zb\stackrel{}{e}_z},(p_z)_0=G_{b\stackrel{}{e}_z}G_{b\stackrel{}{e}_z}.$$
(9)
These three orbits ($`|3/2_\pi ^{}`$, $`|1/2_\pi ^{}`$ and $`|1/2_\sigma ^{}`$) are the basic building blocks for the molecular-orbital structure. Also, $`|3/2_\pi ^{}`$, $`|1/2_\pi ^{}`$ and $`|1/2_\sigma ^{}`$ orbits are introduced by taking the time reversal of $`|3/2_\pi ^{}`$, $`|1/2_\pi ^{}`$, and $`|1/2_\sigma ^{}`$ orbits, respectively.

The Hamiltonian and the effective nucleon-nucleon interaction are the same as in Refs., and parameters of Volkov No.2 for the central part and the G3RS spin-orbit term for the spin-orbit part are determined from the $`\alpha +n`$ and $`\alpha +\alpha `$ scattering phase shifts, and the binding energy of deuteron is also reproduced with these parameters.

## III Results and discussions for bending stability

In the following part, we show the calculated results for the stability of the linear-chain state for various configurations. The isotopes and configurations which we take into account are 12C, 14C$`(3/2_\pi ^{})^2`$ (two $`n`$’s in the $`\pi `$-orbits), 14C$`(1/2_\sigma ^{})^2`$ (two $`n`$’s in the $`\sigma `$-orbits), 16C(($`3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$ (four $`n`$’s in the $`\pi `$-orbits) and 16C(($`3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ (two $`n`$’s in the $`\pi `$-orbits and two $`n`$’s in the $`\sigma `$-orbits). As schematically shown in Fig. 2, two variational paths are introduced corresponding to the breathing-like (Fig. 2 (a)) and the bending-like (Fig. 2 (b)) degrees of freedom. The parameters $`d`$ and $`\theta `$ stand for the $`\alpha `$-$`\alpha `$ distance and the bending angle of the 3$`\alpha `$-core, respectively.

Firstly, we show the $`0^+`$ energy curves for the linear-chain structure against the breathing-path in Fig. 3. It is found that the energy pocket around $`d=3`$ fm becomes deeper as the increase of number of valence neutrons in the $`\pi `$-orbit (12C $``$ 14C$`(3/2_\pi ^{})^2`$ $``$ 16C(($`3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$). Using our framework, the binding energy of one $`\alpha `$-cluster is calculated to be 27.5 MeV, then the $`\alpha `$+$`\alpha `$+$`\alpha `$ +$`n`$-neutrons threshold energy is $`82.5`$ MeV. The 3$`\alpha `$-system (12C) has minimal energy around $`d=3.5`$ fm, however, this is too shallow to conclude the stability of the linear-chain state. On the contrary, in 14C$`(3/2_\pi ^{})^2`$, there appears evident minimal energy around $`d=3`$ fm. The energy ($`82`$ MeV) is lower than 12C by 11 MeV and the energy pocket is much deeper. This energy corresponds to the excitation energy of 18 MeV from the ground state calculated with an equilateral-triangle configuration for the 3$`\alpha `$-core, which is $`101.2`$ MeV. The 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$ configurational state is most stable among states studied and has an energy pocket of $`86`$ MeV, where the $`\alpha `$-$`\alpha `$ distance is $`d=2.5`$ fm, shorter than those for 12C and 14C$`(3/2_\pi ^{})^2`$. Therefore, the $`\pi `$-orbit is found to stabilize the linear-chain structure as the increase of valence neutrons (12C $``$ 14C$`(3/2_\pi ^{})^2`$ $``$ 16C(($`3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$).

Next, we discuss the case where the valence neutrons occupy the $`\sigma `$-orbit. 14C$`(1/2_\sigma ^{})^2`$ has an excitation energy higher by about 14 MeV in comparison with 14C$`(3/2_\pi ^{})^2`$. It is rather surprising that the difference is only 14 MeV in spite of the fact that $`3/2_\pi ^{}`$ has only one node and $`1/2_\sigma ^{}`$ has three nodes. This is because the $`\sigma `$-orbit is along the $`\alpha `$-$`\alpha `$-$`\alpha `$ core: The higher nodal orbits along the symmetry axis become low-lying as a result of the clustering of the core. The $`\sigma `$-orbit enhances the prolonged shape of the 3$`\alpha `$-core, and the optimal $`d`$-value is $``$3.5 fm. However, the 14C$`(1/2_\sigma ^{})^2`$ case has no deep pocket enough to be stabilized against the breathing-like path. When two more valence neutrons occupy the $`\pi `$-orbit, although this minimal energy is higher by 5 MeV than that of 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$, 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ has the minimal $`0^+`$ energy of $`77`$ MeV. The calculated energy pocket is deep enough to guarantee the stability for the breathing-like path. Therefore, it is summarized that the linear-chain structure cannot be stabilized against the breathing-like mode with the $`\sigma `$-orbit only in 14C$`((1/2_\sigma ^{})^2)`$, however, the neutrons in $`\pi `$-orbits increase the binding energy and prevent the break-up of the system in 14C$`((3/2_\pi ^{})^2)`$, 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$, and 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$.

Finally, the stability of these linear-chain states against the bending-like path is examined. The $`0^+`$ energy curves of 12C, 14C$`((3/2_\pi ^{})^2)`$, 14C$`((1/2_\sigma ^{})^2)`$, 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$, and 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ against the $`\theta `$-value are shown in Fig. 4. Except for the case of 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$, the curvature of these states is rather monotonic and the energy minimum does not clearly appear. In 14C, the orthogonality between the linear-chain configuration and low-lying states with equilateral-triangle configuration of 3$`\alpha `$ is taken into account. This effect works but not sufficiently to push up the energy of the state with a finite bending angle. However, as clearly seen in Fig. 4, the 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ case shows a sharp increase of the $`0^+`$ energy as the increase of the bending angle and is found to be stable against the bending-like path. This feature is much different from 12C, 14C$`((3/2_\pi ^{})^2)`$, 14C$`((1/2_\sigma ^{})^2)`$, and 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$ cases. From the analysis above, the linear-chain configuration can be stabilized against the breathing-like path by neutrons in the $`\pi `$-orbit (14C$`((3/2_\pi ^{})^2)`$, 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$, 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$), but 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$ is only the case which is stable also against the bending-like path.

We further discuss the reason for the sharp increase of the $`0^+`$ energy in 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ against the bending-like path. The $`0^+`$ energy increases by 15.7 MeV from $`\theta =0^o`$ to $`\theta =30^o`$, in which the kinetic energy part is 10.3 MeV. To understand the energy increase with the increase of bending angle $`\theta `$ of this case, we calculate and compare the overlap between the wave functions with $`\theta =0^o`$ and $`\theta =30^o`$ for various configurations in Table I, where the results are shown. In 12C, the wave functions with $`\theta =0^o`$ and $`\theta =30^o`$ have the squared overlap of 0.91, and 14C$`((3/2_\pi ^{})^2)`$ has almost the same value. 14C$`((3/2_\sigma ^{})^2)`$ has the value of 0.85, smaller than 14C$`((3/2_\pi ^{})^2)`$ by only 6 $`\%`$, and 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$ has almost the same value as the 14C$`((1/2_\sigma ^{})^2)`$ case. This result shows that the overlaps additionally decreases a little for the $`\sigma `$-orbital neutrons, and also for the $`\pi `$-orbital neutrons as the increase of the valence neutrons. In spite of these, the overlap between $`\theta =0^o`$ and $`\theta =30^o`$ for 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ case shows a significantly large decrease to 0.60. 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ is only the configuration which shows drastic decrease of the overlap between $`\theta =0^o`$ and $`\theta =30^o`$.

As discussed in following part, it can be known that the drastic decrease as the increase of the bending angle is due to the increase of overlap between two neutrons in the $`\pi `$-orbit and two neutrons in the $`\sigma `$\- orbit. When there arises the overlap between them, the overlap component in the total wave function is diminished due to Pauli exclusion principle, that is, so-called Pauli blocking. Therefore, the physical state can be expressed by the modified wave function which is made by subtracting the overlap component from the original wave function. Since the energies of the $`\pi `$\- and the $`\sigma `$\- orbits discussed here are relatively low, the modified wave function involves larger components of higher excitation energy in comparison with the wave function at $`\theta =0`$ free from the Pauli blocking. As a result, the Pauli blocking due to the increase of overlap between valence four neutrons is considered to bring the increase of energy proportional to the decrease of squared overlap. This is a possible explanation for rapid increase of the energy against the bending angle. On this view we can estimate the increase of the energy as follows. Firstly, we compare the cases with Pauli blocking and without Pauli blocking. If there is no Pauli blocking, the squared overlap of 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ between $`\theta =0^o`$ and $`\theta =30^o`$ is estimated as follows:

$$P(^{16}\mathrm{C}_{\pi \sigma })_{NoPauli}=\frac{P(^{14}\mathrm{C}_\pi )\times P(^{14}\mathrm{C}_\sigma )}{P(^{12}\mathrm{C})}=0.85.$$
(10)
Here, $`P(^{14}\mathrm{C}_\pi )`$, $`P(^{14}\mathrm{C}_\sigma )`$, and $`P(^{12}\mathrm{C})`$ are squared overlaps between $`\theta =0^o`$ and $`\theta =30^o`$ of 14C$`((3/2_\pi ^{})^2)`$, 14C$`((1/2_\sigma ^{})^2)`$, and 12C shown in Table I. In the actual squared overlap of 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ shown in Table I, the Pauli blocking effect is automatically included, and the value is $`P(^{16}\mathrm{C}_{\pi \sigma })_{WithPauli}=0.60`$. This difference corresponds to the component excited to higher shells because of the Pauli blocking effect, and this is 25 $`\%`$.

$$P(^{16}\mathrm{C}_{\pi \sigma })_{NoPauli}P(^{16}\mathrm{C}_{\pi \sigma })_{WithPauli}=0.25$$
(11)
Secondly, we estimate the increase of the kinetic energy when two valence neutrons in the $`\pi `$\- or the $`\sigma `$-orbits are excited to higher shells by 25 $`\%`$ due to the Pauli blocking effect.

$$\mathrm{\Delta }E=2\times \mathrm{}\omega \times 0.25=10\mathrm{M}\mathrm{e}\mathrm{V}.$$
(12)
This value is consistent with the increase of the kinetic energy at $`\theta =30^o`$ mentioned above. Therefore, we can conclude that this is one of the most important mechanisms which stabilizes the linear-chain state in 16C.

## IV Summary and discussions

It is summarized that the linear-chain structure of 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ with 3$`\alpha `$ core is the only case to have the simultaneous stabilities for the breathing-like break up path and for the bending-like path among 12C, 14C, and 16C. Other configurations, such as 14C$`((3/2_\pi ^{})^2)`$ and 16C$`((3/2_\pi ^{})^2(1/2_\pi ^{})^2)`$ are stable against the breathing-like path but not stable against the bending-like path. A combination of the $`\pi `$\- and the $`\sigma `$\- orbits occupied by four neutrons plays doubly important roles to make a deep energy pocket for breathing-like path and to prevent the bending-like free motion of the system. The 16C$`((3/2_\pi ^{})^2(1/2_\sigma ^{})^2)`$ configuration forms a rotational band with an energy slope of $`\frac{\mathrm{}^2}{2I}=150`$ keV, and since the two-neutron separation-energy of 16C is experimentally found to be $`5.5`$ MeV, the band head energy is expected to corresponds to around 25 MeV in excitation.

We have shown in 16C that the Pauli-blocking effect among the valence neutrons play an important role for the stability of the linear-chain configuration. This effect is expected to be more important as the neutron number increases. We are interested in 18C the appearance of similar linear-chain structure in low-lying region, where all of the three neutron orbits introduced here are occupied ($`(3/2_\pi ^{})^2(1/2_\pi ^{})^2(1/2_\sigma ^{})^2`$).

Experimentally, the observation of the $`\gamma `$-ray from high-spin states is the most probable way to confirm the presence of rotational band with large moment of inertia. Although the present analysis has been restricted to the $`0^+`$ state, we will investigate also high-spin states, where it is expected that a glue-like role of valence neutrons becomes more important to prevent the breathing-like break-up of the system.

In the present analysis, parameter $`d`$ has been fixed when the stability against the bending motion has been examined, and $`\theta `$ has been fixed to zero when the bending motion has been examined. However, to assure theoretically the stability of these linear-chain states, we intend to investigate the case when two degrees of freedom ($`d`$ and $`\theta `$) simultaneously activate. This shall be examined by superposing states on the energy surface of $`d`$ and $`\theta `$ based on generator coordinate method.

The authors thank members of RI beam science laboratory in RIKEN for discussions and encouragements. One of the authors (N.I) thanks Prof. R. Lovas, Prof. W. von Oertzen, Prof. H. Horiuchi, and Dr. Y. Kanada-En’yo for fruitful discussions.