text:
# The ELODIE survey for northern extra-solar planets III. Three planetary candidates detected with ELODIE Based on observations made with the ELODIE echelle spectrograph mounted on the 1.93–m Telescope at the Observatoire de Haute-Provence (CNRS) and with the CORALIE echelle spectrograph mounted on the 1.2–m Euler Swiss Telescope at ESO–La Silla Observatory. The ELODIE and CORALIE measurements discussed in this paper are only available in electronic form at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsweb.u-strasbg.fr/cgi-bin/qcat?J/A+A/

## 1 Introduction

The ELODIE survey for northern extra-solar planets, a programme aiming at detecting and characterising planetary companions in orbit around solar-type stars of the solar neighbourhood, started in 1994. This survey uses the ELODIE echelle spectrograph (Baranne et al. 1996) mounted on the Cassegrain focus of the 1.93–m Telescope of the Observatoire de Haute-Provence (CNRS, France).

The original Mayor & Queloz (1996) sample consisted of 142 stars, including 51 Peg, the first solar-type star known to host a planetary companion (Mayor & Queloz 1995). A new sample of 330 stars was defined in 1996. The new sample, the observing procedure and the achieved precision are described in details in Perrier et al. (2003, Paper I). The complete list of planets detected with ELODIE can also be found in Paper I.

In this paper, we present two systems hosting long-period planets. We present in Sect. 2.2 our ELODIE radial-velocity measurements for the solar-type stars HD 74156 and 14 Her. Updated velocities and orbits for HD 209458 and 51 Peg are presented in Sect. 3. Section 4 contains the description of the ELODIE radial-velocity data and orbital solutions for Ups And, 55 Cnc, 47 Uma, 70 Vir and HD 187123. Our results are summarized in Sect. 5.

## 2 New ELODIE planet candidates

### 2.1 Stellar characteristics

In this section, we briefly describe the main properties of the stars hosting the new planet candidates detected with ELODIE, namely HD 74156 and 14 Her. The properties of the other stars appearing in this paper are described elsewhere.

The main characteristics of HD 74156 and 14 Her are listed in Table 1. Spectral Type, $`m_\mathrm{V}`$, $`BV`$, $`\pi `$, $`\mu _\alpha \mathrm{cos}(\delta )`$ and $`\mu _\delta `$ are from the HIPPARCOS Catalogue (ESA 1997). The effective temperature $`T_{\mathrm{eff}}`$, the surface gravity $`\mathrm{log}g`$, the metallicity $`[`$Fe/H$`]`$, the microturbulance velocity $`\xi _\mathrm{t}`$ and the stellar masses $`M_{}`$ are taken from Santos et al. (2003). The bolometric correction $`B.C.`$ is computed from the effective temperature with the calibration in Flower (1996).

The projected rotational velocity $`v\mathrm{sin}i`$ is derived from the ELODIE cross-correlation function widths using the calibration in Queloz et al. (1998). We only have a $`v\mathrm{sin}i`$ upper limit for 14 Her because this star does not exhibit any rotational broadening measurable with ELODIE. The $`\mathrm{log}R_{\mathrm{HK}}^{^{}}`$ chromospheric activity indicator for 14 Her is from Butler et al. (2003). No quantitative activity value is available for HD 74156. We do not see any trace of emission in the core of the $`\lambda `$ 3968.5 Å Ca ii H line on our ELODIE coadded spectra. This star is thus probably not very active although moderately fast rotating. $`P_{\mathrm{rot}}`$ and $`age(\mathrm{HK})`$ are derived from the activity indicator using the calibrations in Noyes et al. (1984) and Donahue (1993)111also quoted in Henry et al. (1996), respectively.

HD 74156 and 14 Her are not photometrically variable. Their measured HIPPARCOS scatter is 11 and 9 mmag, respectively. Moreover, 14 Her is flagged as constant star in this catalogue.

We have 3 adaptive optics high-angular resolution images of 14 Her. They were obtained between May 1997 and August 2001 using the PUE’O adaptive optics system mounted on the 3.6–m Canada-France Hawaii Telescope. No visual companion is detected around this star.

### 2.2 Radial-velocity analysis

#### 2.2.1 HD 74156 (HIP 42723)

The detection of a two-planet system orbiting HD 74156 was already announced by our team in April 2001 (April 4th 2001 ESO Press Release). The long-period of the outer required however further observations of this system.

We now have 51 ELODIE and 44 CORALIE radial-velocity measurements of HD 74156. These velocities have been obtained between JD = 2 450 823 (January 1998) and JD = 2 452 910 (September 2003). The mean velocity uncertainty of these measurements is 12.7 m s-1 for the ELODIE data and 8.5 m s-1 for the CORALIE data. The total measurement time span is 2087 d.

The fitted orbital elements for HD 74156 are presented in Table 2. Figure 1 shows the fitted orbits and the residuals around the solution. Assuming a primary mass of 1.27 M (Santos et al. 2003), we compute the planetary minimum masses $`m_\mathrm{b}\mathrm{sin}i`$ = 1.86 MJup and $`m_\mathrm{c}\mathrm{sin}i`$ = 6.17 MJup as well as the semi-major axes $`a_\mathrm{b}`$ = 0.294 AU and $`a_\mathrm{c}`$ = 3.40 AU.

#### 2.2.2 14 Her (HD 145675, HIP 79248, GJ 614)

The detection with ELODIE of a long-period planetary companion to 14 Her was announced during the Protostars & Planets IV conference by M. Mayor 1998, oral contribution, also in Marcy, Cochran, & Mayor 2000. An updated orbital solution was presented by Udry et al. (2000). This detection was recently confirmed by Butler et al. (2003).

The orbital solution of Table 2 was obtained by fitting a Keplerian orbit to the 119 ELODIE radial-velocity measurements. These measurements were obtained between JD =2 449 464 (April 1994) and JD =2 452 857 (August 2003). The time span of these data is though 3392 d. Our velocities have a mean velocity uncertainty of 7.2 m s-1.

The residuals to a simple Keplerian fit are abnormally large: 14 m s-1. We thus computed another solution including an additional linear velocity drift. The fitted drift parameter is highly significant: 3.6 $`\pm `$ 0.3 m s-1 yr-1. The residuals around the Keplerian+drift solution are much smaller: $`\sigma _{\mathrm{O}\mathrm{C}}`$ = 11.3 m s-1. This value is still a bit large compared to our uncertainties.

Figure 2 shows the ELODIE radial-velocity data, the fitted Keplerian + linear drift solution and the residuals to the solution. The residuals of the oldest measurements have clearly a higher scatter value. We computed an orbital solution using only the most recent observations (velocities taken after JD =2 450 500, the date of the set up of the double scrambling device). The solution we obtained in this case was compatible with the solution quoted in Table 2 (slightly longer period $`P`$ =1853 $`\pm `$ 35 d, same eccentricity and semi-amplitude) but with a significantly smaller residuals value: $`\sigma _{\mathrm{O}\mathrm{C}}`$ = 8.6 m s-1 ($`\chi _{\mathrm{red}}^2`$ = 1.54). The linear drift slope of this solution was slightly less significant: 3.9 $`\pm `$ 1.1 m s-1 yr-1 which is not a surprise since only one radial-velocity minimum is present in the most recent data. The temporal coverage of these latter data approximately corresponds to the coverage of the Butler et al. (2003) measurements. We also computed an orbital solution to this velocity set without including the drift parameter. The fitted parameters obtained in this case are in good agreement with the Butler et al. (2003) values but the residuals are slightly larger than for the drift-included solution. This shows that the non detection of the velocity drift by these authors is probably due to their much shorter observing span (1882 d). Moreover, their data only cover one velocity maximum and one minimum.

The computed minimum mass of 14 Her b, assuming $`M_1`$ = 0.9 M (Santos et al. 2003), is 4.74 MJup and its semi-major axis is 2.80 AU.

## 3 Updated orbits of known extra-solar planets

In this section, we present updated orbital solutions for HD 209458 and 51 Peg. These two stars host planetary companions that were detected with ELODIE (Mazeh et al. 2000; Mayor & Queloz 1995). The new orbital solutions are based on larger and more precise radial-velocity sets.

### 3.1 HD 209458 (HIP 108859)

The detection of HD 209458 b was simultaneously announced by Mazeh et al. (2000) and Henry et al. (2000). Since the discovery paper, we have obtained more radial-velocity measurements of HD 209458. Moreover, the numerical template we use for computing the cross-correlation function has been changed in order to reduce the effect of telluric atmospheric lines (Pepe et al. 2002). Another effect of this new template is the change of the radial-velocity zero point.

We now have 46 ELODIE and 141 CORALIE radial-velocity measurements of HD 209458 spanning from JD = 2 450 681 (August 1997) and JD = 2 452 566 (October 2002). The mean velocity uncertainty of these measurements is 11.0 m s-1 with ELODIE and 10.2 m s-1 with CORALIE.

HD 209458 b is known to be a transiting extra-solar planet (Charbonneau et al. 2000). It has also been shown by Queloz et al. (2000) that the passage of the planet in front of the star was also affecting the measured radial-velocity (departure from the Keplerian model). We thus have removed all the measurement taken at or near the phase of the transit (measurements with $`\varphi _{\mathrm{trans}}`$0.025 $``$$`\varphi `$$``$$`\varphi _{\mathrm{trans}}`$ \+ 0.025).

We first tried to fit an orbital solution with the eccentricity as a free parameter. The obtained value was very low: $`e`$ = 0.015 $`\pm `$ 0.019. We then compared the residuals around this solution with the value obtained from a circular solution. We found the same residuals for both solutions. According to the Lucy & Sweeney (1971) test, the probability of significant eccentricity is zero for these data. Monte-Carlo simulations performed with the ORBIT software (Forveille et al. 1999) give the following 3 $`\sigma `$ uncertainty for this parameter: +0.063/$``$0.067. An eccentricity larger than 0.08 can thus be rejected with a very high confidence level.

We present in Table 3 the circular orbital solution fitted to the ELODIE+CORALIE radial velocities.

The velocities and orbital solution are displayed in Fig. 3. This Figure also shows the residuals to the fitted orbit. We note a spectacular improvement of the solution quality since the initial publication (Mazeh et al. 2000). It results from new template used for the cross-correlation (Pepe et al. 2002).

Wittenmyer et al. (2003) have recently presented updated high-precision period and transit time for HD 209458. Their determinations, based on space and ground-based photometric measurements of transits, are: $`P`$ = 3.5247542 $`\pm `$ 0.0000004 d and $`T_{\mathrm{trans}}`$ = 2 452 618.66774 $`\pm `$ 0.000007. Our period is compatible with their value (1.5 $`\sigma `$ away). Our transit time prediction obtained from the fitted time of maximum radial velocity ($`T_{\mathrm{trans}}`$ = 2 452 618.639 $`\pm `$ 0.021) is 1.3 $`\sigma `$ away from their value. We have computed an orbital solution with $`P`$ fixed to the Wittenmyer et al. (2003) value and $`T`$ fixed to the time of maximum velocity computed from their $`T_{\mathrm{trans}}`$. The obtained values for the three remaining free parameters are less than 1 $`\sigma `$ away from the values in Table 3.

With a mass of 1.15 M for HD 209458 (Santos et al. 2003), we compute the planetary minimum mass and semi-major axis: $`m_2\mathrm{sin}i`$ = 0.699 MJup and $`a`$ = 0.048 AU.

### 3.2 51 Peg (HD 217014, HIP 113357, GJ 882)

The detection of the first extra-solar planet around a Solar-type star, 51 Peg b, was published in Mayor & Queloz (1995). This detection was rapidly confirmed by Marcy et al. (1997). Since the discovery, we have gathered many ELODIE radial-velocity measurements of 51 Peg.

The number of measurements for this star now amounts to 153. The time span of these observations is 3278 days (from JD = 2 449 610 (September 1994) to JD = 2 452 888 (September 2003)) and their mean velocity uncertainty is 7.3 m s-1.

The updated orbital solution for 51 Peg is presented in Table 3. As for HD 209458, the Lucy & Sweeney (1971) test shows that a fitted non-zero eccentricity is not significant. We thus fitted a circular orbit. Figure 4 shows the phase-folded velocities and fitted orbit together with the residuals to the solution. The rms of the residuals for this solution is 11.8 m s-1, a value slightly smaller than the 13 m s-1 obtained by Mayor & Queloz (1995).

Assuming a mass of 1.04 M for 51 Peg (Santos et al. 2003), we compute the planetary minimum mass and semi-major axis: $`m_2\mathrm{sin}i`$ = 0.468 MJup and $`a`$ = 0.052 AU.

## 4 Confirmations of planetary companions

In the same way as the presence of the companion to 51 Peg was confirmed by Marcy et al. (1997), we decided to start a follow-up of some of the stars with planetary companions announced by other teams. Confirming and cross-checking the claimed detections is very important, especially when the signal is close to the detection limits.

In this section, we present our radial-velocity measurements for Ups And, 55 Cnc, 47 UMa, 70 Vir and HD 187123. These stars host one or several planetary companions detected by the California & Carnegie Planet Search team. We confirm the orbital solutions for Ups And, 70 Vir and HD 187123. We are not able to fully confirm the published solutions for 55 Cnc and 47 UMa because of insufficient data (too short measurement span or poorer precision).

### 4.1 Ups And (HD 9826, HIP 7513, GJ 61)

The presence of 4.61–d planet in orbit around Ups And was announced by Butler et al. (1997). Two additional longer-period planets were later uncovered (Butler et al. 1999) by the same team in collaboration with the Advanced Fiber-Optic Echelle spectrometer (AFOE, Brown et al. 1994) planet search team (Korzennik et al. 1998).

71 ELODIE radial-velocity measurements have been gathered between JD = 2 450 319 (August 1996) and JD = 2 452 886 (September 2003). The mean velocity uncertainty of these measurements is 8.3 m s-1. We present in Table 4 the 3-planet orbital solution derived from our data. We note some small differences between our solution and the Butler et al. (1999) solution as e.g. a slightly lower period value for Ups And c and a slightly longer period for Ups And d.

The ELODIE radial velocities and orbital solution are shown in Fig. 5. The residuals to the fitted orbit are clearly too large as it can be seen from the $`\chi _{\mathrm{red}}^2`$ value: 4.34. The residuals obtained with the orbital parameters fixed to the Butler et al. (1999) values are even worse: 22.4 m s-1 ($`\chi _{\mathrm{red}}^2`$ = 7.75) We have no explanation for these facts but we are nevertheless able to confirm the presence of the three planets around Ups And.

Adopting the same primary mass as in Butler et al. (1999), $`M_1`$ = 1.3 Msun, we compute the planetary minimum masses: $`m_\mathrm{b}\mathrm{sin}i`$ = 0.75 MJup; $`m_\mathrm{c}\mathrm{sin}i`$ = 2.25 MJup; $`m_\mathrm{d}\mathrm{sin}i`$ = 3.95 MJup. The semi-major axes are: $`a_\mathrm{b}`$ = 0.059 AU; $`a_\mathrm{c}`$ = 0.82 AU; $`a_\mathrm{d}`$ = 2.57 AU.

### 4.2 55 Cnc (HD 75732, HIP 43587, GJ 324 A)

55 Cnc b, a 14.6–d planet was announced by Butler et al. (1997). These authors also quoted the presence of an additional radial-velocity drift superimposed on the periodic signal. Marcy et al. (2002) have recently announced the detection of two additional signals with periods of 44.3 d and 14.7 yr. The interpretation of the 44.3–d signal is uncertain since it corresponds to the rotational period of the star. It thus could be due to stellar intrinsic effects. The 14.7–yr signal is interpreted as the result of the presence of a 4 MJup planet more than 5 AU away from the star.

48 ELODIE radial-velocity measurements of 55 Cnc have been obtained between JD = 2 450 207 (May 1996) and JD = 2 452 779 (May 2003). They have a mean uncertainty of 7.3 m s-1. Our measurement time span (2571 d) is shorter than the announced long period by a factor of 2. Also, the radial-velocity semi-amplitude of the 44.3–d signal is only 13 m s-1. These features make difficult the characterization of the outer planets with the ELODIE measurements but they however allow us to partially confirm the orbital solution presented in Marcy et al. (2002).

We present in Table 4 a 2-Keplerian orbital solution for 55 Cnc. This solution does not account for the intermediate signal. We are able to fit a 44–d orbit to the residuals around this solution but the fit parameters are too poorly constrained to be significant. Thus, they do not appear in this paper. Our data are insufficient for confirming the presence of the intermediate planet. Figure 6 shows the fitted orbits to our velocities as well as the residuals around our solution.

Assuming a primary mass of 0.95 M for 55 Cnc, we compute the minimum masses of the two companions: $`m_\mathrm{b}\mathrm{sin}i`$ = 0.91 MJup and $`m_\mathrm{d}\mathrm{sin}i`$ = 2.89 MJup. The orbital semi-major axes are: $`a_\mathrm{b}`$ = 0.115 AU and $`a_\mathrm{d}`$ = 5.28 AU.

### 4.3 47 UMa (HD 95128, HIP 53721, GJ 407)

The detection of a long-period planet in orbit around 47 UMa was announced in Butler & Marcy (1996). A preliminary ELODIE orbit for this companion was published in Naef et al. (2001). More recently, Fischer et al. (2002) have presented an updated orbital solution including another planetary companion on a longer-period orbit.

We have obtained 44 ELODIE radial-velocity measurements of 47 UMa between JD = 2 450 150 (March 1996) and JD = 2 452 752 (April 2003). The mean velocity uncertainty of these measurements is 7.3 m s-1. We are unable with our data to confirm the presence of the second planetary companion. The solution presented in Table 5 only includes one Keplerian orbit. The residuals to this solution (7.4 m $`s^1`$) are in good agreement with our uncertainties. We do not see any clear velocity trend in these residuals (Fig. 7). With the same primary mass as in Fischer et al. (2002), $`M_1`$ = 1.03 M, we compute the planet minimum mass and its orbital semi-major axis: $`m_2\mathrm{sin}i`$ = 2.76 MJup and $`a`$ = 2.11 AU.

We have tried without success to adjust a 2-Keplerian orbital solution to our velocities, starting from the solution published in Fischer et al. (2002). Moreover, the residuals we obtain by fixing the Fischer et al. (2002) solution to our data are large: $`\sigma _{\mathrm{O}\mathrm{C}}`$ = 13 m s-1. This solution is thus not compatible with the ELODIE data. In a recent paper, Ford (2003) has shown that the orbital solution for 47 UMa was almost unconstrained by the Fischer et al. (2002) velocities. For example, he has obtained a statistically equivalent solution with a period three times longer and an eccentricity eight times larger than the published values. At this point, the ELODIE data are not sufficient for confirming or ruling out the presence of another body in the system. They are well fitted by a single planet orbit. These facts demonstrate that the Fischer et al. (2002) orbit is very preliminary.

### 4.4 70 Vir (HD 117176, HIP 65721, GJ 512.1)

The detection of a massive planet on a moderately eccentric 116–d orbit around 70 Vir was announced in Marcy & Butler (1996). A preliminary ELODIE orbit confirming the Marcy & Butler (1996) solution was presented in Naef et al. (2001).

We now have 35 ELODIE radial-velocity measurements. They have been obtained between JD = 2 450 150 (March 1996) and JD 2 452 777 (May 2003). They have a mean velocity uncertainty of 7.1 m s-1. We present in Table 5 the orbital solution fitted to our velocities. This solution is in full agreement with the Marcy & Butler (1996) solution. It is also much more precise than our preliminary solution (Naef et al. 2001). The improvement is due to the use of the new cross-correlation template reducing the influence of the atmospheric telluric lines (Pepe et al. 2002). With residuals of 6.1 m s-1, it is to date the most precise published ELODIE orbital solution and it thus illustrates well the level of precision presently achieved with this instrument. Figure 8 shows our radial velocities and the fitted solution.

With the primary mass quoted in Marcy & Butler (1996), $`M_1`$ = 0.92 M, we find a minimum mass of 6.56 MJup for 70 Vir b. Its orbital semi-major axis is 0.456 AU.

### 4.5 HD 187123 (HIP 97336)

The presence of a 3.1–d planet in orbit around HD 187123 was announced by Butler et al. (1998). An updated orbital solution including an additional linear velocity drift was presented in Vogt et al. (2000).

We started our ELODIE measurements of HD 187123 shortly after the discovery announcement. We now have 57 radial-velocity measurements of this star spanning from JD 2 451 087 (September 1998) to JD 2 452 889 (September 2003). They have a mean velocity uncertainty of 9.6 m s-1. The orbital solution fitted to this velocities is presented in Table 5. We fitted a circular orbit since our data do not exhibit a significant eccentricity according to the Lucy & Sweeney (1971) test. Our solution also includes an additional linear velocity drift. The derived parameters are in full agreement with those presented in Vogt et al. (2000). We therefore confirm the presence of the short-period planet and the linear drift probably due to another massive body in the system. Figure 9 shows the ELODIE velocities and the fitted orbital solution.

Using the primary published in Santos et al. (2003), $`M_1`$ = 1.05 M, we derive $`m_2\mathrm{sin}i`$ = 0.52 MJup and $`a`$ = 0.042 AU for HD 187123 b.

## 5 Summary

In this paper, we have presented our ELODIE radial-velocity data and derived planetary solutions for HD 74156, 14 Her, HD 209458, 51 Peg, Ups And, 55 Cnc, 47 UMa, 70 Vir and HD 187123, including three new ELODIE planets, updated orbital solutions of previously published candidates and confirmations of detections from other planet search programmes.

The radial-velocity variations detected for HD 74156 are found to be due to the presence of a planetary system consisting of a 1.86 MJup planet at 0.294 AU and a 6.17 MJup planet on a 5.5–yr orbit. For 14 Her, a 4.74 MJup planet at 2.80 AU induces the main detected signal. An additional linear velocity drift revealing the presence of another not yet characterized massive body in the system is detected as well for 14 Her.

We have presented updated more precise orbits for HD 209458 and 51 Peg, two stars hosting short-period planets detected with ELODIE.

We have confirmed the previously published orbital solutions for Ups And, 70 Vir and HD 187123. We have also partially confirmed the solutions for 55 Cnc and 47 UMa. The velocity signals of the inner and the outer planet orbiting 55 Cnc are clearly detected. We have no evidence for the presence of 47 UMa c. Our data is not compatible with the Fischer et al. (2002) solution. However, the radial-velocity semi-amplitude induced by this object (11 m s-1) is small compared to our measurement precision (7.3 m s-1). We thus cannot rule-out this claimed detection.

###### Acknowledgements.

We acknowledge support from the Swiss National Research Found (FNRS), the Geneva University and the French CNRS. We are grateful to the Observatoire de Haute-Provence for the generous time allocation. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France.


format: text

text:
# 1 Cluster growth and corresponding derivative as function of the overall filling factor 𝜈 (a) and of the cluster filling factor 𝜂 (b).

BI-TP 2003/30

Predictions for $`J/\psi `$ Suppression

by Parton Percolation

S. Digala, S. Fortunatoa and H. Satza,b

a: Fakultät für Physik, Universität Bielefeld

Postfach 100 131, D-33501 Bielefeld, Germany

b: Centro de Fisíca das Interacções Fundamentais (CFIF)

Instituto Superior Técnico, Av. Rovisco Pais, P-1049-001 Lisboa, Portugal

Abstract:

Parton percolation provides geometric deconfinement in the pre-equilibrium stage of nuclear collisions. The resulting parton condensate can lead to charmonium suppression. We formulate a local percolation condition viable for non-uniform collision environments and show that it correctly reproduces the suppression observed for $`SU`$ and $`PbPb`$ collisions at the SPS. Using this formulation, we then determine the behavior of $`J/\psi `$ suppression for $`InIn`$ collisions at the SPS and for $`AuAu`$ collisions at RHIC.

1. Introduction

In recent years, the role of the partonic initial stages of high energy nuclear collisions has attracted increasing interest. While the ultimate aim of the experimental program is the production of the quark-gluon plasma predicted by statistical QCD, it is becoming more and more evident that a thermalized medium consisting of quarks and gluons can be produced only if already the initial state of the collision provides the conditions necessary for deconfinement and thermalization. Dilute initial state parton configurations do not lead to the colour connection needed so that partons from different collisions can combine to form a collective medium. On the other hand, if the initial state of primary partons is sufficiently dense, cluster percolation will set in, leading to a condensate of connected and hence interacting partons, which are no longer associated to any “parent” hadrons. Thus parton percolation is a geometric, pre-equilibrium form of deconfinement ; it must occur in any description based on partons of an intrinsic transverse momentum.

Subsequent interactions can thermalize this interconnected partonic system, turning it into a quark-gluon plasma, so that parton percolation constitutes an essential prerequisite for QGP formation. Moreover, pre-equilibrium deconfinement leads to a further important question. Much of the investigation of high energy nuclear collisions is devoted to the search for quark-gluon plasma signatures. However, some of the features observed in nuclear interactions might be determined by the system present before any equilibration has occurred, and they could be independent of a subsequent thermalization. We must therefore investigate observable consequences of parton percolation, and if possible, show how they can be distinguished from those due to QGP formation.

One of the features of particular interest in this connection is the suppression of $`J/\psi `$ production in nuclear collisions, predicted as deconfinement signal for thermal media quite some time ago . If the parton condensate formed through percolation contains partons hard enough to resolve the produced charmonium states, these partons can also dissociate charmonia, so that the onset of $`J/\psi `$ suppression could coincide with that of parton percolation. It was shown in a first study  that for $`PbPb`$ collisions at the CERN-SPS as function of centrality the resulting thresholds appear quite reasonable. A more quantitative comparison, which can then also provide the basis for further predictions, requires a formulation of percolation in a non-uniform environment, as given by the partonic source profile in nuclear collisions, and adapted to the local nature of charmonium probes. The aim of this paper is first to formulate percolation for such conditions, then to compare the results to existing SPS data ($`SU`$ and $`PbPb`$), and finally to present the predictions of this approach for the forthcoming $`J/\psi `$ production experiments at CERN and BNL. All results are obtained by Monte Carlo simulations taking into account the centrality dependence of the collisions. In closing, we shall consider some features which could distinguish $`J/\psi `$ suppression by parton percolation from a suppression in a thermal medium.

2. Local Parton Percolation Conditions

Consider a flat two-dimensional circular surface of radius $`R`$ (the transverse nuclear area), on which $`N`$ small discs of radius $`rR`$ (the transverse partonic size) are randomly distributed, allowing overlap. With increasing density $`nN/\pi R^2`$, clusters of increasing size appear. The crucial feature is that this cluster formation shows critical behaviour: in the limit $`N\mathrm{}`$ and $`R\mathrm{}`$, with $`n`$ finite, the cluster size diverges at a certain critical density $`n=n_c`$. This percolation threshold is given by

$$n_c=\frac{\nu _c}{\pi r^2};$$
(1)
the critical value $`\nu _c1.13`$ of the ‘filling factor’ $`\nu n(r/R)^2`$ is determined by numerical studies. For finite $`N`$ and $`R`$, percolation sets in when the largest cluster spans the entire circular surface from the center to the edge. The resulting cluster growth is illustrated in Fig. 1a for a ratio $`r/R=1/100`$; we show the percolation probability, defined as the relative size of the largest cluster, as function of the filling factor $`\nu `$. Because of overlap, a considerable fraction of the surface is still empty at the percolation point; in fact, at that threshold, only $`1\mathrm{exp}\{\nu _c\}2/3`$ of the surface is covered by discs.

In high energy nuclear collisions, the incident partons are distributed on a plane transverse to the collision axis. However, since they originate from the nucleons within the colliding nuclei, this distribution is highly non-uniform, with more nucleons and hence more partons in the center than towards the edge of the transverse nuclear plane. If the basic surface is not flat, it is still possible to define the percolation threshold as the point at which the cluster size shows singular behavior in the limit of large $`R`$ and $`N`$ . Such a global definition is, however, not always the most useful. Hard probes, such as quarkonia, probe the medium locally and thus test only if it has reached the percolation point and the resulting geometric deconfinement at their location; they cannot register the global features of the entire cluster area. It is thus necessary to define a more local percolation criterium, and this is in fact quite straight-forward.

As mentioned, at the percolation point, $`\mathrm{exp}\{\nu _c\}1/3`$ of the surface remains empty. Hence the disc density in the percolating cluster must be greater than $`(1.5\nu _c/\pi r^2)`$. Numerical studies show that in fact percolation sets in when the density $`m`$ of constituents in the largest cluster reaches the critical value

$$m_c=\frac{\eta _c}{\pi r^2},$$
(2)
with $`\eta _c1.72`$, slightly larger than $`1.5\nu _c`$. This result provides the required local test: if the parton density at a certain point in the transverse nuclear collision plane has reached this level, the medium there belongs to a percolating cluster and hence to a deconfined parton condensate. In Fig. 1b, we illustrate how the percolating cluster grows as function of the cluster filling factor $`\eta =m(r/R)^2`$, again with $`r/R=1/100`$.

Let us recapitulate: there are two equivalent criteria which can be used to specify the onset of percolation. In the ‘global’ definition, it occurs when the total number of partons distributed over the entire transverse nuclear area reaches $`1.13/\pi r^2`$. In the ‘local’ definition, percolation sets in when the parton density in the largest cluster reaches $`1.72/\pi r^2`$. Which of the two is preferable in a given case depends on the physics question addressed, and if we want to study local phenomena such as $`J/\psi `$ suppression, it is clearly the local parton density which is relevant, not the average over the entire transverse area.

We now turn to the implications of the approach to nuclear collisions ; for illustration, we concentrate for the moment on central $`AA`$ interactions at a c.m.s. energy $`\sqrt{s}`$ per nucleon-nucleon collision. The distribution of nucleons in the colliding nuclei is specified by a Glauber calculation using Woods-Saxon nuclear distributions ; this provides the density $`n_s(A)`$ of nucleons in the transverse collision plane. The parton content of a nucleon is given by parton distribution functions $`dN_q(x,Q^2)/dy`$ determined in deep inelastic scattering experiments; here $`x`$ denotes the fraction of the nucleon momentum carried by the parton and $`Q`$ the momentum resolution scale. At central c.m.s. rapidity $`y=0`$, we have $`x=k_T/\sqrt{s}`$, where $`k_T`$ denotes the average transverse momentum of the parton. In nuclear collisions, $`k_T`$ defines the transverse size of the partons and thus also sets the resolution scale, $`k_TQ`$. Using these quantities, we have for the ‘global’ percolation condition in nuclear collisions

$$n_s(A)\left(\frac{dN_q(x,Q_c^2)}{dy}\right)_{x=Q_c/\sqrt{s}}=\frac{\nu _c}{(\pi /Q_c^2)};$$
(3)
it determines for what value of $`A`$ at a given $`\sqrt{s}`$ percolation sets in, and it specifies the value $`Q=Q_c(A,\sqrt{s})`$ at the onset point. Using the nuclear source density $`n_s(A)`$ together with the parton distribution function, we can also calculate the density of partons $`m(A,Q,\sqrt{s})`$ for the largest cluster in the transverse collision plane. The relation

$$m_c(A,Q_c,\sqrt{s})=\frac{\eta _c}{(\pi /Q_c^2)};$$
(4)
then provides the ‘local’ parton percolation condition to specify the onset values of $`A,\sqrt{s}`$ and $`Q`$. These considerations can be extended to non-central collisions in a straight-forward fashion .

It is not a priori evident that the threshold values $`\nu _c`$ and $`\eta _c`$ for the non-uniform distribution specified by the Woods-Saxon distributions remain the same as for a uniform distribution, nor is it clear how large finite size corrections are at given values of $`r/R`$. We have therefore calculated the pseudo-critical values (derivative peak positions) of $`\nu `$ and $`\eta `$ for $`R/r=20,50`$ and $`100`$, using the Woods-Saxon distribution for central $`PbPb`$ collisions to determine the parton density as well as the average overall interaction region. Finite-size scaling techniques from statistical physics then provide the asymptotic limit; in this way, we have verified that $`\nu _c`$ and $`\eta _c`$ for this non-uniform distribution in fact agree with the above mentioned values for the uniform case. Moreover, we found that for the mentioned values of $`R/r`$, the deviations from the asymptotic critical values are less than 2 %.

2. Comparison to Pb-Pb and S-U Data

We now want to study parton percolation for CERN-SPS conditions and compare the resulting consequences for $`J/\psi `$ suppression to the data of the NA38/NA50 collaborations.

The parton content of a nucleon is obtained from GRV94 parton distribution functions , since these are available for the required kinematic range; for more details, see . With this specified, we calculate the size of the largest cluster as well as the overall size of the interaction region for different centralities (the latter is here defined by the ratio of the number to the density of wounded nucleons). For $`PbPb`$ collisions at $`\sqrt{s}=17.4`$ GeV, we show in Fig. 2a the relative size of the largest cluster as function of the filling factor $`\eta `$. It is seen that the cluster size grows rapidly in a narrow band around $`\eta _c`$, and at the percolation onset, the percolating cluster covers about 1/3 of the total collision area. In Fig. 2b, we show the cluster filling factor $`\eta `$ as function of centrality, measured by the number of participating nucleons. The critical value for parton percolation, $`\eta _c=1.72`$, is found to be attained for $`N_{\mathrm{part}}125`$.

The parton distribution functions specify both the number of partons and the resolution scale $`Q_c`$ at the percolation point. For $`PbPb`$ collisions at $`\sqrt{s}=17.4`$ GeV, $`Q_c0.7`$ GeV. The scales of the charmonium states $`\chi _c`$ and $`\psi ^{}`$ , as determined by the inverse of their radii calculated in potential theory, are around 0.6 GeV and 0.5 GeV, respectively. Since the parton condensate can thus resolve these states, we assume, following , that it will dissociate them. To be specific: we assume that all $`\chi _c`$ and $`\psi ^{}`$ states formed inside the percolating cluster disappear; the location of their formation is determined by the collision density as given by the Glauber formulation . The first onset of $`J/\psi `$ suppression in $`PbPb`$ collisions at the SPS should therefore occur at $`N_{\mathrm{part}}125`$, where the $`J/\psi `$’s due to feed-down from $`\chi _c`$ and $`\psi ^{}`$ states in the percolating cluster are eliminated. Directly produced $`J/\psi `$’s survive because of their smaller radii (leading to a scale of 0.9 - 1.0 GeV), as do those coming from excited states outside the percolating cluster. We emphasize that with parton distribution functions and charmonium radii given, the description contains no adjustable parameters.

As noted, the parton condensate formed at the onset of percolation has a resolution scale $`Q_c0.7`$ GeV and thus cannot resolve the ground state charmonium $`J/\psi `$ with its smaller radius. The dissociation of directly produced $`J/\psi `$’s thus requires more central collisions, which lead to a higher parton density and simultaneously to a better resolution, i.e., to an increase of $`Q_c`$. The precise onset point is here more difficult to determine, since it is fixed by the value of the $`J/\psi `$ radius. For $`r_{J/\psi }=0.20`$ fm, we have $`Q1.0`$ GeV and obtain the centrality dependence shown in Fig. 3, indicating that the suppression of directly produced $`J/\psi `$’s starts at $`N_{\mathrm{part}}320`$. However, here a note of warning has to be added. The $`\chi _c`$ and $`\psi ^{}`$ states can be resolved by the parton condensate when it is first formed; the values of their intrinsic scales do not enter explicitly, since they are below the resolution scale $`Q_c`$ of the condensate. The onset point for the suppression of the $`1S`$ ground state is directly determined by the radius of that state, and since the dependence of parton density on centrality (see Fig. 3) is quite flat, small changes in $`r_{J/\psi }`$ lead to large changes in the threshold value. For instance, $`r_{J/\psi }=0.22`$ fm ($`Q0.9`$ GeV) shifts the onset from $`N_{\mathrm{part}}320`$ to about 190, as also shown in Fig. 3.

We now combine both thresholds and use the nuclear collision profile to determine the distribution of charmonium production in the transverse plane , in order to determine which are formed in the cluster and which are not. As function of centrality, this leads to the survival pattern shown in Fig. 4.

It is clear that various theoretical and experimental sources (fluctuations, resolution, binning, etc.) will lead to some smearing of the theoretical survival pattern. To obtain some idea of what this may lead to, we have folded our calculated distribution $`S(N_{\mathrm{part}})`$ with a Gaussian smearing function,

$$S(N_{\mathrm{part}})=𝑑\overline{N}_{\mathrm{part}}S(\overline{N}_{\mathrm{part}})\frac{1}{\sqrt{2\pi }\sigma }\mathrm{exp}\{(N_{\mathrm{part}}\overline{N}_{\mathrm{part}})^2/2\sigma ^2\},$$
(5)
using a relative resolution of $`\sigma /N_{\mathrm{part}}`$ = 15 % at the respective thresholds. The result is also shown in Fig. 4.

The centrality dependence of the $`J/\psi `$ survival has been measured in Pb-Pb collisions at the SPS using two independent variables, the associated transverse energy, $`E_T`$, and the energy of the undeflected nucleons, $`E_{\mathrm{ZDC}}`$. The first effectively measures the number of hadronic secondaries in the transverse plane, the second the number of participating nucleons. We want to study anomalous charmonium suppression, i.e., effects above and beyond the normal suppression already observed in p-A collisions, due to absorption of the nascent charmonium state in normal nuclear matter . To determine this normal suppression, we use a recent analysis  of $`pA`$ and $`SU`$ data by the NA50 group, which finds a somewhat weaker nuclear absorption, $`\sigma _d4.3\pm 0.6`$, than previous studies. Normalizing the Pb-Pb data to the normal nuclear suppression thus specified leads to the points shown in Fig. 5, where they are compared to our parton percolation suppression pattern. We see good agreement between the calculated curve and the data points. We recall, however, that the second threshold value depends crucially on the $`J/\psi `$ radius.

Next we want to check the behaviour predicted for $`SU`$ collisions at the SPS, since this shows at best only very slight anomalous suppression111Initially, it was claimed that the data are fully accounted for by normal nuclear absorption. However, the lower absorption cross sections obtained in the mentioned new p-A studies allow a slight anomalous suppression for the most central $`SU`$ collisions.. In Fig. 7, we show the percolation behaviour, and in Fig. 7, the resulting survival pattern together with the data . The threshold obtained here corresponds to the general onset of percolation; the resolution scale of the direct $`J/\psi `$ cannot be reached in $`SU`$ collisions. The corresponding smeared form of the calculated distribution is also included.

We close this section with some general comments on the suppression of charmonium production. We have here considered parton percolation as mechanism for anomalous $`J/\psi `$ suppression, defined as the reduction of the observed production rate beyond that predicted by normal absorption in standard nuclear matter. In doing so, we neglect other conceivable suppression mechanisms, in particular dissociation by (non-percolating) secondary hadronic or partonic comovers. This is not expected to play a significant role for the suppression of $`1S`$ states, i.e., directly produced $`J/\psi `$’s: the large binding energy (650 MeV) makes a dissociation by confined partons very unlikely . However, since about 40 % of the observed $`J/\psi `$’s come from the decay of higher excited states (feed-down from $`\chi _c`$ and $`\psi ^{}`$ production), we have to consider the fate of these as well. The $`\chi _c`$ ($``$ 250 MeV) is still rather strongly bound, but the much smaller $`\psi ^{}`$ binding energy ($``$ 60 MeV) means that any interaction will cause break-up. It is thus not surprising that the $`\psi ^{}`$ suffer an anomalous suppression even in the most peripheral S-U interactions . This is most likely not related to any parton deconfinement and has to be addressed separately.

3. Predictions for NA60 and RHIC

We now turn to the forthcoming $`InIn`$ collisions at the SPS. The percolation pattern is shown in Fig. 9; the general percolation threshold is reached at $`N_{\mathrm{part}}140`$. In contrast to $`PbPb`$ collisions, the threshold for directly produced $`J/\psi `$’s here remains above the produced density even for the most central collisions (but recall the caveat concerning its dependence on the $`J/\psi `$ radius). The resulting survival pattern for $`InIn`$ collisions at the SPS is given in Fig. 9, including also the smeared form. We see that the most central 30 % of the collisions show a suppression of the $`J/\psi `$’s produced by feed-down from higher excited charmonium states.

Finally we address the suppression pattern for $`AuAu`$ collisions at RHIC, with the percolation pattern shown in Fig. 11. Here the increased parton density shifts the onset of percolation to a higher resolution scale, so that from the threshold on, all charmonium states are suppressed. This leads to a single step suppression pattern, starting at $`N_{\mathrm{part}}90`$; it is given in Fig. 11. Recall that we are here specifying the predicted anomalous suppression, beyond the normal suppression resulting from (pre-resonance) absorption in standard nuclear matter. For various reasons, this normal suppression may depend on the collision energy and can thus be quite different from that observed at SPS energy. To determine the anomalous suppression at RHIC, it is therefore crucial to first measure the normal absorption in $`pA`$ (or $`dA`$) collisions. The threshold for $`AuAu`$ collisions as obtained here occurs for considerably more central collisions than found in previous studies , in which the dependence of parton distribution functions on the resolution scale had not taken into account.

In summary: if the origin of the two-step suppression pattern observed in the $`PbPb`$ SPS collisions is indeed the onset of parton percolation at different resolution scales, as discussed, then we expect a single step pattern both in $`InIn`$ collisions at the SPS and in $`AuAu`$ collisions at RHIC. In the former, the step is due to the suppression of the higher excited states only, with directly produced $`J/\psi `$’s surviving; in the latter, there is a common threshold for all charmonium states and hence only one step.

4. Discussion

We have presented the $`J/\psi `$ suppression pattern as obtained from parton percolation under the assumption of complete charmonium dissociation within a percolating cluster of sufficient resolution scale, full survival in media below percolation density or of insufficient resolution scale. In this closing section, we want to summarize the uncertainties of such an approach and compare its outcome to that of a thermal description.

The onset of parton percolation is specified directly through the parton distribution functions (PDF’s) as obtained from an analysis of deep inelastic lepton-hadron data. Uncertainties in the PDF’s produce some uncertainty in the precise values of the percolation threshold, but not in the onset as such. As mentioned, percolation must occur in any approach based on partons with an intrinsic transverse momentum. The assumption tested by comparing our results to $`J/\psi `$ data is that for charmonium suppression in nuclear collisions, the primary parton density is the crucial variable and the formation of a percolating cluster the relevant mechanism.

As noted, the onset of the suppression of directly produced $`J/\psi `$’s at SPS energy is quite sensitive to the choice of charmonium parameters (and on that of the parton distribution functions, for which slight changes would also produce considerable threshold shifts). On the other hand, the onset of the first step in SPS $`PbPb`$ pattern, the predicted single step in SPS $`InIn`$ collisions, and the predicted single step in RHIC $`AuAu`$ collisions should be quite robust. This should help in identifying the origin of the second step in SPS $`PbPb`$ collisions. It has been suggested  that this drop is largely due to transverse energy fluctuations for the most central collisions; this would imply that it should be absent in zero-degree calorimeter studies, where (see Fig. 5) it seems to still occur to some extent . A further test should come from the SPS $`InIn`$ and the RHIC $`AuAu`$ data, where the existence of a second step could not be attributed to parton percolation.

Finally we turn briefly to the question of the distinguishing the initial state charmonium suppression presented here from a suppression occurring at a later thermalized stage. We had here assumed that the primary parton configuration determines the pattern of $`J/\psi `$ suppression. What would one expect if the initial parton configuration leads to subsequent thermalization, with the formation of a more or less equilibrated quark-gluon plasma, and it is this QGP which causes the observed $`J/\psi `$ suppression? The most recent (quenched) lattice QCD studies  indicate that in a hot thermal medium of quarks and gluons, higher excited charmonium states are dissociated at or below $`T_c`$. In contrast, the ground state $`J/\psi `$ survives to considerably higher temperatures; present calculations still show a clear $`J/\psi `$ signal at least up to $`T=1.5T_c`$. This would imply that thermal dissociation removes the feed-down contributions from $`\psi ^{}`$ and $`\chi _c`$ states for energy densities around or below 1 GeV/fm3, while energy densities of more than 6–8 GeV/fm3 are required for the melting of directly produced $`J/\psi `$’s. It seems difficult to accommodate the existing SPS data with these thresholds, while the pattern resulting from parton percolation can do so quite naturally. For $`AuAu`$ collisions at RHIC, thermal melting would suppress the higher state feed-down contributions for practically all centralities, while the directly produced $`J/\psi `$’s disappear at best for the most central collisions. We thus conclude that the charmonium suppression thresholds from initial state and from thermal suppression are quite different also at much higher energy.

Acknowledgements

It is a pleasure to thank F. Karsch, E. Laermann, C. Lourenço and P. Petreczky for stimulating remarks and discussions.


text:
# Inferences from Bose-Einstein correlations in multiple particle production processes Supported in part by the KBN grant 2P03B 093 22

## 1 Introduction

An important reason for studying Bose-Einstein correlations (BEC) in multiple particle production processes is that they yield information about the interaction regions i.e. about the regions, where the hadrons are produced. There is much to learn about such regions. Let us consider two examples. In high energy $`e^+e^{}`$ annihilations the first stages of the process can be described perturbatively. The two primary leptons merge into a heavy photon, or a $`Z^0`$, this particle decays into a $`q\overline{q}`$ pair, the quark and the antiquark radiate gluons, the gluons split into $`q\overline{q}`$ pairs or into gluons, and so on. Then, however, something much more complicated happens: the swarm of parton gets converted into hadrons. This is a nonperturbative process, known as hadronization, which is poorly understood. One would like to know the size and shape of the region where the hadrons are produced, the time interval between the first annihilation and the end of the hadronization process, as well as the time interval between the beginning and the end of hadronization. As a more complicated example let us consider central heavy ion collisions at high energy. In this case the initial color fields are so strong that perturbative methods may be unreliable. It is usually assumed that the initial partons, whatever their production mechanism, rapidly thermalize. From that time on they can be considered as a fluid (a kind of liquid or a kind of gas), which expands according to the laws of hydrodynamics or according to some simplified version of the Boltzmann equation. Finally the hadrons are produced in a hadronization process which may be different from that in $`e^+e^{}`$ annihilations. Here, besides the hadronization, one is interested in the flows of the fluid, in its equation of state and its phase transitions, if any. The study of Bose-Einstein correlations among the final hadrons supplies tentative answers to all such questions. The problem is, however, how reliable these answers are? Till a few years ago the mood was optimistic (cf. e.g. ). The recent results from RHIC, and also some results from LEP, however, have been so unexpected and so puzzling (references can be traced e.g. from the recent review ) that perhaps some important ideas are still missing. In the present paper we review some of the most popular ways of getting from the experimental data to the physical conclusions. It will be seen that much remains to be clarified there.

The key object in the transition from the experimental data to the physical conclusions is usually a single particle density matrix $`\rho (\text{p}_1,\text{p}_2)`$. There are approaches where this is not the case, e.g. the string model (cf.  and references quoted there), but most models can be formulated in terms of this matrix . In this section we will present its relation to the data. The discussion is grossly simplified (cf. e.g. ), but it contains the main ideas. The relation between the density matrix and the features of the interaction region is the main subject of this paper and is discussed in the following sections. Let us stress that the density matrix $`\rho (\text{p}_1,\text{p}_2)`$ is an auxiliary construct, which does not have to coincide with the actual single particle density matrix for the mesons of a given type in the final state.

The diagonal elements of the density matrix are related, as usual, to the single particle momentum distribution

$$\frac{dN}{d^3p}\rho (\text{p},\text{p}).$$
(1)
The tilde in most models means equality up to a constant factor. Only in the GGLP model it is a more complicated operation. $`d^3p`$ may denote the infinitesimal volume in momentum space, or the covariant infinitesimal volume i.e. the volume in momentum space divided by the corresponding energy $`E(\text{p})`$. The two-particle distribution for identical bosons is

$$\frac{dN}{d^3p_1d^3p_2}\rho (\text{p}_1,\text{p}_1)\rho (\text{p}_2,\text{p}_2)+|\rho (\text{p}_1,\text{p}_2)|^2.$$
(2)
For uncorrelated, distinguishable particles the second term on the right-hand side would be absent. It results from the symmetrization of the product of the single particle density matrices of the two particles and thus reflects the BEC. Usually one considers the correlation function obtained by dividing the two-particle distribution by the product of the corresponding single particle distributions

$$C(\text{p}_1,\text{p}_2)1+\frac{|\rho (\text{p}_1,\text{p}_2)|^2}{\rho (\text{p}_1,\text{p}_1)\rho (\text{p}_2,\text{p}_2)}.$$
(3)
In practice the procedure is usually much more complicated than presented here, but the strategy is as described: the data are used to get the correlation function which is simply related to the single particle density matrix.

## 2 GGLP approach

The first model of BEC in multiple particle production was proposed by the Goldhabers Lee and Pais  (GGLP). As a realistic description of BEC in multiple particle production processes this model is outdated, but it is still the best introduction to the subject. GGLP introduced a density $`\rho (\text{x})`$ of particle sources. They considered identical pions and for definiteness also we will call the identical bosons pions, though the analysis can be applied to any identical bosons (or even to identical fermions, when the sign between the two terms in formula (3) is changed from plus to minus). GGLP made the following two assumptions: firstly, that all the pions are produced simultaneously and instantaneously; secondly, that the production is completely incoherent – the pions produced in two different space point do not interfere. Under these assumptions the single particle density operator is

$$\widehat{\rho }=d^3x|\text{x}\rho (\text{x})\text{x}|.$$
(4)
The corresponding density matrix reads

$$\rho (\text{p}_1,\text{p}_2)=\text{p}_1|\widehat{\rho }|\text{p}_2d^3xe^{i\text{qx}}\rho (\text{x}),$$
(5)
where $`\text{q}=\text{p}_1\text{p}_2`$. If this were the density matrix which can be obtained from formula (3), it would be a very nice result. The density of sources $`\rho (\text{x})`$ could then be obtained unambiguously just by inverting the Fourier transform. E.g. for

$$\rho (p_1,p_2)e^{\frac{1}{2}R^2\text{q}^2},$$
(6)
where $`R^2>0`$ is a constant, one would obtain

$$\rho (\text{x})e^{\frac{\text{x}^2}{2R^2}}.$$
(7)
Unfortunately, the density matrix (5) is untenable as a density matrix proportional to the true single particle density matrix, because it gives a single particle momentum distribution constant in all momentum space. GGLP introduced, therefore, a projection on the states allowed by energy-momentum conservation, i.e. they used the density matrix as calculated here to find the momentum distribution for all the particles present in the final state, multiplied it by the delta function of energy-momentum conservation and integrated over the (covariant) momentum space of all the particles except two identical pions for the two particle distribution and over all the momenta except one for the single particle distribution. They got results in qualitative agreement with experiment. The method, however, was cumbersome – no one did calculations for more than six particles in the final state – and was in violent disagreement with experiment at the quantitative level .

## 3 Kopylov and Podgoretskii model

Very significant progress was obtained by Kopylov and Podgoretskii (cf.  and references contained there). In their model the density operator in the Schrödinger picture is

$$\widehat{\rho }=d^4x_se^{iH_0(t_st)}|\psi _s\rho (x_s)\psi _s|e^{H_0(t_st)},$$
(8)
where

$$\text{x}|\psi _s=\psi (\text{x}x_s);$$
(9)
or equivalently

$$\text{p}|\psi _sd^3xe^{i\text{px}}\psi (xx_s)=e^{i\text{px}_s}A(\text{p}),$$
(10)
where

$$A(\text{p})=d^3xe^{i\text{px}}\psi (\text{x}).$$
(11)
In these formulae $`x_s`$ is a point in space-time, as well as a label which defines unambiguously a source. It should be related to the position of the corresponding source in space-time, but the actual relation is a matter of choice. E.g. it could be the point where the source got created, or its mean position in space-time. Each source is labelled by a space-time point, but the particle produced by a source is not localized in a point. Its wave function in space is $`\psi (\text{x}\text{x}_s)`$ with an additional time dependent phase which is zero at $`t=t_s`$. All these functions are related by shifts in space-time. The wave function in momentum space is proportional to $`A(\text{p})`$, thus all the sources yields particles with the same momentum distribution. $`H_0`$ is the free particle Hamiltonian and we are interested in times $`t`$ larger than the latest time $`t_s`$, thus the particle evolves freely.

The corresponding density matrix,

$$\rho (\text{p}_1,\text{p}_2)A(\text{p}_1)A^{}(\text{p}_2)d^4x_xe^{iqs}\rho (x_s),$$
(12)
yields the single particle distribution,

$$\frac{dN}{d^3p}|A(\text{p})|^2,$$
(13)
which can be made to agree with any experimental distribution by a suitable choice of $`A(\text{p})`$. Thus, there is no obvious need to introduce the energy-momentum conservation constraint. Kopylov and Podgoretskii assumed that for final states with not too few particles in the interesting momentum region the energy-momentum conservation does not affect significantly the one- and two-body distributions, and that consequently their matrix (12) can be assumed equal, up to a normalizing factor, to the actual single particle density matrix. This assumption makes it easy to get results for high multiplicity exclusive channels and for high energy inclusive processes. Both were beyond the reach of the GGLP method. Besides the satisfactory formula for the single particle distribution given above, the model gives for the correlation function

$$C(\text{p}_1,\text{p}_2)=1+d^4x_se^{iqx_s}\rho (x_s).$$
(14)
Here the factors $`A(p)`$ cancel and the result is as for point sources. Note, however, that because of the four-fold integration this relation cannot be inverted to give $`\rho (x_s)`$, when the correlation function is known. Moreover, a priori the correlation function could depend on the vector q and on the vector $`\text{K}=(\text{p}_1+\text{p}_2)/2`$, while the right hand side depends only on the four-vector $`q`$. The reason is that the factor $`A(\text{p})`$ is the same for every source and cancels. Physically this means that there is no correlation between the momentum of the pion and the position of the source. Since such correlation follow from almost every model and since experimentally it is not true, that the correlation function depends only on $`q`$ (cf. e.g.  and references quoted there) this is a serious weakness of the model.

## 4 Yano and Koonin method

Another attempt to go beyond the GGLP approximations is due to Yano and Koonin . In our notation their key formula is

$$|\rho (\text{p}_1,\text{p}_2)|^2=Red^4x_1D(x_1,\text{p}_1)e^{iqx_1}d^4x_2D(x_2,\text{p}_2)e^{iqx_2}.$$
(15)
Putting

$$D(x,\text{p})=\delta (t)\rho (\text{x})$$
(16)
one recovers the GGLP model. Putting

$$D(x,\text{p})=|A(\text{p})|^2\rho (x),$$
(17)
which corresponds to the assumptions made in , one reproduces the results of Kopylov and Podgoretskii. It is not clear, whether this approach can be given sense in the framework of quantum mechanics, when there are position-momentum correlations, i.e. when the $`D`$ function does not factorize into a momentum dependent and a space dependent factor.

## 5 Covariant current formalism

The covariant current formalism (cf.  and references given there) can be derived from the model of Kopylov and Podgoretskii introducing the following two generalizations. The label characterizing the source is changed from $`x_s`$ to $`x_s,p_s`$. Thus, the source is characterized by its position in space-time and by its four momentum. Correspondingly, the density of sources $`\rho (x_s)`$ gets replaced by $`\rho (x_s,p_s)`$. Such labels are not subject to the Heisenberg uncertainty principle. E.g. a one-dimensional harmonic oscillator may have $`x=0`$ and $`p=0`$. The universal, momentum dependent function $`A(\text{p})`$ gets replaced by the source dependent $`j(\frac{p_sp}{m_s})`$, where $`m_s`$ is the mass and $`p_s`$ the momentum of the source. In the Kopylov Podgoretskii model all the sources had the same momentum distribution in some overall reference frame, e.g. in the center of mass frame of the collision. In the covariant current formalism each source has the same momentum distribution, when considered in its rest frame. Another way of formulating this model is to replace the incoherent sources by incoherent wave packets. Then $`p_s`$ and $`x_s`$ characterize the wave packet, and since the waves in the packet describe pions, it is natural to put $`m_s=m_\pi `$. In this model, assuming a fixed mass for all the sources,

$$\rho (p_1,p_2)=d^4x_sd^3p_s\rho (x_s,p_s)e^{iqs_s}j(\frac{p_sp_1}{m_s})j^{}(\frac{p_sp_2}{m_s}).$$
(18)
From the point of view of model builders, a nice feature of this approach is that one can assume a classical motion of the source, e.g. $`\text{x}_s=\text{x}_s(t_s)`$ and $`\text{p}_s=\text{p}_s(t_s)`$, and still have a formula which is consistent with quantum mechanics. In the covariant current formalism position-momentum correlations are naturally included.

## 6 Emission function method

The method of emission functions is very popular nowadays (cf. e.g. ). The emission function is built by analogy with the Wigner function , . The Wigner function $`W`$ is related to the single particle density matrix in the momentum representation by the formula

$$\rho (\text{p}_1,\text{p}_2)=d^3Xe^{i\text{qX}}W(\text{X},\text{K}),$$
(19)
Where $`\text{X}=(\text{x}_1+\text{x}_2)/2`$, and as usual $`\text{q}=\text{p}_1\text{p}_2`$ and $`\text{K}=(\text{p}_1+\text{p}_2)/2`$. The Wigner function is well-defined as a Fourier transform of the density matrix in momentum representation. On the other hand, as a function of X and K it may play the role of a phase space distribution. When studying the interactions regions, it is important to combine information about the space and the momentum distributions. Kopylov, Podgoretskii, and followers proposed to use the position of the source and the momentum of the particle. These two vectors can be measured simultaneously and their joint distribution gives an idea about the phase space distribution of particles. The Wigner function gives another ersatz phase space distribution. The trouble is, however, that the Wigner function refers to a given moment of time, while the hadrons are produced during a time interval. The emission function $`S`$ is supposed to improve on that and is related to density matrix by the formula

$$\rho (\text{p}_1,\text{p}_2)=d^4Xe^{iqX}S(X,K).$$
(20)
The usual strategy is to find from some classical or quasi-classical argument a phase space distribution as function of time, and then to interpret it as the emission function. Once this is done, one can, by a four-fold integration, obtain the density matrix and compare it with experiment. This comparison may be used to fix the free parameters of the model.

Let us note some difficulties of this approach. Given a density matrix there is an infinity of very different emission functions, which can reproduce it. For instance perfect agreement is obtained for

$$S(X,K)=\delta (X_0)W(\text{X},\text{K}).$$
(21)
This, however, corresponds to simultaneous and instant production of all the hadrons, which is not a very plausible scenario. When particle production is an incoherent sum of production amplitudes at various moments of time, the emission function $`S(\text{X},t,K)`$ can be related to the Wigner function of the particles produced at time $`t`$. In the general case, however, the relation of the emission function to a Wigner function is hardly visible . Thus, the formula can be used to eliminate wrong models rather than to prove that a model is implied by the data. This is not necessarily bad. Eliminating the unacceptable values of the parameters of a model, one learns which are the good ones. Nevertheless, one must always keep in mind that a completely different model can give equally good, or better results. One should also keep in mind that X and K are not really position and momentum, but only half sums of the corresponding arguments of the density matrix. When the contributions of the sources are strongly smeared in coordinate space or in momentum space, the difference may be very significant.


text:
# In flight calibration of the ISGRI camera

## 1 Introduction

The ISGRI low energy camera (15 keV-1 MeV) of the IBIS imager (Ubertini et al. (2003)) on-board INTEGRAL (Winkler et al. (2003)) is an array of 128 by 128 Cadmium Telluride (CdTe) pixels $`4\times 4\times 2`$ $`\mathrm{mm}^3`$ in size. They are grouped in 8 Modular Detector Units (MDU) which are assembled on an aluminium structure. For more details on the design and expected performances, see Lebrun et al. (2003). CdTe has a large photoelectric cross section, due to its high Z. However, ballistic losses cause an incomplete collection of the charges created during the gamma-ray interaction and thus produce a large continuum up to the photopeak energy. Nevertheless since the charge loss is proportional to the transit time which governs the pulse rise time, the deposited energy can be corrected using a measurement of rise time. For this reason, the ISGRI ASIC has been designed to read both the pulse height and rise time of 4 pixels, and a charge correction algorithm has been developed to provide good energy resolution over the whole energy range (Lebrun et al. (1996)). The values obtained on ground vary from 8% at 60 keV to 4% at 511 keV.

The in-flight calibrations have been performed using scientific data from Cygnus X-1, Crab, Galactic Center and empty field observations, as well as calibration unit data (in coincidence with the $`{}_{}{}^{22}\mathrm{Na}`$ source on board also called S2 data, see Ubertini et al. (2003) and Bird et al. (2003)). The aim is to:

* provide a good algorithm to detect and suppress bad and noisy pixels
* update the pulse height and rise time gains and offsets determined on ground using in-flight background and calibration source lines
* check the validity and update the charge loss correction algorithm
* determine the background structure and provide a correction for deconvolution

Calibration of the ISGRI point spread function (PSF) (Gros et al. (2003)) as well as the instrument response matrix modeling and in-flight validation (Laurent et al. (2003)) are presented elsewhere.

## 2 Pixel behaviour

When ISGRI was first switched on, the pixels were randomly bursting with events having a null rise time. The origin of these events is still unclear. They represent around 15% of the total ISGRI count rate, and since they are useless, a patch was applied to the on-board software to reject low rise time events. Moreover, noisy pixels are cut on-board by the NPHS (Noisy Pixels Handling System), which detects high number of triggers in a single pixel and in a module. Because of the zero rise time events, it was necessary to lower the NPHS sensitivity in order to limit the number of pixels switched off. For more details on NPHS and the zero rise time events, see also Lebrun et al. (2003).

To limit the number of noisy pixels, a ground based algorithm has been developed to update the thresholds of each pixel depending on its individual spectrum after every revolution. Noisy pixels have a strong low energy peak but have a correct high energy spectrum. A new threshold is computed in order to remove the whole peak. Bad pixels usually have a low efficiency and a very hard spectrum with a large bump between 50 and few hundred keV. They are selected using an hardness ratio and switched off. All the pixels, identified as noisy in a previous revolution, have their threshold lowered by one step ($``$ 1 keV).

The introduction of the aforementioned algorithm has strongly increased the number of pixels switched off and is now steady at around 420, i.e. less than 3% of the camera. In the same time the number of low efficiency pixels has decreased, as can be seen on figure 1. The number of noisy pixels found in each revolution is $``$250. The total fraction of pixels switched off or noisy is therefore $``$4.5%, which is better than the 5% specification.

Some other unexpected features of the CdTe pixels response in flight can be seen on the background biparametric diagram on figure 5. For instance, a large horizontal structure at high rise time can be seen spanning over the whole energy range. Its origin is still unclear, though likely due to cosmic rays since it was barely visible on ground tests. It does not strongly limit the sensitivity since there are not many events at such high rise times. An on board upper limit has been applied to remove those background events and lower the telemetry load.

## 3 Energy calibration

The gains and offsets of rise time and pulse height have been calibrated on the ground at various biases and temperatures using several sources. A look-up table (LUT1) valid at $`0^{}\mathrm{C}`$ and 100V was produced, and a dependency on temperature and bias was found. This is detailed in the next section and in Lebrun et al. (2003).

In order to check the validity of the LUT1, the background lines were used to calibrate the in-flight gains and offsets. It is important to note that the calibration of rise time is especially difficult because there are no special features in the rise time distribution. Therefore, it is necessary to build a biparametric diagram for the 511 keV line for each pixel, but this requires a much larger amount of data than is currently available.

The calibration data spectrum, shown on figure 2 is characterised by several lines: the fluorescence lines of tungsten (59.3 keV and 67.2 keV) and lead (75 keV and 84.9 keV), and also the 511 keV line a portion of which is due to $`{}_{}{}^{22}\mathrm{Na}`$ calibration source (see Bird et al. (2003)). These is also a large peak at low energies (around 24 keV) whose origin is not completly understood. The fluorescence lines of Cd (23.1 and 26.9 keV) and Te (27.5 and 31 keV) certainly contribute. Besides, the irradiation of CdTe creates a large number of radioactive isotopes among which $`{}_{}{}^{111\mathrm{m}}\mathrm{Cd}`$ could contribute (Murakami et al. (2003)). It has a life time of 45 minutes, and decay with the emission of a 150 keV line (which is marginaly detected because of the reduced peak efficiency), and also lines between 23 and 26 keV.

### 3.1 Thermal effects

The temperature dependence of gains and offsets has been evaluated on ground for temperatures ranging from $`10^{}\mathrm{C}`$ to $`20^{}\mathrm{C}`$. It was found that the offsets are roughly constant, but that both pulse height and rise time gains strongly depend on temperature. They were fitted with the following expressions:

$$\frac{g_{\mathrm{PH}}}{g_{\mathrm{PH}}(273\mathrm{K},100\mathrm{V})}=\left(\frac{T}{273\mathrm{K}}\right)^{1.11}\left(\frac{U}{100\mathrm{V}}\right)^{0.08}$$

$$\frac{g_{\mathrm{RT}}}{g_{\mathrm{RT}}(273\mathrm{K},100\mathrm{V})}=\left(\frac{T}{273\mathrm{K}}\right)^{0.52}\left(\frac{U}{100\mathrm{V}}\right)^{0.58}$$

$`g_{\mathrm{PH}}`$ and $`g_{\mathrm{RT}}`$ being expressed respectively in keV per channel and $`\mu \mathrm{s}`$ per channel.

The temperature variations can be as large as 10K in a single revolution, if the pointing is changed or some PICSIT modules are switched off in the radiation belts. This can induce gain variations of the order of 5% which can be strong limiting factors for the energy resolution especially at high energies where the relative resolution is significantly lower. Thermal sensors are placed on the side of every ISGRI module as shown on figure 3 and allow us to correct for temperature effects on a module per module basis.

The pulse height gains and offsets were then calibrated using the 59.3 keV fluorescence line, and the 511 keV line in S2 data. A systematic bias was observed between the LUT1 and the coefficients obtained. This bias can be interpreted in terms of temperature variations in the camera, for the on-ground calibration was not done in vacuum; the thermal equilibrium is therefore very different. Using the previous formula, an effective temperature map was deduced, and is shown in figure 3. Once the temperature is taken into account, the pulse height calibration errors are very small. Note also that this correction is the only one applied to the rise time gains and offsets.

### 3.2 Charge loss corrections, LUT2

Charge losses have been measured on ground using various radioactive sources: $`{}_{}{}^{22}\mathrm{Na}`$, $`{}_{}{}^{139}\mathrm{Ce}`$, $`{}_{}{}^{54}\mathrm{Mn}`$, $`{}_{}{}^{57}\mathrm{Co}`$, $`{}_{}{}^{133}\mathrm{Ba}`$, etc. The pulse height – rise time relations obtained for various energies were interpolated to produce a general charge correction look-up table (LUT2). This is detailed in Lebrun et al. (2003).

The in-flight bias voltage and temperature conditions being different from those used for the LUT2 production, a simple correction was applied to it in order to recover good performances for the 59.3 and 511 keV lines: a multiplicative factor and an offset have been applied to the LUT2 so that both lines match those observed in the biparametric diagram of the flight calibration data (see figure 5).

### 3.3 Energy resolution

After correcting for temperature effects, we evaluated the energy resolution for the 59.3 keV and 511 keV lines using S2 data, and studied its variation with rise time. The resolution obtained for low rise time is coherent with expectations, but it increases faster with rise time as is shown on figure 4. At 511 keV, the integrated width (up to 4 $`\mu s`$) is of the order of 25 keV (FWHM) resulting in a 4.9% resolution, which is larger than what was expected from ground measurements. This is primarly caused by the rapid degradation of resolution with rise time: from 3.5% to more than 20% at $`4\mu \mathrm{s}`$. Poor calibration of rise time gain and offset is mainly responsible for the degradation, and an effort is needed in order to improve these performances. This effect is less important at lower energies where the charge losses are less important. At 59.3 keV, the integrated FWHM is 5 keV resulting in a 8.4% resolution. In this case, the variation with rise time is small, since charges losses are low, see figure 4.

## 4 Instrumental background

The background shape and spectrum is a key issue in optimizing the sensitivity of a coded mask imaging system like IBIS. The various contributions to the background include: galactic and extragalactic diffuse gamma-ray emission that dominates in the low energy range, nuclear lines produced in the external payload and spacecraft mass distribution, as well as scattering of cosmic rays. In the standard analysis (Goldwurm et al. (2003)) a background model is subtracted from the shadowgram before the deconvolution is performed. This background model is produced using empty fields observations or by summing shadowgrams obtained in dithering mode. For more details on the background non-uniformity in general, see Natalucci et al. (2003).

We have studied the background spatial and spectral behaviour for two different VETO configurations: using only lateral anticoincidence (veto LAT) up to revolution 33, and with lateral and bottom coincidence (veto ALL) for subsequent revolutions. To build a background model, we have used two empty field staring observations (rev 14 and 24) and a part of revolution 38 as well as some high latitude observations where the source contribution is weak (veto ALL). Images were produced with these data in various energy bands correcting for bad and noisy pixels as well as pixel efficiency in each science window. These images are then used as background model, leaving the normalization parameter free, because the background structures are supposed to be constant, but the global flux is expected to be variable.

Figure 6 (left column) shows the background model obtained in several standard imaging energy bands in the vetoALL configuration. The significance mosaic maps obtained with the standard analysis for three science windows with (figure 6, middle) and without (figure 6, right) correction are also shown. In the low energy regime ($`<`$40 keV), a rather flat structure is observed in the camera. The shadow of the aluminium frame is visible on the first pixel row in the upper and lower modules, suggesting an important contribution of diffuse gamma-rays to the background. From 40 to 100 keV, the background is dominated by the contribution of tungsten and lead fluorescence lines coming from the hopper and the mask. Since the incidence angle on the camera can be large, the aluminium frame’s shadow is much broader, thus creating convex module structures. At higher energies (see figure 6, 100-200 and 200-400 keV bands), the background is more intense on the sides of the camera, especially on the lower side closest to SPI, due to scattering in passive materials and possible leaks in the shielding. The pixels around the aluminium frame also show a strong flux enhancement for the same reasons, in this case creating a concave module structure. These structures produce very strong lines in the deconvolved images. In the vetoLAT configuration, the high energy regime is also characterized by a large bump in the whole detector producing a saddle shape in the resulting images.

Subtracting the background model from the shadowgrams yields flat sky maps, as shown on the right column of figure 6 which was obtained using 3 Galactic center observations science windows. The bright lines, due to background enhancement on the detector sides, seen above 100 keV are very well suppressed, but some errors remain although at a much smaller level. However, the background variations with time are large and the model used does not correct properly all observations, depending on latitude, cosmic ray fluxes, and status of the camera. The background structure variations are under study in order to quantify these as a function of time, latitude and cosmic ray flux.

A final point is that background substracted images have to be corrected for efficiency non-uniformity due to for instance gamma-ray scattering (see Natalucci et al. (2003)). This, however, is a second order effect, and it is not taken into account at the moment.

## 5 Conclusion

ISGRI is the first large CdTe gamma-ray imager used in orbit, and despite a few unexpected features like zero rise time events, it performs very well with only 4.5% noisy or disabled pixels. The energy calibration has been performed using the on board calibration source and background lines. Thermal effects are at the origin of the largest difference between ground and in-flight data. Correcting for these effects yields good spectral performances close to the expectations with 8.4% at 59.3 keV and 4.9% at 511 keV. The resolution in the high energy band is broader than before launch because of residual rise time gains uncertainties. Handling of these errors requires a larger amount of calibration data than what is available today.

Background non-uniformity in the camera has also been studied, and a simple energy dependent model has been produced using empty field observations. It provides rather good flat fielded images, even above 100 keV where the background non uniformity is largest. Unfortunately, background variations with time are large, and further studies are required in order to attain the full sensitivity.

###### Acknowledgements.

R.T. acknowledges financial support from CNES. This work is based on observations with INTEGRAL, an ESA project with instruments and science data centre funded by ESA member states (especially the PI countries: Denmark, France, Germany, Italy, Switzerland, Spain), Czech Republic and Poland, and with the participation of Russia and the USA.


text:
# Ultra peripheral heavy ion collisions and the energy dependence of the nuclear radius

## Acknowledgments

We are grateful to Y. Hama for discussions. This research was supported by the Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) (AAN) and by Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP) (CGR).


text:
# WIMP direct detection and halo structure

## 1. WIMPs a (very) brief introduction

Any stable weakly interacting massive particle (WIMP) in thermal equilibrium in the early universe will generically have an interesting present day density, $`\mathrm{\Omega }_{\mathrm{WIMP}}𝒪(\mathrm{\Omega }_{\mathrm{CDM}})0.3`$. Furthermore supersymmetry provides a natural WIMP candidate, the lightest supersymmetric particle, the neutralino. There are basically two methods of detecting WIMPs: indirect detection, which involves detecting the products of WIMP annihilation ($`\gamma `$, $`\nu `$, $`\overline{p}`$, $`e^+`$), and direct detection, which involves detecting the energy deposited in a detector due to elastic scattering of WIMPs on the detector nuclei. I will focus on WIMP direct detection. For a review of particle dark matter see e.g. Bergström (2000).

## 2. Direct detection signals

Direct detection experiments are just reaching the sensitivity required to detect WIMPs. The expected event rates are very small ( $`𝒪(10^510)`$ counts $`\mathrm{kg}^1\mathrm{day}^1`$) and distinguishing a putative Weakly Interacting Massive Particle (WIMP) signal from backgrounds, such as neutrons from cosmic-ray induced muons or natural radioactivity, is crucial. The event rate depends on the velocity of the detector relative to the Galactic rest frame and the Earth’s motion (as shown in Fig. 1) provides two potential WIMP smoking guns. Firstly the event rate is direction dependent, being greatest in the forward direction (Spergel 1988). Secondly the Earth’s velocity, and hence the event rate, varies annually due to the Earth’s orbit about the Sun (Drukier, Freese & Spergel 1986). If the local WIMP velocity distribution is isotropic then the annual modulation is roughly sinusoidal with a maximum in early June (when the Earth’s speed with respect to the Galactic rest frame is largest) and amplitude of order a few per-cent.

There are currently more than 20 WIMP direct detection experiments being carried out around the world. I will focus on those currently producing the most interesting results. The DAMA collaboration, using NaI with an exposure of 108 000 kg-days at Gran Sasso, have detected an annual modulation with the properties described above which they interpret as a WIMP signal (Bernabei et al. 2003). Assuming a standard halo model with a Maxwellian velocity distribution with dispersion $`\sigma =270\mathrm{k}\mathrm{m}\mathrm{s}^1`$ (corresponding to an asymptotic halo circular velocity of $`v_\mathrm{c}=220\mathrm{km}\mathrm{s}^1`$), they find a best fit WIMP mass $`m_\chi 50`$ GeV and scattering cross-section $`\zeta \sigma 7\times 10^{46}\mathrm{m}^2`$111$`\zeta =\rho _\chi /(0.3\mathrm{GeV}\mathrm{cm}^3)`$ parameterizes the uncertainty in the local WIMP density $`\rho _\chi `$.. Also of interest are three other experiments using different targets and a different strategy: Cryogenic Dark Matter Search (Ge target, 28 kd-days exposure at the Stanford Underground Facility (Akerib et al. 2003)), Edelweiss (Ge, 32 kg-days at Modane (Benoit et al. 2002)) and Zeplin I (liquid Xe, 230 kg-days at Boulby (Barton et al. 2002)). With their smaller exposures these experiments aim to constrain the mean WIMP scattering rate, rather than attempting to detect the annual modulation signal. The full DAMA allowed region of WIMP mass cross-section parameter space is shown in fig.2, along with the exclusion limits from the CDMS, Edelweiss and Zeplin I experiments. The allowed region and the exclusion limits are all calculated assuming the standard halo model as described above and, taken at face value, the allowed region appears to be incompatible with the exclusion limits.

## 3. Local velocity distribution

The differential elastic scattering rate depends on the local WIMP density, $`\rho _\chi `$, and the normalised WIMP speed distribution, in the rest frame of the detector, $`f_v`$, as

$$\frac{\mathrm{d}R}{\mathrm{d}E}\rho _\chi _{v_{\mathrm{min}}}^{\mathrm{}}\frac{f_v}{v}dv,$$

where $`v_{\mathrm{min}}`$ is the minimum WIMP velocity which can kinematically produce a nuclear recoil of energy $`E`$.

As discussed above, data analyses usually assume that WIMPs have an isotropic Maxwellian velocity distribution (i.e. that the Milky Way halo is an isothermal sphere), but observed and simulated halos are triaxial, anisotropic and contain substructure. Exclusion limits, which depend on the time averaged speed distribution, vary by of order tens of per-cent when triaxial and anisotropic halo models (with parameters chosen to match the properties of observed and simulated halos) are considered, with the shift being experiment dependent (Green 2002).

The annual modulation signal is far more sensitive to the WIMP velocity distribution, and Belli et al. (2002) found that considering “non-standard” halo models led to a large increase in the size of the region of $`m_\chi \sigma `$ parameter space consistent with the DAMA annual modulation signal. Accurate calculation of the shape and phase of the annual modulation signal requires all three components of the Earth’s velocity with respect to the Sun, and also the Sun’s motion with respect to the Local Standard of Rest, to be taken into account (Green 2003). Significantly, if the velocity distribution is not close to isotropic then the phase and shape of the annual modulation change and become incompatible with the DAMA annual modulation signal (Copi & Krauss 2003, Fornengo & Scopel 2003, Green 2003).

## 4. (very) Small scale structure

The velocity distributions of analytic halo models are derived via the Jean’s equations, which assume that the phase space distribution function has reached a steady state. In Cold Dark Matter (CDM) cosmologies structure forms hierarchically, with galaxy halos forming from the merger and accretion of smaller subhalos (which themselves formed from even smaller subhalos), and the local velocity distribution may not have reached a steady state. Helmi, White & Springel (2001) examined the velocity distribution of particles within a 4 kpc box located 8kpc from the centre of a simulated Milky Way like halo and found that, apart from a stream of fast moving particles from a late accreting subhalo, the velocity distribution was well approximated by a multi-variate gaussian. The WIMP direct detection event rate, however, depends on the dark matter distribution on sub-mpc scales, many orders of magnitude smaller than the scales probed by even the highest resolution galaxy simulations. Moore et al. (2001) argued that the local velocity distribution will depend sensitively on the structure and merger history of the first halos to form, while Stiff & Widrow (2003) have used a novel reverse technique to probe the velocity distribution at a single point in a simulation, finding that it appears to consist of a series of discrete peaks.

The power spectrum on very small scales is an essential input for any attempt to study the fate of the first subhalos to form and the dark matter distribution on very small scales. For neutralino CDM collisional damping and free-streaming erase power on comoving scales smaller than of order a pc (Hofmann, Schwarz & Stöcker 2001, Green, Hofmann & Schwarz 2003), and the first subhalos to form have mass of order $`10^6M_{}`$ (Green et al. 2003), 12 orders of magnitude lighter than the smallest subhalos resolved in galaxy simulations.

Streams of WIMPs with small velocity dispersion (from late accreting subhalos which pass through the solar neighbourhood or from the remnants of the first dense subhalos to form, should they survive to the present day) will produce steps in the differential event rate, the position and amplitude of which vary annually.

## 5. Summary

Direct detection of WIMPs would confirm the existence of Cold Dark Matter (and probe particle physics beyond the standard model). Accurate astrophysical input (not just the local WIMP velocity distribution, but also the motion of the detector with respect to the Galactic rest frame) is required when calculating the WIMP annual modulation signal. Analyzing data assuming a sinusoidal modulation with fixed phase could lead to erroneous constraints on, or best fit values, for the WIMP mass and cross-section, even worse a WIMP signal could be overlooked. On the other hand using unrealistic halo models or parameter values could lead to overly restrictive exclusion limits or a misleadingly large range of allowed values of the WIMP mass and cross-section. Finally if WIMPs are directly detected then we will be able to probe the local velocity distribution and perhaps learn about the (sub-)structure of the Milky Way halo.

### Acknowledgments.

I would like to thank my collaborators on the neutralino power spectrum on sub-galactic scales, Stefan Hofmann and Dominik Schwarz.

## References

Akerib, D. S. et al. 2003, hep-ex/0306001
Barton, J. C. et al. 2002, proceedings of the Fourth International Workshop on the Identification of Dark Matter, ed. N. J. C. Spooner & V. Kudryavtsev (World Scientific), 302
Belli, P., et al. 2002, Phys.Rev.D66, 043503
Benoit, A. et al. 2002, Phys.Lett.B545, 43
Bergström, L. 2000, Prog.Phys.63, 793
Bernabei, R. et al. 2003, Riv.N.Cim23, 1
Copi, C. J. & Krauss, L. M. 2003, Phys.Rev.D67, 103507
Drukier, A. K., Freese, K. & Spergel, D. N. 1986, Phys.Rev.D33, 3495
Fornengo, N. & Scopel, S. 2003, astro-ph/0301132
Green, A. M. 2002, Phys.Rev.D66, 083003
Green, A. M. 2003, Phys.Rev.D68, 023004
Green, A. M., Hofmann, S. & Schwarz, D. 2003, astro-ph/0309621
Helmi, A., White, S. D. M. & Springel, V. 2002, Phys.Rev.D66, 063502
Hofmann, S., Schwarz, D. J. & Stöcker, H. 2001, Phys.Rev.D64, 083507
Moore, B. et al. 2001, Phys.Rev.D64, 063508
Spergel, D. N. 1988, Phys.Rev.D37, 1353
Stiff, D. & Widrow, L. M. 2003, Phys.Rev.Lett90, 211301