2 Outline of methodology

2.1 Outline of the SCIDAR technique

The SCIDAR technique is based upon calculating the average correlation of a large number of short exposure images of the scintillation pattern within a telescope pupil produced by a binary star.

  

Figure 1: A simplified one-dimensional representation of the SCIDAR technique for a single turbulent layer at height $h_{\rm L}$. a) A large number of short exposure images of the scintillation pattern due to a binary star (of angular separation $\theta$) are collected. b) The average auto-correlation of these data is calculated, from which the altitude of the turbulent layer above the telescope pupil is given simply by $h_{\rm L}=r_{\rm L}/\theta$

  

Figure 2: Schematic diagram of the optical system for SCIDAR. a) shows the optical system for pupil plane SCIDAR. In this mode the telescope entrance pupil is imaged by the telescope mirror system L₁ and the field lens L₂. The focal length of the lens L₂ is chosen to give the correct image size at the detector plane. The image is effectively at infinity so that the final pupil image is formed at approximately the focal distance f₂ behind L₂. b) shows the optical system for generalised SCIDAR. At the CCD plane this system forms a pupil image from a plane located the required distance behind the telescope mirror ($\ell_1$). For this case then, the telescope mirrors (L₁) form a virtual image of the plane of focus (P₁) at a distance $\ell'_1$ behind the telescope mirror. The new field lens L'₂ re-images this virtual object to form the pupil image at the CCD plane. For example, for the JKT at La Palma; the telescope has a 1 m pupil with a F/15 focal ratio. A 159 mm field lens L₂ is required to give the correct image size for pupil plane SCIDAR. For Generalised SCIDAR the field lens L₂ is replaced by a lens with a focal length of 46 mm (L'₂), still located at 159 mm from the CCD plane. A virtual image of a plane at 3.5 km is formed $\ell_2$ $\simeq$ 64 mm in front of L'₂ and this virtual object is imaged by L'₂ onto the CCD detector

Figure 1 shows that in a simplified one-dimensional representation of the observing procedure for a binary object observed at zenith (a), a peak in the correlation plane occurs at a pupil separation $r_{\rm L}$, (b), and thus by simple geometry, denotes the presence of a single turbulent layer at an altitude $h_{\rm L}=r_{\rm L}/\theta$ above the telescope pupil. This simple situation can be generalised to the case where multiple turbulent layers exist, each one having an associated altitude calculable from the peak position within the average correlation. Also for the general case where the binary star is observed at some zenith angle z, then the geometrical conversion above becomes $h=r/(\theta {\rm sec}(z))$.

In the conventional operation of SCIDAR, the conjugate plane of the imaging system is placed at the altitude of the pupil of the telescope - referred to here as pupil-plane SCIDAR. However, low altitude turbulence cannot be determined in this mode since the scintillation requires a large distance over which to develop a signal strong enough to be detected. Consequently the correlation signal close to zero altitude is occupied only by a central core corresponding to the sum of all other higher altitude peaks. This is shown diagrammatically as the central peak in Fig. 1a. Recently, however, a solution to this problem was proposed by Fuchs et al. (1994) in which one changes the conjugate height of the imaging system from that of the telescope pupil, to one some distance below it. In effect, this creates an extra distance for the scintillation pattern to develop to a detectable level and also shifts the central core to below the pupil (zero) altitude. Experiments performed by us in this generalised SCIDAR mode appear to confirm the validity of this approach and the results are discussed later in Sect. 4.

2.2 Experimental setup

The optical system used to acquire the SCIDAR data is shown schematically in Fig. 2. The imaging system consists of a GenIII image intensifier lens-coupled to a high speed 128$\times$128 pixel CCD camera (DALSA CA-D1-0128A). The lens-coupling demagnifies, by a factor of approximately 5.2, the intensifier output phosphor image onto the 2 mm$\times$2 mm (16 $\mu$m pixel) CCD chip. In Fig. 2a, the pupil plane mode, the telescope entrance pupil is imaged onto the intensifier by the field lens (L₂). The focal length of this lens is chosen such that the final CCD pixel size corresponds to approximately 1 cm$\times$1 cm at the telescope entrance pupil. This choice of pixel sampling is made as a compromise between the received intensity per pixel and the detectability of atmospheric scintillation (see Sect. 3.1).

  

Figure 3: Example average scintillation pattern image obtained in the generalised SCIDAR mode at the JKT, in La Palma from 18000 1.6 ms exposure frames. The scale is in relative intensity units

For example, for the Jacobus Kapteyn Telescope (JKT) in La Palma, the telescope has a 1m pupil with an F/15 focal ratio and therefore an ideal focal length for the field lens (L₂), to give this required sampling, would be 159 mm. This would produce a 10.6 mm diameter pupil image at the intensifier, or equivalently a 2.04 mm pupil image on the CCD.

For the generalised mode shown in Fig. 2b, the field lens (L₂) is replaced by one of a different focal length, L'₂. The separation of field lens position to intensifier is unchanged however and thus the pupil plane sampling remains unchanged also. Using the thin lens equation one can calculate the required focal length for this generalised field lens for any chosen height of generalisation. For our measurements we are primarily concerned with ascertaining the strength and velocity of all turbulence present above the telescope altitude by conjugating to altitudes below that of the telescope pupil. Therefore for example, again for the JKT, in order to conjugate to a distance of 3.5 km, say, then a field lens of 46mm focal length is required.

The digital output from the DALSA Camera is acquired by a Bitflow Raptor PCI frame capture board controlled by a fast PC and the data is stored directly onto Digital Linear Tape. The frame exposure times normally used are 1.6 ms and 2.7 ms, with frame acquisition rates of 305 and 188 frames per second respectively.

Figure 3 shows an example of the average scintillation pattern of 18000 frames obtained using this system on the JKT whilst running in the generalised SCIDAR mode. The sheared image of the telescope pupil is clearly visible.

2.3 Data analysis

2.3.1 Determination of C_n²(h)

Prior to calculation of the C_n²(h) profile, each $i^{\rm th}$scintillation pattern image, I_i(x,y), is first mean-normalised to obtain the relative intensity,

$${I'}_{i}(x,y) = \frac{I_{i}(x,y) - < I(x,y) \gt}{< I(x,y) \gt},$$ (9)

(where < I(x,y) > is the average image) and the average auto-correlation of all such images calculated. To prevent the suppression of signal within the correlation plane for large separations within the pupil, this is then divided by the auto-correlation of an artificial pupil function (i.e. the optical transfer function of the pupil). This step is particularly important if using data collected from a telescope with a large obstruction ratio, as is the case for the example shown in Fig. 3 (see Sect. 3.1). The resulting function thus represents the two-dimensional intensity covariance, $C_I({\underline r})$.

  

Figure 4: Example a) spatial and b) temporal correlation plots generated from the same data set obtained from observations of the binary star $\gamma$ Aries (separation 7.8'' and components of equal magnitude 4.8 mv)

  

Figure 5: a) An example 2-d average cross-correlation plot and b) the corresponding Cn2(h) profile calculated using a maximum entropy based method. Note that the apparent single turbulent layer centred at 7.75 km in a) has been resolved into two layers in b). The resolution of the inversion calculation is 300 m

Figure 4a shows an example of this obtained using the generalised SCIDAR mode at the JKT for an effective conjugate distance of 2.5 km, which clearly indicates the presence of at least two distinct turbulent layers. Each of these can be identified by their individual ``triple-peak" pattern in a direction parallel to that of the binary star; comprising of a central peak contribution plus a (symmetrically-reflected) peak at a distance from the centre corresponding to the height of the layer above the altitude of conjugation via $h=r/(\theta {\rm sec}(z))$.

The finite width of these peaks results directly from the Fresnel length associated with each turbulent layer, given by $\Delta h_{\rm F}(h) = 0.7 \sqrt{\lambda h}/\theta$ (Vernin & Azouit 1983). This in turn effectively limits the degree of resolution in C_n²(h) attainable from the SCIDAR data. For the data collected thus far, $\Delta h_{\rm F}(h)$ ranges from around 700 m (for a pupil-conjugate distance of say 3 km) to around 2000 m (for the maximum altitudes considered). Assuming a pupil plane sampling, $\Delta r$, of 1 cm pix^-1 and binary star separations between 4-9'', then the altitude sampling available, given by $\Delta h = \Delta r/(\theta {\rm sec}(z))$, corresponds to around 200-500 m, and thus represents a reasonable sampling rate of the Fresnel length at all altitudes.

A slice taken through the correlation plane, parallel to the direction of the binary star, $C_{\rm I}(r,\theta)$, can then be expressed as (Vernin & Roddier 1975),

$$\begin{eqnarray} C_{\rm I}(r,\theta) & = & \left( \frac{B_{1}^{2} + B_{2}^{2}}{(... ...B_{2})^{2}} \right)[ A_{\rm I}(r,\theta) + A_{\rm I}(r,-\theta) ]\end{eqnarray}$$ (10)

where B₁ and B₂ are the relative brightnesses of the binary star components, $A_{\rm I}(r,0)$ represents the central (single star) correlation peak and where $A_{\rm I}(r,\pm\theta)$ contains the information that is required regarding the strength and position of the turbulent layers.

The log amplitude covariance of the scintillation fluctuations then follows as,

$$A_{\chi}(r,\theta) = \frac{ \ln [A_{\rm I}(r,\theta) + 1]}{4},$$ (11)

and finally, if one can assume both the Rytov approximation (i.e. that $A_{\chi}(0) = \sigma^{2}_{\chi} < 0.3$, where $\sigma^{2}_{\chi}$ is the variance of the intensity fluctuations within a single pixel) and Kolmogorov statistics for the optical turbulence fluctuations, then it can be shown that this function can be expressed in terms of C_n²(h) as (Tyler 1992),

$$A_{\chi}(r,\theta) = \frac{8.16 k^2}{4\pi} \int_{0}^{\infty}{ C_{n}^{2}(h) h^{5/6} F(Q) {\rm d}h } ,$$ (12)

where,

$$F(Q) = \int_{0}^{\infty}{ K^{-8/3} J_{0}(QK) [1 - \cos(K^2)] {\rm d}K },$$ (13)

and

$$Q = \left[ \frac{k}{h {\rm sec}(z)} \right]^{1/2} \mid r - \theta h {\rm sec}(z) \mid. \vspace{0.5cm}$$ (14)

Therefore, given the known quantities (13) and (14) one is required to calculate C_n²(h) from Eq. (12) - a classic inversion problem (specifically, a Fredholm integral of the first kind).

The particular method one might employ for this problem (and also to some extent the confidence of the solution) depends very much upon the C_n²(h) altitude sampling that is required. Tests performed upon a number of high resolution C_n²(h) profiles obtained by balloon-borne instruments (Vernin & Muñoz-Tuñôn 1992) have shown that the main parameters which determine AO anisoplanatic performance can be calculated with reasonable accuracy provided that the altitude sampling, $\Delta h$, is less than 1km. For the data collected thus far, the altitude resolution available is between around 200-500 m, and therefore comfortably within this constraint. However it has been found that a simple least-squares inversion of Eq. (12) for these data is incapable of producing a reliable solution unless one rectangularises the two-dimensional integral equation such that multiple pupil samples map onto each altitude sample and $\Delta r$ remains fixed at 1 cm pix^-1 but $\Delta h$ is extended to 2-3 km.

Although, ultimately, it is the Fresnel length that limits the altitude resolution attainable, in an effort to realise the altitude sampling that is made available by that within the pupil, a number of iterative inversion methods have been explored. One that has proved successful is based upon a maximum entropy routine which also allows for the controlled propagation of noise between the data and solution (Skilling & Bryan 1984). In practice, we first calculate a crude (large $\Delta h$) least-squares inversion of C_n²(h) using a rectangularised version of Eq. (12) and then use an interpolated version of this as the first guess estimate for the iterative routine. As an aid to deciding when the routine should be terminated, the following RMS difference error metric is calculated following each iteration:

$${\cal E}_{\rm A} = \sqrt{ \frac{< (A_{\rm data}(r,\theta)-A_... ... best-fit}(r,\theta))^2\gt} {< A^2_{\rm data}(r,\theta) \gt} },$$ (15)

where, $A_{\rm data}(r,\theta)$ is the observed one-dimensional slice through the correlation plane and $A_{\rm best-fit}(r,\theta)$ is the function calculated using the current (best-fit) solution to C_n²(h) from Eq. (12).

In Fig. 5, (a) shows the correctly scaled correlation slice, $A_{\rm data}(r,\theta)$, of Fig. 4a, and in (b) the corresponding C_n²(h) profile obtained. Also shown in Fig. 5a, for comparison is $A_{\rm best-fit}(r,\theta)$ which gives some indication of the accuracy of the solution obtained where, in this case, ${\cal E}_{\rm A} = 4.3\%$. It is interesting to note that for the particular example shown here the apparent single correlation peak at a height of 7.5 km in Fig. 5a has been resolved into two distinct turbulent layers within the corresponding C_n²(h) profile of Fig. 5b at heights of 7.2 km and 8.1 km.

The apparent width of the dominant turbulent layers shown here probably does not reflect the true physical depth over which the optical turbulence exists. It is more likely that this apparent spread of C_n²(h) values is due to the effective resolution imposed by the Fresnel length, as discussed in Sect. 2.3.1. Therefore, although it appears that one is oversampling each peak in the C_n²(h) profile, the suspicion is that in fact the layers of turbulence themselves could be substantially thinner.

Confidence in the inversion methods used has been obtained by testing them on realistic simulated SCIDAR data generated using a Fresnel wavefront propagation code developed within our group (Adcock et al. 1996).

  

Figure 6: a) Simulated spatial correlation data along with b) a comparison of the model input Cn2(h) profile (bar-plot) and that obtained using our integral inversion method

Figure 6 gives one such result, which shows a typical simulated correlation image used in these tests and compares the true C_n²(h) input data with the solution obtained. Although the blurring effect of the Fresnel length is clearly apparent, the integrated values of C_n²(h) for each layer agree to within a few percent.

2.3.2 Determination of V(h)

Prior to calculation of the velocity profile, V(h), the same initial processing steps are taken as for the C_n²(h) except that now the average cross -correlation of consecutive mean-normalised images is calculated, each separated by a time interval $\Delta t$. The resulting two-dimensional function thus represents the average spatio-temporal correlation of the turbulence, the sampling interval of which is given simply by $\Delta V = \Delta r/\Delta t$.

Figure 4b shows such a plot generated from the same data as in Fig. 4a, with $\Delta t$ = 23 ms (corresponding to $\Delta V$=0.42 ms^-1). The zero velocity position within this plot is denoted at the centre by a cross. In comparing this with the spatial-correlation plot of Fig. 4a one can see that, as indicated by the the C_n²(h) of Fig. 5b, the apparent single peak at the greatest distance from the central core has been resolved into two layers at heights of 7.2 km and 8.1 km with associated velocities of 14.6 ms^-1 (wind direction Position Angle, PA = 103$^{\circ}$) and 9.8 ms^-1 (PA = 130$^{\circ}$) respectively. It is also worth noting that the correlation peak closest to the central core (corresponding to a height close to the telescope pupil) has shifted only very slightly from its central position in the spatial-correlation (V = 0.01 ms^{- 1}). This is indicative perhaps of the presence of almost static turbulence within the telescope dome since the ground wind speed recorded at this time was 13 ms^-1.

The algorithm to allow automative calculation of V(h) is based upon the assumption that the optical turbulence is highly discretised (the justification for this is clarified later in Sect. 4). Having first determined the altitudes of each discrete layer and the orientation of the binary star observed, a template is then constructed for each layer, each one representing their individual spatial correlation signal. The cross-correlation of each of these templates with the spatio-temporal plot then allows the speed and direction of each layer to be identified.