text:
# Optimization with numerical simulations of the conjugate altitudes of deformable mirrors in an MCAO system

## 1 Introduction

Adaptive Optics (AO) was first proposed by Babcock (1953) as a means to compensate for the loss of spatial resolution of telescopes due to the disturbing effects caused by atmospheric turbulence (Roddier 1981) and thereby obtain diffraction limited resolution. A major limitation of AO is the small size of the Field of View (FoV) over which the correction is effective, referred to as the isoplanatic patch. Under good seeing conditions a typical value for the isoplanatic radius is just a few (2-3) arcsec in the visible ($`0.5\mu m`$) and about 12-18 arcsec in the near infrared ($`2.2\mu m`$). The concept of MCAO as a technique to deliver large corrected FoVs was proposed by Dicke (1975) and later developed by Beckers (1988) with the main goal of increasing the corrected FoV. It did not however receive much attention at that time as the AO community focused on conceptual problems such as obtaining adequate sky coverage with and without the use of Laser Guide Stars (LGS) (Sandler et al. 1994; Le Louarn et al. 1998) and the problems associated with LGS-based AO systems, in particular tip-tilt indetermination and the cone effect. During the early ’90s the concept of atmospheric turbulence reconstruction in 3-D (i.e. tomography) was proposed as a solution to overcome the problems associated with sky coverage and the cone effect in LGS-based AO systems (Tallon & Foy 1990). Shortly afterwards, the theoretical foundations of MCAO were developed (see e.g. Johnston & Welsh 1994; Ellerbroek 1994). In the last few years it has been recognised that tomography and MCAO offer the possibility of attaining large corrected FoVs and there is a large on-going effort in MCAO on both the theoretical side (see e.g. Ragazzoni et al. 1999; Fusco et al. 2000; Ragazzoni et al. 2000; Rigaut et al. 2000; Fusco et al. 2001) and using numerical simulation (see e.g. Flicker et al. 2000; Carbillet et al. 2001a; Conan et al. 2001; Bello et al. 2001; Ellerbroek 2001; Tordi et al. 2001; Tokovinin et al. 2001; Femenía et al. 2001; Ellerbroek & Rigaut 2001; Le Louarn 2002) to assess the performance of MCAO systems.

Conventional AO systems employ a single Deformable Mirror (DM) usually conjugated to the telescope entrance pupil. MCAO requires at least a second DM conjugated to some altitude above the pupil. The addition of a second DM represents a modest increase in complexity (of course the wave-front sensor also has to be modified). Such a dual-conjugate mode is considered as a short term upgrade to the AO system being designed for the 10m GTC telescope (Devaney et al. 2002) at the Observatorio del Roque de los Muchachos (ORM, Spain). In this work we present numerical simulations of a dual-mirror modal MCAO system for the case of a 10-m telescope. The main goal of our simulations is to determine the optimal conjugation altitude of the high-altitude DM (from now on referred to as $`DM_2`$) assuming that the low-altitude DM ($`DM_1`$) is conjugated to ground level. Our numerical simulations investigate the dependence of the optimal high-altitude conjugation altitude on the nature of the GS constellation (either NGS or NGS plus LGS constellations are considered), the geometry and the number of reference sources in the GS constellation. These simulation results will be compared with theoretical predictions based on the estimates of Tokovinin & Viard (2001) for the optimal conjugation altitude of $`DM_2`$.

## 2 MCAO model

### 2.1 Modal tomography approach

In modal tomography (Lawrence & Chow 1984; Ragazzoni et al. 1999) one seeks to estimate the coefficients of Zernike expansions of the wavefronts over metapupils defined at the conjugate altitudes of the DMs. The metapupils are circles enclosing the footprints of the beams from the guide stars. If $`\stackrel{}{b}_{(\mathrm{g})}`$ denotes the Zernike polynomial decomposition of the incoming wavefront when looking at guide star $`g`$ in the direction $`\stackrel{}{\theta }_{(\mathrm{g})}`$ (i.e. $`W(\stackrel{}{x};\stackrel{}{\theta }_{(\mathrm{g})})=_ib_{i,(\mathrm{g})}Z_\mathrm{i}(\stackrel{}{x}/R)`$, with R the pupil radius) then the Zernike polynomial decomposition of the wavefront defined on each of the metapupils to be reconstructed at heights $`h_l`$ ($`l=1\mathrm{}L`$) is $`\left\{\stackrel{}{a}^{(l)}\right\}_{l=1\mathrm{}L}`$ which is related to $`\left\{\stackrel{}{b}_{(g)}\right\}_{g=1\mathrm{}G}`$ by :

$$\stackrel{}{b}^\mathrm{T}=𝐏\stackrel{}{a}^\mathrm{T},$$
(1)
where the subscript T denotes transpose and the different terms in Eq. (1) are:

$`\stackrel{}{b}`$ $`=`$ $`\left[\begin{array}{cccc}\stackrel{}{b}_{(1)}\left(\stackrel{}{\theta }_{(1)}\right)& \stackrel{}{b}_{(2)}\left(\stackrel{}{\theta }_{(2)}\right)& \mathrm{}& \stackrel{}{b}_{(\mathrm{G})}\left(\stackrel{}{\theta }_{(\mathrm{G})}\right)\end{array}\right],`$ (3)
$`𝐏`$ $`=`$ $`\left[\begin{array}{cccc}𝐏(\stackrel{}{\theta }_{(1)};h_1)& 𝐏(\stackrel{}{\theta }_{(1)};h_2)& \mathrm{}& 𝐏(\stackrel{}{\theta }_{(1)};h_\mathrm{L})\\ 𝐏(\stackrel{}{\theta }_{(2)};h_1)& 𝐏(\stackrel{}{\theta }_{(2)};h_2)& \mathrm{}& 𝐏(\stackrel{}{\theta }_{(2)};h_\mathrm{L})\\ \mathrm{}& \mathrm{}& \mathrm{}& \mathrm{}\\ 𝐏(\stackrel{}{\theta }_{(\mathrm{G})};h_1)& 𝐏(\stackrel{}{\theta }_{(\mathrm{G})};h_2)& \mathrm{}& 𝐏(\stackrel{}{\theta }_{(\mathrm{G})};h_\mathrm{L})\end{array}\right],`$ (8)
$`\stackrel{}{a}`$ $`=`$ $`\left[\begin{array}{cccc}\stackrel{}{a}^{(1)}& \stackrel{}{a}^{(2)}& \mathrm{}& \stackrel{}{a}^{(L)}\end{array}\right].`$ (10)
and in Eq. (8) each of the $`𝐏(\stackrel{}{\theta }_{(\mathrm{g})};h_\mathrm{l})`$ is a Zernike projection matrix relating the Zernike decomposition coefficients of the wavefront defined on the metapupil ($`\stackrel{}{a}^{(\mathrm{l})}`$) and $`\stackrel{}{b}_{(\mathrm{g})}^{(\mathrm{l})}(\stackrel{}{\theta }_{(\mathrm{g})};h_\mathrm{l})`$ are the Zernike decomposition coefficients of the wavefront defined on the intersection of the metapupil and the propagation path at the $`l`$th layer when looking toward $`\stackrel{}{\theta }_{(\mathrm{g})}`$:

$$\stackrel{}{b}_{(\mathrm{g})}^{(\mathrm{l})}(\stackrel{}{\theta }_{(\mathrm{g})};h_\mathrm{l})=𝐏(\stackrel{}{\theta }_{(\mathrm{g})};h_\mathrm{l})\stackrel{}{a}^{(\mathrm{l})}.$$
(11)
The existence of such $`𝐏(\stackrel{}{\theta }_{(\mathrm{g})};h_\mathrm{l})`$ matrices is demonstrated in Ragazzoni et al. (1999) and an analytical expression for their computation is provided in Femenía & Brusa (2003) and depends exclusively on the geometry of the problem (i.e. on the metapupil radius, radius of the circular region defined by the intersection of the metapupil and the propagation path and the relative position of the circular region with respect to the metapupil center). The required $`\left\{\stackrel{}{a}^{(l)}\right\}_{l=1\mathrm{}L}`$ coefficients may be obtained from Eq. (1) by a proper inversion method (e.g. SVD or least square estimator, optimal reconstructor or maximum-a-priori reconstructor).

### 2.2 Implementation in the simulation

The previous description of modal tomography does not give any indication as to how the Zernike polynomial coefficients $`\stackrel{}{b}`$ of the GS wavefronts are retrieved from real data. In a real wavefront sensor the data will be wavefront local slope measurements (as in the case of Shack-Hartmann and pyramid wavefront sensors) or wavefront local curvature values (as in the case of the curvature sensor). We require a reconstruction process from real measurements to the metapupil Zernike coefficients $`\stackrel{}{a}`$ without having to first reconstruct the GS Zernike coefficients $`\stackrel{}{b}`$. This is achieved by means of the *interaction matrix* $`𝐌`$ of the MCAO system linking wavefront sensor outputs ($`\stackrel{}{s}`$) and metapupil Zernike coefficients: $`\stackrel{}{s}=𝐌\stackrel{}{a}`$. While the $`𝐏`$ matrix in Sect. 2.1 is obtained from purely geometric considerations, building the interaction matrix $`𝐌`$ requires explicit knowledge of the WFS such as the extended size of the GS spot as seen by the WFS, pixelization effects in the WFS or the fact that the WFS sees only a limited FoV. This leads us to use interaction matrices generated by the simulation in which we sequentially place different Zernike polynomials on each layer and register the local tilts produced on each subaperture of each Shack-Hartmann Sensor (SHS) in our MCAO system. In the simulations using LGSs we consider that there is a global tilt subtraction for the wavefront coming from each LGS.

## 3 Description of the software employed: CAOS

The main software tool used is the IDL-based Code for Adaptive Optics Systems (CAOS, Carbillet et al. 2001b, 1999). The main feature of CAOS is its modular structure: each elementary physical process is modeled by a specific set of routines (i.e. a module) , e.g. the simulation of atmospheric turbulence or the propagation of the wavefront through turbulent layers. This modular design makes CAOS a very versatile tool as it can be extended by the user with new modules as required. The other main feature of the CAOS software is its graphic user interface, the “Application Builder” (Fini et al. 2000), which helps to design simulation projects in a very straightforward way even for users with no knowledge of IDL. After the graphical design, the Application Builder automatically generates an IDL code which can be run as the desired simulation. In order to simulate MCAO systems new modules have been developed in addition to the modules in the standard CAOS library. A dedicated module conducts the calibration of the MCAO system yielding the interaction matrix: at each layer to be reconstructed we place a single polynomial $`Z_\mathrm{j}`$ at a time and perform the geometrical propagation (i.e near-field approximation) of the wavefront error corresponding to this $`Z_\mathrm{j}`$ along the line of sight of each Guide Star (GS) down to the SH sensors. The vector of x- and y-slopes from all the SH sensors gives us the column of the matrix associated with the j-th Zernike polynomial on the considered layer; this process is then repeated for all the modes to be reconstructed on each layer. The reconstruction process is done with a dedicated module which uses a Least Square Estimator (LSE) and by inverting with Singular Value Decomposition (SVD) the linear relation $`\stackrel{}{s}=𝐌\stackrel{}{a}`$ between the measured slopes ($`\stackrel{}{s}`$) and the coefficients ($`\stackrel{}{a}`$) on each reconstructed metapupil. The multi-mirror reconstructor is simulated assuming perfect Zernike DMs in the sense that any mirror is able to generate a pure Zernike mode up to a given mode $`j`$ and allows different degrees of correction to be used at different layers.

## 4 Simulation description & simplifications

### 4.1 Generation of Guide Stars (GS)

NGS are simulated as point-like sources at an infinite distance from the telescope and sampling a cylindrical turbulence volume. LGS are simulated as extended sources at a finite distance so that they sample a cone shaped turbulence volume. The width of the LGS is evaluated taking into account the broadening due to diffraction by a launching telescope of 0.5 m diameter, the upward propagation of the laser beam, the perspective elongation of the LGS beacon as seen from the observing telescope and the geometry of the beacon size averaged over the range gate (for Rayleigh LGS) or width of the Na-layer (for Na LGS) using the expression presented in Ellerbroek & Rigaut (2001):

$$\theta _\mathrm{B}^2=\left(\frac{\lambda }{d_\mathrm{T}}\right)^2+\left(\frac{\lambda }{r_0}\right)^2+\frac{1}{3}\left(\frac{L\delta h}{h^2}\right)^2+\left(\frac{d_\mathrm{T}\delta h}{2\sqrt{3}h^2}\right)^2$$
(12)
where $`\lambda `$ is the wavelength of observation ($`\lambda =589.3\mathrm{nm}`$ for Na LGS, for Rayleigh LGS we consider $`\lambda =355\mathrm{nm}`$), $`d_\mathrm{T}`$ the launching telescope diameter, $`r_0`$ the Fried parameter at the WFS operation wavelength, $`L`$ the distance between the launching telescope and the SHS subaperture, $`h`$ the altitude of the LGS and $`\delta h`$ the “range gate”. By “range gate” we mean the vertical distance of the LGS on which we integrate while it is propagating upward. Such a value can be fixed by us when considering Rayleigh LGS (by means of e.g. a Pockels cell). In the case of a Na LGS it is not possible to fix a range gate since the Sodium layer drifts in altitude and thus in the above expression we must use the full width of the Sodium layer weighted by the fact that the sodium layer has a Gaussian density profile so that an appropriate value in the above expression is $`\delta h=5`$ Km. Then the LGSs are projected onto the WFS so that they look like circular Gaussian spots of width equal to $`\theta _\mathrm{B}`$. This neglects the differential elongation due to the finite distance from the telescope to the LGS but on the other hand such an approach is pessimistic since the LGS is not circular but elongated and in practice this involves a loss of resolution only along the elongation direction.

### 4.2 Generation of atmospheric turbulence

For the purpose of our simulations we have constructed a 7-layer profile summarized in Table 1. This profile was built using the measurements of balloon flights (Vernin & Muñoz-Tuñón 1994) launched from the Observatorio del Roque de los Muchachos (ORM). Data from the 6 available balloon flights were used to determine an averaged profile in the following way; first all profiles from the balloon measurements were normalized to have the same total $`r_0`$ as the annual mean (i.e. $`r_0(500\mu \mathrm{m})=0.15\mathrm{m}`$) at ORM and then a continuous average profile was generated. Then, following Ellerbroek et al. (1994), a 7-layer turbulence profile was obtained by performing a discrete fit to the continuous average profile that matches the first 14 moments of the continuous turbulence profile.

Atmospheric turbulent phase screen generation was performed by summing a finite set of Zernike polynomials weighted by a set of coefficients drawn from a Gaussian multivariate distribution of zero mean and covariance matrix as given in Noll (1976), i.e. assuming the turbulence follows Kolmogorov statistics. The truncation of the infinite series implies missing high spatial frequency components. This is easily obviated by truncating at a large enough mode $`j_{\mathrm{max}}`$. The total Strehl ratio due to all missing modes on all layers is estimated to amount to about 10% and this reduction is applied to the iso-Strehl maps. A marginal effect is due to the actual statistics followed by the Zernike coefficients on portions of the phase screens the same size as the telescope pupil. This problem can also be circumvented by the use of a large enough number of Zernike modes when generating the turbulence phase screens.

### 4.3 Propagation through turbulence and wavefront sensing

Propagation of light through atmospheric turbulence is performed by adding linearly the phase fluctuation produced by each turbulent layer which assumes the near-field approximation (Roddier 1981) and corresponds to neglecting diffraction effects in the propagation of the phase of the complex field between layers and from layers to the ground (a similar approximation is implicit when assuming thin layers in Sect. 4.2: diffraction effects between the top and the bottom of the layer are neglected). Ellerbroek (2001) presents results using a code implementing wave-optics for phase propagation and makes a comparison of the results obtained under the near-field approximation for the case of an MCAO system at the Gemini telescope. He concludes that including Fresnel-diffraction in the propagation of light involves a phase rms of only a few nm (¡20 nm) with respect to the results obtained with the geometric propagation. Another effect related to not including Fresnel diffraction is that we neglect the effect of scintillation on wavefront sensing. CAOS implements a wave-optics code to simulate the WFS. However the problems of scintillation are not handled since they do not appear in the near-field approximation. Our simulations assume Shack-Hartmann WFS for which the performance dependence on the scintillation has been studied in Voitsekhovich et al. (2001) in the case of open-loop measurements. A direct extrapolation of the results of Voitsekhovich et al. (2001) indicates that for our case with $`r_0(500\mathrm{nm})=0.15`$ m scintillation effects in the WFS will involve errors $`<`$ 8% in the estimation of the high order reconstructed modes. A much more important simplification is that of no noise in the WFS procedure as we are not concerned at this stage with problems associated with sky coverage and by not considering noise we are in the limit of very bright objects.

### 4.4 Least Square Estimator (LSE) reconstructor

Several works (Fusco et al. 1999; Brusa et al. 2000; Flicker et al. 2000; Fusco et al. 2001; Conan et al. 2001) have focused on estimators for MCAO which assume a real-time knowledge of the statistics of phase perturbations and the turbulence vertical distribution to derive optimal reconstructors operating in open-loop: the so-called Maximum a posteriori (MAP) estimators. In conditions of low Signal-to-Noise Ratio (SNR) it has been shown (Flicker et al. 2000; Conan et al. 2001; Le Louarn 2002) that these optimal reconstructors perform better than the SVD-based reconstruction while yielding similar results in the case of high SNR values. Since we work in the limit of very bright objects (see Sect. 4.3) nothing justifies using MAP estimators here except for the advantage of including knowledge of the vertical distribution of turbulence. On the other hand, it is clear that a MCAO system will be operated in closed-loop which imposes no problem for the implementation of the LSE reconstructor while, to our knowledge, statistics of the residual signal on the WFS have not been derived analytically to be implemented in the design of the MAP reconstructor.

### 4.5 Perfect Zernike deformable mirrors

The DMs considered in these simulations are assumed to be able to reproduce exactly the Zernike modes being reconstructed, that is, if the command is that of generating $`a_\mathrm{j}Z_\mathrm{j}`$ (with $`\mathrm{j}=2,\mathrm{},153`$) on $`DM_1`$ then this is exactly what we get (idem with $`DM_2`$ but $`\mathrm{j}=2,\mathrm{},91`$) . This assumption neglects the influence functions of real DMs which in addition to the commanded $`Z_\mathrm{j}`$ mode will probably introduce other modes at a much lower level. Realistic influence functions will be considered in future work using more detailed numerical simulations.

### 4.6 Open-loop operation.

Our simulation assumes open-loop operation. In doing this we can easily perform a statistical analysis by generating independent input turbulent wavefronts and retrieving the response of the MCAO system. Open-loop operation neglects degradation of the MCAO system performance due to time-delay errors while on the other hand it gives an underestimation of the performance of the system as compared to closed-loop operation since the WFS is more affected by problems due to the restricted FoV in the SHS and the extended size of the GS.

## 5 Simulation results

In the case of single-conjugate AO with moderate-good correction, the performance may be characterised adequately by a single parameter, namely the on-axis Strehl ratio. The Strehl ratio can be maximised by choosing a threshold value in the SVD reconstruction; singular values lower that the threshold are set to zero and this eliminates modes to which the system is insensitive. In the case of an MCAO system, the performance to be determined is the Strehl ratio over the field of view. We find in our simulations that the Strehl ratio maps are sensitive to the SVD threshold (for example, see Fig. 2). Since the SVD threshold cannot be selected *a priori* we need criteria for selecting the optimal threshold values for each configuration. In order to characterize SR maps, we extract two parameters which quantify the average Strehl ratio over the field and the amount of variation of the Strehl ratio over the field. For a given Strehl map we determine the value of the Strehl ratio averaged on concentric annular radii (the first annulus is the inner circle of radius $`7\mathrm{}`$ and concentric annuli have an outer radius $`2.5\mathrm{}`$ larger than the inner radius) and the root mean square deviation of Strehl in each annulus with respect to the averaged Strehl ratio on the circle of radius equal to the annulus outer radius. As shown in the example in Fig. 1 we obtain two families of curves: the first set of curves shows for each SVD threshold value the evolution of Strehl averaged on annuli as a function of the annulus outer radius while the second family of curves shows the evolution of Strehl $`rms`$. This process is conducted on contour plots such as those depicted in Fig. 2. To produce these contour plots we interpolate the reconstructed FoV sampled at 15x15 locations in the simulation to a finer grid. Due to this interpolation when we get close to the edges we observe a large gradient in S. To avoid a significant bias in the SVD selection process we generate curves of averaged Strehl ratio and Strehl $`rms`$ on concentric annuli up to a maximum outer radius of about 85% of the whole reconstructed FoV radius. Once we have generated the above sets of curves for each GS configuration and $`DM_2`$ conjugation altitude the SVD threshold value selection is based on (i) having flat curves of averaged Strehl within concentric annuli (in general we observe curves which are mostly monotonically increasing or decreasing, so that we select those curves with the lowest peak-to-valley values and/or the lowest gradient on annulus averaged Strehl), (ii) having the smallest Strehl $`rms`$ and (iii) in those situations where several threshold values yield similar results based on the previous two points we select the largest SVD threshold value. This is justified in terms of SNR considerations.

### 5.1 Performance using Natural Guide Stars

Fig. 2 shows K-band Strehl maps for the case of having three natural guide stars evenly spaced on a circumference of radius $`30\mathrm{}`$ in order to reconstruct a circular field of view of $`1\mathrm{}`$ diameter for several SVD threshold values. The panels correspond to setting different singular value thresholds; it can be seen that the performance is fairly constant until a large number (15) of SVD modes are rejected, some of them with important contributions of absolute Zernike tilts on each of the two reconstructed layers. When only a few SVD modes are rejected (i.e. from THR=$`7\times 10^5`$ to THR=$`3\times 10^4`$ involving 2 to 7 SVD modes rejected) the Strehl ratio peaks at the positions of the guide stars as would be expected considering the least-squares reconstruction process employed while the Strehl ratio across the reconstructed FoV is relatively uniform when rejecting 8 to 10 SVD modes (THR=$`4\times 10^4`$ to THR=$`10^3`$). This is a general behavior observed in our simulations: uniformity of the Strehl ratio across the reconstructed FoV is achieved at the cost of lowering the peak Strehl ratio at the GS positions by an adequate SVD filtering. Fig. 3 shows the best SVD Strehl maps when the three natural guide stars are on a circumference of radius of $`45\mathrm{}`$, corresponding to a reconstructed FoV of $`1.5\mathrm{}`$ diameter, and $`DM_2`$ conjugated to 6.1 and 10 Km above the telescope pupil. We observe that with $`DM_2`$ conjugated to 10 Km the average Strehl ratio is significantly smaller than when $`DM_2`$ is conjugated to 6.1 Km while the variation of the Strehl ratio across the field is more pronounced for $`DM_2`$ conjugated to 10 Km. Making measurements on the maps it is found that the average Strehl ratio drops from 0.43 to 0.20 and the normalised Strehl $`rms`$ increases from 16% to 48% (see also Fig. 4). The Strehl $`rms`$ value when $`DM_2`$ is conjugated to 10 Km (and to 13 Km) is due to such a strong center-to-edge Strehl variation that there is little sense in talking of an average reconstructed Strehl ratio for this scenario.

By varying the $`DM_2`$ conjugation altitude in our simulations we obtained that for constellations of 3 NGS the optimal conjugation altitudes are approximately 8.5 Km and 6 Km for reconstructed FoVs of $`1\mathrm{}`$ diameter and $`1.5\mathrm{}`$ diameter, respectively. These optimal altitudes are surprisingly low when we consider the vertical distribution of turbulence in the 7-layer model in Table 1. We would expect the pupil-conjugate DM to correct the turbulence in the four lower layers (up to 5 Km), and the second DM would therefore be positioned to correct the turbulence in the three upper layers (10-18 Km). Tokovinin & Viard (2001) published expressions allowing optimal conjugation altitudes to be calculated for up to three mirrors for a given turbulence profile. Applying these expressions to our case (i.e. two DMs, balloon 7-layer model), we obtain an optimal $`DM_2`$ altitude of 13 Km. We believe the reason for such a disagreement has to do with the strong assumptions in the aforesaid theoretical work: the authors assume that whatever the altitudes of conjugation for the DMs, tomography is able to perfectly reconstruct the wavefront at those altitudes so that the only effect considered is the correction of a continuous turbulence profile by a finite number of $`DM`$s. In fact, the further from the telescope we conjugate $`DM_2`$, the less overlap between the GS footprints at the metapupil at the $`DM_2`$ conjugation altitude. This effect will tend to push down the optimal conjugation altitude for $`DM_2`$. This statement is supported by the observation that the optimal conjugation altitude is lower for the case of a FoV of $`1.5\mathrm{}`$ than for the case of a FoV of $`1\mathrm{}`$. In order to further test this hypothesis, we carried out simulations using five natural guide stars and keeping the rest of the parameters of the simulation as in the case of 3 NGS. Increasing the number of NGS in our constellation from 3 to 5 has the effect of increasing the overlap between their footprints and thus we decrease the importance of voids in the metapupil for the wavefront reconstruction at that range. A comparison between the results achieved with constellations of 3 NGS and constellations of 5 NGS is given in the left top and left bottom panels in Fig. 4 where we plot the average Strehl ratio on the reconstructed FoV as a function of the $`DM_2`$ conjugation altitude and the Strehl ratio $`rms`$ variation, respectively (the right top and right bottom panels of Fig. 4 show the results obtained when considering different LGS constellations plus an on-axis NGS and which will be the subject of the following sections). We can see in Fig. 4 that when increasing the number of NGSs from 3 to 5 NGS the conjugation altitude of $`DM_2`$ also increases and this effect can only be due to a better sampling (and reconstruction) of the turbulence at the height at which $`DM_2`$ is conjugated. From these results we conclude that the effect neglected by Tokovinin & Viard (2001) is very important when determining the best $`DM_2`$ conjugation altitude.

We also compared two different geometries employing 5 NGS: placing all NGS on the vertices of a pentagon or putting 4 of them on the corners of a square and one at the center. It turns out that it is the latter configuration which yields better results (see top left panel in Fig. 4 and iso-Strehl contour plots in Fig. 5) despite the fact that it leaves unsampled a slightly bigger portion of the metapupil. This suggests that not only the fractional sampling of the metapupil is important but also the distance of points within the reconstructed FoV to the GS. Although we have not conducted any test to check our explanation of this effect we believe that this may be due to the LSE reconstruction in which by definition the best reconstructed positions are those corresponding to the GS positions while for a MAP reconstructor we may encounter a different behavior since it is built by optimizing the reconstruction at a set of positions which do not have to be coincident with those of the GS.

### 5.2 Performance using Sodium Laser Guide Stars

The probability of finding constellations of 3 (not to mention 5) bright natural guide stars arranged to nicely span fields of $`1\mathrm{\hspace{0.17em}1.5}\mathrm{}`$ is very low. As an indication, using the Guide Star Catalogue-II111Copyright Association of Universities for Research in Astronomy, Inc. Further details and access to database at http://www-gsss.stsci.edu/gsc/gsc2/gsc22\_release\_notes.htm we estimate that the probability of finding three stars brighter than $`m_\mathrm{R}=16`$, which is an optimistic value for limiting guide star magnitude, is approximately 3% for a field of $`2\mathrm{}`$ at galactic coordinates $`(l,b)=(90\mathrm{°},40\mathrm{°})`$, and this value will be further reduced when requiring the 3 NGS to be properly arranged. In view of these considerations it will be necessary to employ laser guide stars (LGS) in order to provide reasonable sky coverage for star-oriented MCAO systems, at least on 8-10m telescopes.

It has already been pointed out by different authors (see e.g. Ellerbroek & Rigaut 2001) that the indetermination of the tip-tilt of individual LGS leads to tip-tilt anisoplanatism in an LGS-based MCAO system. The tip-tilt anisoplanatism may be solved in a natural way by detecting the tip-tilt of several NGS in the FoV (Rigaut et al. 2000). In the Gemini-South MCAO system it is planned to reconstruct the tip-tilt using three dim (i.e. $`m_\mathrm{R}1920`$) natural guide stars, and the sky coverage in the K band is predicted (Rigaut et al. 2000) to be 15% and 80% at Galactic pole and $`b=30\mathrm{°}`$, respectively. Alternative schemes can be derived by considering that tip-tilt anisoplanatism mainly arises from turbulence-induced quadratic wavefront errors (i.e defocus and astigmatism) at altitudes which cannot be determined from the multiple-LGS measurements. If quadratic measurements are made simultaneously with guide stars at different ranges, then it may be possible to correct the tip-tilt anisoplanatism. One approach along these lines is to make wavefront measurements on a natural guide star simultaneously with the laser guide stars. We determine the performance of such an approach here for fields of view of diameter $`1\mathrm{}`$ and $`1.5\mathrm{}`$. In the case of the smaller field of view we employ four Na LGS on a circumference spanning the field of view plus one on-axis NGS. In the case of the $`1.5\mathrm{}`$ field of view we consider two Na LGS configurations; (i) four LGS spanning the field of view plus one on-axis NGS (ii) five LGS on a circumference spanning the field of view plus one on-axis NGS. The circumferences on which the LGS are placed are somewhat larger than the fields of view so that the metapupil sizes are the same when using cone shaped beams from the LGS as when using cylindrical-shaped beams from NGS placed at the edge of the field of view. Fig. 6 shows the best Strehl maps for the $`1\mathrm{}`$ FoV (left panel) and $`1.5\mathrm{}`$ FoV (right panel), employing constellations of 4 off-axis LGS+ on-axis NGS and 5 off-axis LGS + on-axis NGS, respectively . For the $`1\mathrm{}`$ FoV case the average Strehl ratio is similar to that obtained on the same field of view using 3 NGS although it exhibits a larger peak-to-valley value. The $`rms`$ Strehl variation is similar to the 3 NGS case but it does not rise steeply when $`DM_2`$ is conjugate to an altitude higher than the optimal which for the LGS case is approximately 10 Km while in the case of NGS constellation it was about 8.5 Km. In the case of the $`1.5\mathrm{}`$ FoV the average Strehl ratio is almost the same for both 4LGS and 5LGS configurations. The optimal altitude occurs in the range 8-10 Km. If $`DM_2`$ is conjugate to an altitude higher than the optimal then the average Strehl ratio is higher and the rms variation is lower when 5 LGS are employed rather than 4 LGS.

### 5.3 Hybrid Rayleigh-Na system performance

The natural guide star employed to measure quadratic modes has to be bright, thereby again compromising sky coverage. An alternative approach to solving the tip-tilt anisoplanatism caused by the quadratic modes consists of employing both Rayleigh and Na LGS as well as a single natural guide star which is only used to determine global tip-tilt (Femenía et al. 2001; Ellerbroek & Rigaut 2001; Femenía 2003). Here we apply the algorithm described in Femenía (2003) whose starting point is to realize that unresolved second order modes in LGS-based tomography can be reconstructed using an LGS constellation made at several altitudes (e.g. Rayleigh + Sodium LGSs, or Rayleigh LGSs at different altitudes) as shown by Brusa et al. (2000). The reconstruction process consists of two stages in which barycenter tip-tilt is measured from a NGS and in parallel we perform tomography with the LGS constellation placed at at least two different altitudes. From the LGS tomography we obtain a reconstruction of all modes starting from the Zernike second-order modes and this information is used together with the NGS barycenter tip-tilt to obtain an estimate of the global Zernike tilts in the entire reconstructed FoV. The hybrid technique may be implemented in different ways, in particular Femenía (2003) considers a single on-axis Rayleigh LGS on which to conduct high-order wave-front sensing while Ellerbroek & Rigaut (2001) consider several off-axis Rayleigh LGS coupled to low-order WFS. Other parameters to be optimized include the height of the Rayleigh LGS and the wavelength for the Rayleigh LGS. Given the number of new parameters to be investigated the results we present here are by no means exhaustive but only preliminary to give an idea of the potential of the technique. We assume two configurations of 4 Na LGS and 5 Na LGS evenly arranged on a circumference of $`56.45\mathrm{}`$ radius (so that the reconstructed FoV is $`1.5\mathrm{}`$) plus an on-axis Rayleigh LGS. The Na LGSs are observed at $`589\mathrm{nm}`$ and the Rayleigh LGS at $`355\mathrm{nm}`$ with $`16\times 16`$ subaperture SHSs. The rest of the simulation parameters (atmospheric profile, number of CCD pixels and SHS FoV, no noise in the SHS, number of reconstructed modes, etc.) are the same as in the simulations with NGS and LGS+on-axis NGS. The top panels a) & b) in Fig. 7 shows the reconstruction across the $`1.5\mathrm{}`$ FoV for the cases of 4 Na LGS and 5 Na LGS, respectively.

On comparing the results from the constellation using 4 Na LGS against the results using 5 Na LGS we see a limited gain. This limited performance gain is a consequence of the limited reconstruction along the on-axis direction which in turn is used to reconstruct the Zernike tip-tilt modes across the entire reconstructed FoV. This statement is supported by the levels of reconstruction attained if we assume perfect reconstruction of the Zernike tip-tilt modes across the entire FoV as depicted in panel c) in Fig. 7 for the case of 4 Na LGS at $`56\mathrm{}`$ and in panel d) in Fig. 7 for the case of 5 Na LGS at $`56\mathrm{}`$; in both cases also considering on-axis Rayleigh LGS at 30 Km and on-axis Tip-tilt NGS. When comparing in Fig. 7 the bottom panels (perfect correction of Zernike tip-tilt across the entire FoV) against the top panels (actual correction of Zernike tilts) we verify that in the latter figures there is a small amount of tilt anisoplanatism caused by a worse correction of higher-order modes (especially second-order modes) in the on-axis direction. This in turn translates into tilt anisoplanatism when using those higher-order modes to reconstruct off-axis Zernike tilts. This situation also appeared when dealing with 5 NGS where we observed that the best configuration is that corresponding to 4 NGS on the corners of a square plus one on-axis NGS rather than all 5 NGS on the vertices of a pentagon. Based on this conjecture it would be more effective having 4 Na LGS distributed on the vertices of a square plus one on-axis Na LGS than having all 5 Na LGS distributed on the vertices of a pentagon.

## 6 Choosing the $`DM_2`$ conjugation altitude

As seen in Fig. 4 the optimal $`DM_2`$conjugation altitude depends very strongly on the GS nature and configuration. In the case of observing 3NGS on the vertices of equilateral triangles on a circumference of radius $`30\mathrm{}`$ or $`45\mathrm{}`$ the optimal altitudes are about 8.5 and 6.0 Km, respectively. When considering 5 NGSs the optimal $`DM_2`$ conjugation altitude is 10 Km. When hybrid LGS constellations and a single tip-tilt NGS are used aiming at a reconstructed FoV of $`1.5\mathrm{}`$ the best $`DM_2`$ conjugation altitude occurs in the range  Km with little variation of performance within that range leaving the final decision to opto-mechanical considerations. In order to find an optimal $`DM_2`$ conjugation altitude that best suits all possible GS configurations cases in this work we consider the following figure of merit:

$$\chi ^2(h)=\underset{i=1}{\overset{9}{}}\left(\frac{S_{\mathrm{max},\mathrm{i}}S_\mathrm{i}(h)}{S_{\mathrm{max},\mathrm{i}}}\right)^2w_\mathrm{i},$$

where $`h`$ denotes a generic $`DM_2`$ conjugation altitude, the index $`\mathrm{i}=1,\mathrm{},9`$ runs over all the GS configurations in this work, $`w_\mathrm{i}`$ weights different GS configurations, $`S_{\mathrm{max},\mathrm{i}}`$ is the maximum Strehl ratio achieved with the $`i`$th GS configuration and $`S_\mathrm{i}(h)`$ is the achieved Strehl ratio with the $`i`$th GS configuration when conjugating $`DM_2`$ to altitude $`h`$ above the telescope entrance pupil. As the $`\chi ^2`$ mimization process will depend on $`\left\{w_\mathrm{i}\right\}`$ we have considered the sets of weights listed in Table 2: the first weighting set (i.e. column $`w_1`$in Table 2) gives equal importance to all GS configurations, weighting set $`\left\{w_2\right\}`$ is focused on avoiding system complexity and thus assigns higher weight to simpler systems while weighting sets $`\left\{w_3\right\}`$ and $`\left\{w_4\right\}`$ are derived from simple sky coverage considerations based on star counts from the Guide Star Catalogue-2 at $`(l,b)=(90\mathrm{°},40\mathrm{°})`$ and assuming that for the NGS we require at least $`m_\mathrm{R}16`$ when it is used for high-order mode reconstruction and $`m_\mathrm{R}18.5`$ when it is only used for tip-tilt reconstruction (i.e. as in the case of the hybrid laser configurations). The difference between weighting sets $`\left\{w_3\right\}`$ and $`\left\{w_4\right\}`$ is that the former is proportional to the probability of finding the required number of NGSs in the GS configuration and the second is proportional to the square root of that probability (i.e. $`w_{4,\mathrm{i}}=\sqrt{w_{3,\mathrm{i}}}`$ where index $`i`$ runs over the different GS configurations under study). This choice for $`\left\{w_3\right\}`$and $`\left\{w_4\right\}`$ is very qualitative but given the small difference between the results obtained with them (see Table 3) we are confident that a detailed sky coverage analysis is likely to differ little from the results obtained here. The $`\chi ^2(h)`$ minimization process yields an “optimal” $`DM_2`$ conjugation altitude as a function of the weighting set. The first and second columns in Table 3 give the $`DM_2`$ conjugation altitude for which a given GS configuration attains its maximum Strehl ratio and the Strehl ratio itself. Columns 3 to 6 give Strehl ratios for each GS configuration when the $`DM_2`$ is conjugated to the joint optimal conjugation altitude according to the weighing set under consideration (i.e $`S_{\mathrm{opt},1}`$corresponds to using $`w_1`$); the value of such optimal conjugation altitude for each weighting set being displayed in parenthesis below each column title. The last two columns in Table 3 give the Strehl ratio and fractional decrease with respect to $`S_{\mathrm{max}}`$, respectively, when $`DM_2`$ is conjugated to 8 Km.

## 7 Conclusions

Modal MCAO in the “classical” star-oriented fashion holds the promise to reconstruct FoVs as large as $`1.5\mathrm{}`$ in diameter. Ideally one would like to rely on NGS constellations but the chances of finding a suitable one is extremely low even for the simplest case of 3NGS. Thus if one wants to build a star-oriented MCAO system which can be used over a significant fraction of the sky one must necessarily consider LGS-based MCAO systems. However, one then has to face the tip-tilt problem and the indetermination of second-order modes. The simplest approach to LGS-based MCAO systems is that of relying on a bright enough NGS on which one could perform not only tip-tilt & second-order sensing but also high-order WFS. This yields poor results and the initially surprising result that performance does not improve significantly when increasing the number of Na LGSs from 4 to 5. In any case, whether we use 4 Na LGS or 5 Na LGS the reconstructed FoV is strongly affected by tip-tilt anisoplanatism and the reconstruction is extremely sensitive to the SVD filtering process which in practice will be determined by the SNR in the WFS measurements. We also remark that such a naive approach also has small sky coverage as it requires a star which is bright enough to obtain high-order wavefront measurements. However, with hybrid LGS constellations (Sect. 5.3) it may be possible to obtain good levels of reconstruction. Taking the sky coverage values obtained with the three tip-tilt NGS scheme considered by Gemini (see Sect. 5.2) and assuming equivalent magnitude requirements for the tip-tilt NGS with hybrid LGS constellations it is expected that this approach can be applied to a larger fraction of the sky as only a single tip-tilt NGS is required. There is a lot of room for possible improvements to the technique used (Femenía 2003).

The optimal $`DM_2`$ conjugation altitude depends strongly on the guide star constellation used. When considering NGS constellations our MCAO simulations give results in disagreement with the theoretical work by Tokovinin & Viard (2001) which would place $`DM_2`$ at about 13 Km above telescope pupil (Devaney et al. 2002). The reason for such disagreement is believed to be caused by the strong assumptions in these theoretical works resulting in overestimating the $`DM_2`$ conjugation altitude. From our study based on detailed numerical simulations we consider that a suitable $`DM_2`$ conjugation altitude for a 10-m class telescope is around 8 Km and thus we are considering such an option for the design of the AO system at the GTC telescope which considers an initial AO system with a single corrector but upgradeable to a dual-conjugate system. It has been argued that the aforementioned theoretical works would be applicable to the case of extremely large telescopes (30-100 m diameter) as the main limiting assumption (i.e. full metapupil coverage by overlapping the GS beams at the reconstructed layers) is relaxed and it should be possible to fulfill the requirement of perfect wavefront reconstruction at each metapupil. However, one should consider that the larger the telescope diameter, the larger the number of modes required in order to reach the same Strehl Ratio as with a smaller telescope diameter. It follows from this consideration that with extremely large telescopes, although the metapupil coverage is nearly complete, the importance of small gaps in the metapupil coverage is much more important, this implying that the assumptions in the above theoretical works may not be realistic even for extremely large telescopes.

Our last conclusion regards a rule of thumb indicated by our simulations: given a fixed number of GS (either LGS or NGS) the best MCAO performance occurs for the constellation configuration that minimizes the distance from any point in the FoV to a GS position. As we have seen although 5 GS on the vertices of a pentagon give a larger covered surface of the metapupil than 4 GS on the square vertices plus one on-axis GS, the latter configuration yields the better reconstruction.

We would like to end with a comment on the turbulence profile used in our simulations. It has been derived from only 6 turbulence profiles at the ORM obtained with balloon flights which is a rather small sample of turbulence profiles. A larger set of turbulence profiles would be desirable in order to obtain a more significant average profile at ORM. In any case, all the simulations were also conducted on a modified version of the average profile and the optimal conjugation altitude of $`DM_2`$ remained essentially the same.

###### Acknowledgements.

The Guide Star Catalogue-II is a joint project of the Space Telescope Science Institute and 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. The participation of the Osservatorio Astronomico di Torino is supported by the Italian Council for Research in Astronomy. Additional support is provided by European Southern Observatory, Space Telescope European Coordinating Facility, the International GEMINI project and the European Space Agency Astrophysics Division.


format: text

text:
# The Bose-Einstein correlation function 𝐶₂(𝑄) from a Quantum Field Theory point of view

\[

## Abstract

We show that a recently proposed derivation of Bose-Einstein correlations (BEC) by means of a specific version of thermal Quantum Field Theory (QFT), supplemented by operator-field evolution of the Langevin type, allows for a deeper understanding of the possible coherent behaviour of the emitting source and a clear identification of the origin of the observed shape of the BEC function $`C_2(Q)`$. Previous conjectures in this matter obtained by other approaches are confirmed and have received complementary explanation.

\]

In this work we would like to focus attention on two specific features of Bose-Einstein correlations (BEC) clearly visible when BEC are presented in the language of some specific (thermal) version of Quantum Field Theory (QFT) supplemented by the operator-field evolution of Langevin type proposed recently . Because the importance of BEC and their present experimental and theoretical status are widely known and well documented (see, for example,  and references therein), we shall not repeat it here. The features we shall discuss are: $`(i)`$ how the possible coherence of the hadronizing system (modelled here by some external stationary force occurring in the Langevin equations describing hadronization process ) influences the BEC function $`C_2(Q)`$ and $`(ii)`$ what is the true origin of the experimentally observed $`Q`$-dependence of the $`C_2(Q)`$ correlation function in the approach used here .

In what concerns the first point we have obtained identical expression as derived in  by means of a quantum statistical (QS) approach with novel interpretation of the chaoticity parameter $`p`$ introduced there. On the other hand, our results (and that of ) differ from the formula recently obtained in . We shall argue that the origin of this difference lies in different ways of introducing the concept of coherence in both approaches.

In the second point we demonstrate that in order to obtain a given (experimentally observed) shape of the BEC correlation function $`C_2(Q)`$ (i.e., its $`Q`$-dependence, where $`Q=|k_\mu k_\mu ^{}|=\sqrt{(k_\mu k_\mu ^{})^2}`$) one has to account somehow for the finiteness of the space-time region of the particle production (i.e., of the hadronizing source). In our approach this means the necessity of smearing out of some generalized functions (delta functions: $`\delta (Q_\mu =k_\mu k_\mu ^{})`$) appearing in the definition of thermal averages of some operators, which is characteristic feature of QFT approach used here. The freedom in using different types of smearing functions to perform such a procedure allows us to account for all possible different shapes of hadronizing sources apparently observed by experiment. .

Referring to  for details, let us recapitulate here the main points of our approach. The collision process produces a lot of particles out of which we select one (we assume for simplicity that we are dealing only with identical bosons) and describe it by operator $`b(\stackrel{}{k},t)`$ (the notation is the usual one: $`b(\stackrel{}{k},t)`$ is an annihilation operator, $`\stackrel{}{k}`$ is $`3`$-momentum and $`t`$ is a real time). The rest of the particles are then assumed to form a kind of heat bath, which remains in equilibrium characterized by a temperature $`T=1/\beta `$ (which will be one of our parameters). All averages $`(\mathrm{})`$ are therefore thermal averages of the type: $`(\mathrm{})=Tr\left[(\mathrm{})e^{\beta H}\right]/Tr\left(e^{\beta H}\right)`$. We shall also allow for some external (to the above heat bath) influence to our system. Therefore we shall represent the operator $`b(\stackrel{}{k},t)`$ as consisting of a part corresponding to the action of the heat bath, $`a(\stackrel{}{k},t)`$, and also of a part describing action of these external factors, $`R(\stackrel{}{k},t)`$:

$$b(\stackrel{}{k},t)=a(\stackrel{}{k},t)+R(\stackrel{}{k},t).$$
(1)
The time evolution of such a system is then assumed to be given by a Langevin equation

$$i_tb(\stackrel{}{k},t)=F(\stackrel{}{k},t)A(\stackrel{}{k},t)+P$$
(2)
(and a similar conjugate equation for $`b^+(\stackrel{}{k},t)`$). These equations are supposed to model all aspects of the hadronization process. The combination $`F(\stackrel{}{k},t)A(\stackrel{}{k},t)`$ represents the so called Langevin force and is therefore responsible for the internal dynamics of hadronization in the following manner: $`A`$ is related to stochastic dissipative forces and is given by

$$A(\stackrel{}{k},t)=_{\mathrm{}}^+\mathrm{}𝑑\tau K(\stackrel{}{k},t\tau )b(\stackrel{}{k},\tau ),$$
(3)
with the operator $`K(\stackrel{}{k},t)`$ being a random evolution field operator describing the random noise and satisfying the usual correlation-fluctuation relation for the Gaussian noise . The operator $`F(\stackrel{}{k},t)`$ describes the influence of heat bath,

$$F(\stackrel{}{k},t)=_{\mathrm{}}^+\mathrm{}\frac{d\omega }{2\pi }\psi (k_\mu )\widehat{c}(k_\mu )e^{i\omega t}.$$
(4)
Our heat bath is represented by an ensemble of damped oscillators, each described by operator $`\widehat{c}(k_\mu )`$ such that $`[\widehat{c}(k_\mu ),\widehat{c}^+(k_\mu ^{})]=\delta ^4(k_\mu k_\mu ^{})`$, and characterized by some function $`\psi (k_\mu )`$ . Finally, the constant term $`P`$ (representing external source term in Langevin equation) denotes the possible influence of some external force (assumed here to be constant in time). This force would result, for example, in a strong ordering of phases, leading therefore to the coherence effect in the sense discussed in . Out of many details (for which we refer to ) what is important in our case is the fact that the $`2`$-particle correlation function for like-charge particles, is defined as $`(k_\mu =(\omega =k^0,k_j)`$):

$`C_2(Q)`$ $`=`$ $`\xi (N){\displaystyle \frac{\stackrel{~}{f}(k_\mu ,k_\mu ^{})}{\stackrel{~}{f}(k_\mu )\stackrel{~}{f}(k_\mu ^{})}}`$ (5)
$`=`$ $`\xi (N)\left[1+D(k_\mu ,k_\mu ^{})\right],`$ (6)
where $`\stackrel{~}{f}(k_\mu ,k_\mu ^{})=\stackrel{~}{b}^+(k_\mu )\stackrel{~}{b}^+(k_\mu ^{})\stackrel{~}{b}(k_\mu )\stackrel{~}{b}(k_\mu ^{})`$ and $`\stackrel{~}{f}(k_\mu )=\stackrel{~}{b}^+(k_\mu )\stackrel{~}{b}(k_\mu )`$ are the corresponding thermal statistical averages (in which temperature $`T`$ enters as a parameter) with $`\stackrel{~}{b}(k_\mu )=\stackrel{~}{a}(k_\mu )+\stackrel{~}{R}(k_\mu )`$ being the corresponding Fourier transformed stationary solution of eq. (2). As shown in  (notice that operators $`\stackrel{~}{R}(k_\mu )`$ by definition commute with themselves and with any other operator considered here):

$`\stackrel{~}{f}(k_\mu ,k_\mu ^{})`$ $`=`$ $`\stackrel{~}{f}(k_\mu )\stackrel{~}{f}(k_\mu ^{})+`$ (7)
$`+`$ $`\stackrel{~}{a}^+(k_\mu )\stackrel{~}{a}(k_\mu ^{})\stackrel{~}{a}^+(k_\mu ^{})\stackrel{~}{a}(k_\mu )+`$ (8)
$`+`$ $`\stackrel{~}{a}^+(k_\mu )\stackrel{~}{a}(k_\mu ^{})\stackrel{~}{R}^+(k_\mu ^{})\stackrel{~}{R}(k_\mu )+`$ (9)
$`+`$ $`\stackrel{~}{a}^+(k_\mu ^{})\stackrel{~}{a}(k_\mu )\stackrel{~}{R}^+(k_\mu )\stackrel{~}{R}(k_\mu ^{}),`$ (10)
$`\stackrel{~}{f}(k_\mu )`$ $`=`$ $`\stackrel{~}{a}^+(k_\mu )\stackrel{~}{a}(k_\mu )+|\stackrel{~}{R}(k_\mu )|^2.`$ (11)
This defines $`D(k_\mu ,k_\mu ^{})=\stackrel{~}{f}(k_\mu ,k_\mu ^{})/[\stackrel{~}{f}(k_\mu )\stackrel{~}{f}(k_\mu ^{})]1`$ in (6) in terms of the operators $`\stackrel{~}{a}(k_\mu )`$ and $`\stackrel{~}{R}(k_\mu )`$, which in our case are equal to:

$$\stackrel{~}{a}(k_\mu )=\frac{\stackrel{~}{F}(k_\mu )}{\stackrel{~}{K}(k_\mu )\omega }\mathrm{and}\stackrel{~}{R}(k_\mu )=\frac{P}{\stackrel{~}{K}(k_\mu )\omega }.$$
(12)
The multiplicity $`N`$ depending factor $`\xi `$ is in our case equal to $`\xi (N)=N^2/N(N1)`$. This means therefore that the correlation function $`C_2(Q)`$, as defined by eq. (6), is essentially given in terms of $`P`$ and the two following thermal averages for the $`F(\stackrel{}{k},t)`$ operators:

$`F^+(\stackrel{}{k},t)F(\stackrel{}{k}^{},t^{})`$ $`=`$ $`\delta ^3(\stackrel{}{k}\stackrel{}{k}^{})`$ (13)
$``$ $`{\displaystyle \frac{d\omega }{2\pi }\left|\psi \right|^2n(\omega )e^{+i\omega (tt^{})}},`$ (14)
$`F(\stackrel{}{k},t)F^+(\stackrel{}{k}^{},t^{})`$ $`=`$ $`\delta ^3(\stackrel{}{k}\stackrel{}{k}^{})`$ (15)
$``$ $`{\displaystyle \frac{d\omega }{2\pi }\left|\psi \right|^2[1+n(\omega )]e^{i\omega (tt^{})}}`$ (16)
where $`n(\omega )=\left\{\mathrm{exp}\left[(\omega \mu )\beta \right]1\right\}^1`$ is the number of (by assumption - only bosonic in our case) damped oscillators of energy $`\omega `$ in our reservoir characterized by parameters $`\mu `$ (chemical potential) and inverse temperature $`\beta =1/T`$ (both being free parameters). Notice that with only delta functions present in (14) one would have a situation in which our hadronizing system would be described by some kind of white noise only. The integrals multiplying these delta functions and depending on $`(a)`$ momentum characteristic of our heat bath $`\psi (k_\mu )`$ (representing in our case, by definition, the hadronizing system) and $`(b)`$ assumed bosonic statistics of produced secondaries resulting in factors $`n(\omega )`$ and $`1+n(\omega )`$, respectively, bring the description of our system closer to reality.

It should be stressed at this point that, contrary to the majority of discussions of BEC , we are working here directly in phase space , so far no space-time considerations were used. It is easy to realize now that the existence of BEC, i.e., the fact that $`C_2(Q)>1`$, is strictly connected with nonzero values of the thermal averages (14). However, in the form presented there, they differ from zero only at one point, namely for $`Q=0`$ (i.e., for $`k_\mu =k_\mu ^{}`$). Actually, this is the price one pays for the QFT assumptions tacitly made here, namely for the infinite spatial extension and for the uniformity of our reservoir. But we know from the experiment  that $`C_2(Q)`$ reaches its maximum at $`Q=0`$ and falls down towards its asymptotic value of $`C_2=1`$ at large of $`Q`$ (actually already at $`Q1`$ GeV/c). To reproduce the same behaviour by means of our approach here, one has to replace delta functions in eq. (14) by functions with supports larger than limited to a one point only. This means that such functions should not be infinite at $`Q_\mu =k_\mu k_\mu ^{}=0`$ but remain more or less sharply peaked at this point, otherwise remaining finite and falling to zero at small, but finite, values of $`|Q_\mu |`$ (actually the same as those at which $`C_2(Q)`$ reaches unity):

$$\delta (k_\mu k_\mu ^{})\mathrm{\Omega }_0\sqrt{\mathrm{\Omega }(q=Qr)}.$$
(18)
Here $`\mathrm{\Omega }_0`$ has the same dimension as the $`\delta `$ function (actually, it is nothing else but $`4`$-dimensional volume restricting the space-time region of particle production) and $`\mathrm{\Omega }(q)`$ is a dimensionless smearing function which contains the $`q`$-dependence we shall be interested in here. In this way we are tacitly introducing a new parameter, $`r_\mu `$, a $`4`$-vector such that $`\sqrt{(r_\mu )^2}`$ has dimension of length and which makes the product $`Qr=Q_\mu r_\mu =q`$ dimensionless. This defines the region of nonvanishing density of oscillators $`\widehat{c}`$, which we shall identify with the space-time extensions of the hadronizing source. The expression (18) has to be understood in a symbolic sense, i.e., that $`\mathrm{\Omega }(Qr)`$ is a function which in the limit of $`r\mathrm{}`$ becomes strictly a $`\delta `$ function. Making such replacement in eq. (14) one must also decide how to accordingly adjust $`n(\omega )`$ occurring there because now, in general, $`\omega \omega ^{}`$. In what follows we shall simply replace $`n(\omega )n(\overline{\omega })`$ with $`\overline{\omega }=(\omega +\omega ^{})/2`$ (which, for classical particles would mean that $`n(\omega )\sqrt{n(\omega )n(\omega ^{})}`$).

In such way $`r`$ becomes new (and from the phenomenological point of view also the most important) parameter entering here together with the whole function $`\mathrm{\Omega }(Qr)`$, to be deduced from comparison with experimental data . With such a replacement one now has

$$D(k_\mu ,k_\mu ^{})=\frac{\sqrt{\stackrel{~}{\mathrm{\Omega }}(q)}}{(1+\alpha )(1+\alpha ^{})}\left[\sqrt{\stackrel{~}{\mathrm{\Omega }}(q)}+2\sqrt{\alpha \alpha ^{}}\right]$$
(19)
where

$$\stackrel{~}{\mathrm{\Omega }}(q)=\gamma \mathrm{\Omega }(q),\gamma =\frac{n^2(\overline{\omega })}{n(\omega )n(\omega ^{})},\alpha \frac{P^2}{|\psi (k_\mu )|^2n(\omega )},$$
(20)
with $`n(\omega )`$ the same as defined above. The parameter $`\alpha `$ is another very important parameter, which summarizes our knowledge of other than space-time characteristics of the hadronizing source (given by $`\mathrm{\Omega }(q)`$ introduced above). In particular it contains the external static force $`P`$ present in the evolution equation (2). It is combined (in multiplicative way) with information on the momentum dependence of the reservoir (via $`|\psi (k_\mu )|^2`$) and on the single particle distributions of the produced particles (via $`n(\omega =\mu _T\mathrm{cosh}y)`$ where $`\mu _T`$ and $`y`$ are, respectively, the transverse mass and rapidity). Notice that $`\alpha >0`$ only when $`P0`$. Actually, for $`\alpha =0`$ one has

$$1<C_2(Q)<1+\gamma \mathrm{\Omega }(Qr),$$
(21)
i.e., it is contained between limits corresponding to very large (lower limit) and very small (upper limit) values of $`P`$. Because of this $`\alpha `$ plays the role of the coherence parameter . For $`\gamma 1`$, neglecting the possible energy-momentum dependence of $`\alpha `$ and assuming that $`\alpha ^{}=\alpha `$ one gets the expression

$$C_2(Q)=1+\frac{2\alpha }{(1+\alpha )^2}\sqrt{\mathrm{\Omega }(q)}+\frac{1}{(1+\alpha )^2}\mathrm{\Omega }(q),$$
(22)
which is formally identical with what has been obtained in  by means of QS approach. It has precisely the same form, consisting two $`Q`$dependent terms containing the information on the shape of the source, one being the square of the other, each multiplied by some combination of the chaoticity parameter $`p=1/(1+\alpha )`$ (however, in  $`p`$ is defined as the ratio of the mean multiplicity of particles produced by the so called chaotic component of the source to the mean total multiplicity, $`p=N_{ch}/N`$). In fact, because in general $`\alpha \alpha ^{}`$ (due to the fact that $`\omega \omega ^{}`$ and therefore the number of states, identified here with the number of particles with given energy, $`n(\omega )`$, are also different) one should rather use the general form (6) for $`C_2`$ with details given by (19) and (20) and with $`\alpha `$ depending on such characteristics of the production process as temperature $`T`$ and chemical potential $`\mu `$ occurring in definition of $`n(\omega )`$.

Notice that eq. (22) differs from the usual empirical parameterization of $`C_2(Q)`$ ,

$$C_2(Q)=1+\lambda \mathrm{\Omega }(Qr),$$
(23)
with $`0<\lambda <1`$ being a free parameter adjusting the observed value of $`C_2(Q=0)`$, which is customary called ”incoherence”, and with $`\mathrm{\Omega }(Qr)`$ represented usually as Gaussian. Recently eq. (23) has found strong theoretical support expressed in great detail in . The natural question arises: which of the two formulas presented here is correct? The answer is: both are right in their own way. This is because each of them is based on different ways of defining coherence of the source. In  one uses the notion of coherently and chaotically produced particles or, in other words, one divides hadronizing source into coherent and chaotic subsources. In  one introduces instead the notion of partially coherent fields representing produced particles, i.e., one has only one source, which produces partially coherent fields. Our approach is similar as we describe our particle by operator $`b(\stackrel{}{k},t)`$, which consists of two parts, cf. eq. (1), one of which depends on the external static force $`P`$. The action of this force is to order phases of particles in our source (represented by the heat bath). The strength of this ordering depends on the value of the external force $`P`$. In any case, for $`P0`$, it demonstrates itself as a partial coherence .

Let us return to the problem of $`Q`$-dependence of BEC. One more remark is in order here. The problem with the $`\delta (k_\mu k_\mu ^{})`$ function encountered in two particle distributions does not exist in the single particle distributions, which are in our case given by eq. (11) and which can be written as $`\stackrel{~}{f}(k_\mu )\stackrel{~}{a}^+(k_\mu )\stackrel{~}{a}(k_\mu )+|\stackrel{~}{R}(k_\mu )|^2(1+\alpha )\stackrel{~}{a}^+(k_\mu )\stackrel{~}{a}(k_\mu )`$ . To be more precise

$$\stackrel{~}{f}(k_\mu )=(1+\alpha )\mathrm{\Xi }(k_\mu ,k_\mu ),$$
(24)
where $`\mathrm{\Xi }(k_\mu ,k_\mu )`$ is one-particle distribution function for the ”free” (undistorted) operator $`\stackrel{~}{a}(k_\mu )`$ equal to

$$\mathrm{\Xi }(k_\mu ,k_\mu )=\mathrm{\Omega }_0\left|\frac{\psi (k_\mu )}{\stackrel{~}{K}(k_\mu )\omega }\right|^2n(\omega ).$$
(25)
Notice that the actual shape of $`\stackrel{~}{f}(k_\mu )`$ is dictated both by $`n(\omega )=n(\omega ;T,\mu )`$ (calculated for fixed temperature $`T`$ and chemical potential $`\mu `$ at energy $`\widehat{\omega }`$ as given by the Fourier transform of random field operator $`\stackrel{~}{K}`$ and by shape of the reservoir in the momentum space provided by $`\psi (k_\mu )`$) and by external force $`P`$ in parameter $`\alpha `$. They are both unknown, but because these details do not enter the BEC function $`C_2(Q)`$, we shall not pursue this problem further. What is important for us at the moment is that both the coherent and the incoherent part of the source have the same energy-momentum dependence (whereas in other approaches mentioned here they were usually assumed to be different). On the other hand it is clear from (24) that $`N=N_{ch}+N_{coh}`$ (where $`N_{ch}`$ and $`N_{coh}`$ denote multiplicities of particles produced chaotically and coherently, respectively) therefore justifying definition of chaoticity $`p`$ mentioned above.

For an illustration we plotted in Fig. 1 the correlation function $`C_2(Q)`$ as given by eq. (22) for different choices of $`\mathrm{\Omega }(q)`$ corresponding to different hadronizing sources discussed in the literature  (here $`r=|r_\mu |=\sqrt{r_\mu r_\mu }`$):

* Gaussian: $`\mathrm{\Omega }(q)=\mathrm{exp}\left(Q^2r^2\right)`$;
* exponential: $`\mathrm{\Omega }(q)=\mathrm{exp}(Qr)`$;
* Lorentzian: $`\mathrm{\Omega }(q)=1/\left(1+Qr\right)^2`$;
* given by Bessel function :
$`\mathrm{\Omega }(q)=\left[J_1(Qr)/(Qr)\right]^2`$ .

All curves are drawn for the same values of the size parameter $`r=1`$ fm and assuming for simplicity constant and equal values of $`\alpha `$ and $`\alpha ^{}`$ parameters (i.e., using simplified eq. (22)), which have been put equal $`\alpha =0.2`$ here, just for illustrational purpose (it corresponds to $`p=0.8`$ in ). Fig. 2 shows in detail (using Gaussian shape of $`\mathrm{\Omega }(q)`$ function) the dependence of $`C_2(Q)`$ on different values of $`\alpha =0,0.25,1,4`$ (again, used in the same approximate way as before and corresponding to $`p=1.,0.8,0.5,0.2`$) and compare it to the case when the second term in eq. (22) is neglected, as is the case in majority of phenomenological fits to data.

To summarize: using a specific version of QFT supplemented by Langevin evolution equation (2) to describe hadronization process  we have derived the usual BEC correlation function in the form explicitly showing the origin of both the so called coherence (and how it influences the structure of BEC) and the $`Q`$-dependence of BEC represented by correlation function $`C_2(Q)`$. The dynamical source of coherence is identified in our case with the existence of a constant external term $`P`$ in the Langevin equation. Its influence turns out to be identical with the one obtained before in the QS approach  and is described by eq. (19). Its action is to order phases of the produced secondaries. Therefore for $`P\mathrm{}`$ we have all phases aligned in the same way and $`C_2(Q)=1`$. This is because both here and in  the coherence has already been introduced on the level of a hadronizing source, as property of fields (in ) or operators describing produced particles. Dividing instead the hadronizing source itself into coherent and chaotic subsources leads to results obtained in  and given by eq. (23). The controversy between results given by  and  is therefore explained: both approaches are right, one should only remember that they use different descriptions of the notion of coherence. It is therefore up to the experiment to decide which proposition is followed by nature: the simpler formula (23) or rather the more involved (6) together with (19). From Fig. 2 one can see that differences between both forms are clearly visible, especially for larger values of coherence $`\alpha `$, i.e., for lower chaoticity parameter $`p`$.

From our presentation it is also clear that the form of $`C_2`$ reflects distributions of the space-time separation between the two observed particles rather than the distribution of their separate production points  (i.e., it is Fourier transform of two-particle density profile of the hadronizing source, $`\rho (r_1,r_2)=\rho (r_1r_2)`$, without approximating it by the product of single-particle densities, as in ).

Finally, we would like to stress that our discussion is so far limited to only a single type of secondaries being produced. It is also aimed at a description of hadronization understood as kinetic freeze-out in some more detailed approaches. So far we were not interested in the other (highly model dependent) details of the particle production process. This is enough to obtain our general goals, i.e., to explain the possible dynamical origin of coherence in BEC and the origin of the specific shape of the correlation $`C_2(Q)`$ functions as seen from the QFT perspective. Actually, our source of coherence should be regarded as being only one possibility , the others were discussed in detail in . It is then plausible that in general description of the BEC effect they should be somehow combined, especially if experimental data would indicate such necessity. But to do so our approach should first be generalized to allowing, as is the case in , for production of different types of secondaries and allow also for resonance production and final state interactions (both of strong and Coulomb origin). This is, however, outside the scope of the present paper.

We would like to acknowledge support obtained from the Bogolyubov-Infeld program in JINR and partial support of the Polish State Committee for Scientific Research (KBN), grants 621/E-78/SPUB/CERN/P-03/DZ4/99 and 2P03B05724.


text:
# Coulomb dissociation of 8B and the low-energy cross section of the 7Be(p,𝛾)8B solar fusion reaction

(June 3, 2003)

## Abstract

An exclusive measurement of the Coulomb breakup of 8B into 7Be+p at 254 $`A`$ MeV allowed to study the angular correlations of the breakup particles. These correlations demonstrate clearly that E1 multipolarity dominates and that E2 multipolarity can be neglected. By using a simple single-particle model for 8B and treating the breakup in first-order perturbation theory, we extract a zero-energy $`S`$ factor of $`S_{17}(0)=18.6\pm 1.2\pm 1.0`$ eV b, where the first error is experimental and the second one reflects the theoretical uncertainty in the extrapolation.

Exciting new results Ahm01  from the Sudbury Neutrino Observatory (SNO) have proven for the first time that the measured high-energy neutrino flux from the Sun agrees well with the one calculated from standard solar models BPB01 ; Bru99  if non-electron flavour neutrinos are taken into account. This again focusses attention onto the 7Be(p,$`\gamma `$)8B reaction which provides almost exclusively the high-energy neutrinos measured in the SNO experiment. Their flux depends linearily on the 7Be(p,$`\gamma `$)8B cross section at solar energies. Very recently, the latter has been redetermined by new high-precision direct measurements Ham01 ; Str01 ; Jun02 ; Bab03  and extrapolated to zero energy with the help of a theoretical model Des94 . The resulting zero-energy astrophysical $`S`$ factors, $`S_{17}(0)`$, however, do not always agree within their quoted errors: Hammache et al. Ham01  found $`S_{17}(0)=18.8\pm 1.7`$ eV b, in agreement with other direct-capture data Fil83 ; Vau70 ; Str01 . In contrast, Junghans et al. Jun02  report a considerably larger value, $`S_{17}(0)=22.3\pm 0.7\pm 0.5`$ eV b. The very recent result of Baby et al. Bab03  also favours a rather large value of $`S_{17}(0)=21.2\pm 0.7`$ eV b.

In view of their importance for astro- and elementary-particle physics, these conflicting results should be verified and cross-checked by other, indirect measurements that have different systematic errors. One possibility is Coulomb dissociation (CD) of 8B in the electromagnetic field of a high-$`Z`$ nucleus. Such measurements have been performed at low Kol01 , intermediate Kik98 ; Dav01 , and high energies Iwa99 . Alternatively, $`S_{17}`$(0) can also be calculated from asymptotic normalization coefficients (ANC) which in turn are determined in low-energy proton-transfer or in proton-removal reactions Azh01 ; Tra01 ; Cor03 .

In the present Letter, we focus on a crucial question that must be answered if one wants to use the CD method to derive a precise value for $`S_{17}`$(0). The astrophysical $`S`$ factors of the 7Be(p,$`\gamma `$) reaction can only be calculated reliably from the energy-differential CD cross sections if the electromagnetic multipole components relevant for direct capture and the time-reversed process have the same strength. In low-energy proton capture the E1 contribution by far dominates the cross section. While E1 is the dominant multipolarity also in CD, one can show easily that the equivalent photon field emitted from a high-$`Z`$ target nucleus contains a strong E2 component. This is particularly true for CD at low energies. At higher energies (see Ref. Iwa99 ) the relative amount of E2 multipolarity is expected to be reduced, but may still be substantial enough to affect the final result. To remove this ambiguity, it is indispensable to either determine the E1/E2 ratio in CD experimentally, or to extract $`S_{17}`$ with such cuts that any E2 contribution is negligible.

Experimental limits for a possible E2 contribution were extracted in the work of Kikuchi et al. Kik98  and Iwasa et al. Iwa99 . Both papers found negligible E2 contributions. Recently, Davids et al. have reported positive experimental evidence for a finite E2 contribution in CD of 8B, mainly from the analysis of inclusive longitudinal momentum ($`p_{||}`$) spectra of 7Be fragments measured at 44 and 81 $`A`$ MeV Dav01 . The asymmetries in the $`p_{||}`$ spectra were interpreted to be due to E1-E2 interference in terms of first-order perturbation-theory Esb96 .

In order to resolve these discrepancies, we decided to perform an exclusive CD experiment at high energy (254 $`A`$ MeV) at the kaon spectrometer KaoS at GSI Sen93  with the aim to measure quantities that should be sensitive to contributions of E2 multipolarity, namely the angular correlations of the 8B-breakup particles, proton and 7Be. Experimentally, this requires high-resolution measurements of the positions and angles of the incident 8B beam as well as those of the breakup fragments. The 8B secondary beam was produced at the SIS/FRS radioactive beam facility at GSI Gei92  by fragmenting a 350 $`A`$ MeV 12C beam in a 8 g/cm2 Be target and separating it from contaminant ions in a 1.4 g/cm2 wedge-shaped Al degrader placed in the FRS intermediate focal plane. Typical 8B beam intensities in front of KaoS were $`5\times 10^4`$ per 4 sec spill; the only contaminant consisted of about 20% 7Be ions which could be identified event by event with the help of a time-of-flight measurement.

Positions and angles of the secondary beam incident on the Pb breakup target were measured with the help of two parallel-plate avalanche counters (PPAC) located at 308.5 cm and 71 cm upstream from the target, respectively. The detectors, which were designed at RIKEN Kum01 , had areas of $`10\times 10\mathrm{cm}^2`$ and allowed to track the incident 8B beam with about 90% efficiency and with position and angular resolutions of 1.3 mm and 1 mrad, respectively. Downstream from the Pb target (which consisted of 50 mg/cm2 208Pb enriched to 99.0$`\pm `$0.1%), the angles and positions as well as the energy losses of the outgoing particles were measured with two pairs of Si strip detectors (300 $`\mu `$m thick, 100 $`\mu `$m pitch) located at distances of about 14 cm and 31 cm. Proton and 7Be momenta were analyzed with the KaoS spectrometer which was set up almost identical

to our previous experiment Iwa99 , except for a newly constructed plastic-scintillator wall near the KaoS focal plane with 30 elements (each 7 cm wide and 2 cm thick) used for trigger purposes.

The coincident p and 7Be signals resulting from breakup in the 208Pb target were identified by reconstructing their vertex at the target, this removed all breakup events in layers of matter other than the target. The measured momentum vectors of the outgoing p and 7Be particles allowed to construct the invariant-mass spectrum of the excited 8B system prior to breakup.

Fig. 1 shows the coordinate systems used. The angle $`\theta _8`$ is the laboratory scattering angle of 8B relative to the incoming 8B beam. The polar angles, $`\theta _{cm}`$, and the azimuthal angles, $`\varphi _{cm}`$, of the breakup protons are measured in the rest frame of the 8B\* system, as shown in Fig. 1. In the same way, one can calculate e.g. the transverse proton momentum vector in the reaction plane ($`p_t^{in}`$).

In the following we will present some angular distributions of the emitted proton in the frame of the 8B system that can be shown to be sensitive to an E2 amplitude in CD. To interprete the measured distributions we need guidance by a theoretical model. To this end, we have performed standard first-order perturbation-theory (PT) calculations of the CD process in the semi-classical approach Ber94 ; Typ02 , using a simple Woods-Saxon potential model for 8B. The potential depth for the ground state of 8B was adjusted to match the proton binding energy; the potential depths of the scattering states were fitted to the scattering lengths of the 7Li+n mirror system Koe83 . We used a radius parameter of $`r_0=1.25`$ fm and a diffuseness of $`a=0.65`$ fm. For channel spin $`I=2`$ (the dominant contribution) we obtained a potential depth of $`V_2=52.60`$ MeV. The resulting scattering length for this channel of $`a_{02}^{theo}=8`$ fm agrees well with the recently measured value of $`a_{02}^{exp}=7\pm 3`$ fm (Angulo et al. Ang03 ).

To take into account absorption due to nuclear overlap in CD, we have introduced a diffuse absorptive nuclear potential with a depth of 20 MeV and a radius of 9.91 fm, i.e. the sum of the projectile and target radii. This choice reproduces well the integral $`\theta _8`$ angular distribution. Technically, the results of the PT calculations were returned as a statistically-distributed ensemble of $`\mathrm{500\hspace{0.17em}000}`$ CD-“events” that were analyzed in the same way as the experimental data, thus imposing the experimental cuts.

We first present in Fig. 2 the distribution of $`p_t^{in}`$ for three different upper limits in $`\theta _8`$, 0.62, 1.0, and 2.5. In classical Rutherford scattering, this corresponds to impact parameters of 30 fm, 18.5 fm, and 7 fm, respectively. Relative energies between p and 7Be up to 1.5 MeV were selected. The experimental data for all three $`\theta _8`$-cuts can be reproduced well by a PT calculation that includes only E1 multipolarity (full histograms in Fig. 2, the theoretical curves were normalized individually to the data points). If E1-plus-E2 multipolarity is used in the PT calculation, the different impact-parameter dependences of E1 and E2 multipolarity lead to markedly different shapes for the different $`\theta _8`$-cuts (dashed histograms in Fig. 2). The latter distributions are, however, in clear disagreement with our data points.

Fig. 3 depicts the experimental $`\varphi _{cm}`$ and $`\theta _{cm}`$ distributions for three different $`E_{rel}`$ bins,

as indicated in the figure. A “safe” $`\theta _8`$ limit of 1 was chosen. As expected, these distributions are mostly isotropic at low $`E_{rel}`$ (indicative of $`s`$-waves) and become increasingly anisotropic for larger values (contributions from $`d`$-waves). For the $`\varphi _{cm}`$ distributions, which are most sensistive to E2 admixtures, the calculations for pure E1 multipolarity clearly fit best; inclusion of an E2 component shifts the maxima away from 90 and 270 with increasing $`E_{rel}`$, while at the same time the anisotropy is reduced. Similar conclusions can be drawn from the bottom part of Fig. 3, where the proton polar angular ($`\theta _{cm}`$) distributions are shown. The low-$`E_{rel}`$ bins show little sensitivity to E2 multipolarity, whereas inclusion of E2 leads to a marked discrepancy near $`\mathrm{cos}(\theta _{cm})=1`$ for the highest $`E_{rel}`$ bin. More detailed calculations show that at most E2 amplitudes of $``$0.3 times the theoretical one from our simple model are simultaneously compatible with all our measured observables. Since this would correspond to E2-contributions to the cross sections of less than 1%, much less than the errors of the data points, we neglect the effect of E2 multipolarity. This is in line with conclusions drawn by Kikuchi et al. Kik98  and by Iwasa et al. Iwa99  from their respective $`\theta _8`$ distributions (which are, however, less sensitive to a small E2 component than the present angular correlations). Our findings contradict the conclusions of Davids et al. Dav01  that a substantial E2 cross section has to be subtracted from the total measured CD cross section.

Our results allow to interpret the relative-energy distributions of the breakup particles in an easy way. In the following, we have restricted the angles $`\theta _8`$ to values below 0.62 to ensure both dominance of CD and reduction of the effect of any possible E2 contribution. The data are compared to a simulation with GEANT that includes two electromagnetic multipole components: a resonant M1 contribution located at $`E_{rel}`$=0.63 MeV with resonance parameters taken from Filippone et al. Fil83 , and the non-resonant E1 contribution from our theoretical model as described above. The latter was scaled by a normalization factor of 0.79. Note that we have added to the GEANT simulation a contribution that feeds the first excited state in 7Be at 429 keV using the measurements of Kikuchi et al. Kik98 .

Subtracting the small M1 contribution (that affects only a narrow $`E_{rel}`$ region around the resonance), the remaining $`d\sigma /dE_{rel}`$ distribution can be converted to the E1 astrophysical $`S`$ factor $`S_{17}(E_{rel})`$.

The resulting $`S_{17}`$ factors (averaged over $`E_{rel}`$ bins 0.2 to 0.3 MeV wide) are visualized in Fig. 4. The error bars do not include a common systematic error of 5.6%. The top panel (a) compares our results to those of other CD experiments Kik98 ; Iwa99 ; Dav01  (the data of Ref.Dav01  represent their E1-$`S_{17}`$ factors after subtraction of the E2-contribution). At low $`E_{rel}`$, the CD $`S`$ factors are in good agreement, though the Davids et al. Dav01  data are systematically lower. The bottom panel (b) compares our data to those of the 7Be(p,$`\gamma `$)8B measurements where the authors have subtracted the contribution from the M1 resonance (Refs. Ham01 ; Jun02 ; Bab03 ). At low energies the (p,$`\gamma `$) data of Refs. Ham01 ; Bab03  and ours are in good agreement, whereas the Seattle data Jun02  deviate considerably. The opposite behaviour is noted above the M1 resonance: our data and those of Refs. Jun02 ; Bab03  match excellently, whereas the other (p,$`\gamma `$) experiments Ham01 ; Str01 ; Vau70 ; Fil83  consistently report lower values. We want to emphasize the remarkably good agreement of our CD data up to 1.1 MeV with the most recent direct-proton-capture experiment where an ion-implanted 7Be target was used Bab03 .

To extrapolate to zero energy, all recent (p,$`\gamma `$) experiments have chosen the cluster model of Descouvemont and Baye Des94 . When we fit our data points up to $`E_{rel}`$ = 1.5 MeV to this model and add in quadrature a common systematic error of 5.6%, we obtain $`S_{17}(0)=20.8\pm 1.3`$ eV b (dashed lines in Fig. 4). Restricting the fit to energies below 0.6 MeV, where the model-dependence has been shown to be weaker Jen98 , $`S_{17}(0)=19.6\pm 1.4`$ eV b is obtained. Our potential model, however, reproduces the data over the entire energy range up to 1.5 MeV, yielding $`S_{17}(0)=18.6\pm 1.2`$ eV b (full lines in Fig. 4). It is interesting to note that a fit of the Baby et al. (p,$`\gamma `$) data to our model yields practically the same result, $`S_{17}(0)=18.1\pm 0.3`$ eV b. Clearly, still more high-precision experimental data are needed to resolve the discrepancies between the experimental data sets and to pin down the correct theoretical extrapolation of the measured data to solar energy. In the meantime, an additional “extrapolation error” of $`\pm `$ 1.0 eV b seems appropriate.

We conclude that Coulomb dissociation has been proven to be a valuable method to provide a rather precise value for the low-energy 7Be(p,$`\gamma `$) cross section. Since in CD all energy bins are measured simultaneously, CD provides a reliable measurement of the shape of the $`S_{17}`$ distribution. By setting tight constraints to the scattering angle $`\theta _8`$ and analyzing p-7Be angular correlations, a significant contribution from E2 multipolarity can be excluded. Small modifications of the Woods-Saxon potential parameters allow to reproduce the data in first-order perturbation theory with remarkable accuracy up to about $`E_{rel}=1.5`$ MeV.

The authors wish to thank K.-H. Behr, K. Burkard, and A. Brünle for technical assistance. Vivid discussions with B. Davids, P. Descouvemont, M. Hass, and A. Junghans are gratefully acknowledged.