text:
# Photon and Pion Production at High 𝑝_𝑇

## Abstract

We present a study of high-$`p_T`$ photon and pion production in hadronic interactions, focusing on a comparison of the yields with expectations from next-to-leading order perturbative QCD (NLO pQCD). We examine the impact of phenomenological models of $`k_T`$ smearing (which approximate effects of additional soft-gluon emission) on absolute predictions for photon and pion production and their ratio.

Single and double direct-photon production in hadronic collisions at high transverse momenta ($`p_T`$) have long been viewed as an ideal testing ground for the formalism of pQCD. A reliable theoretical description of the direct-photon process is of special importance because of its sensitivity to the gluon distribution in a proton through the quark–gluon scattering subprocess ($`gq\gamma q`$). The gluon distribution, $`G(x)`$, is relatively well constrained for $`x<0.1`$, but much less so at larger $`x`$ huston-uncertainty . In principle, fixed-target direct-photon production can constrain $`G(x)`$ at large $`x`$, and such data have therefore been incorporated in several modern global parton distribution function (PDF) analyses cteq4 ; grv92 ; mrst .

However, both the completeness of the NLO description of the direct-photon process, as well as the consistency of results from different experiments, have been questioned mrst ; baerreno ; huston-discrepancy ; E706-kt ; ktprd ; aurenche-dp ; aurenche-pi0 ; kimber . The inclusive production of hadrons provides a further means of testing the predictions of the NLO pQCD formalism. Deviations have been observed between measured inclusive direct-photon and pion cross sections and NLO pQCD calculations. Examples of such discrepancies are shown in Fig. 1 where ratios of data to theory are displayed as a function of $`x_T=2p_T/\sqrt{s}`$ for photon and pion data. (Unless otherwise indicated, all NLO calculations aurenche-nlo ; qiu ; gordon ; aversa ; bailey  in this paper use a single scale of $`\mu =p_T/2`$, CTEQ4M PDFs cteq4 , and BKK fragmentation functions for pions BKK .) It has been suggested that part of the deviations from theory for both photons and pions can be ascribed to higher-order effects of initial-state soft-gluon radiation huston-discrepancy ; E706-kt ; ktprd .

Given the scatter of the data shown in Figs. 1, it may be instructive to consider measurements of the $`\gamma /\pi ^0`$ ratio over a wide range of $`\sqrt{s}`$ osaka . Both experimental and theoretical uncertainties tend to cancel in such a ratio, and the ratio should also be less sensitive to incomplete treatment of gluon radiation.

A compilation of comparisons between data and theory, shown for simplicity without their uncertainties, is presented in Fig. 2. In this figure, the ratios of data to theory for the $`\gamma `$ to $`\pi ^0`$ measurements have been approximated as a constant value at high-$`p_T`$, and the results plotted as a function of $`\sqrt{s}`$ (see osaka  for details). The figure suggests an energy dependence in the ratio of data to theory for $`\gamma /\pi ^0`$ production. There is, however, an indication of substantial differences between the experiments at low $`\sqrt{s}`$ (where the observed $`\gamma /\pi ^0`$ is smallest), which makes it difficult to quantify this trend. Recognizing the presence of these differences is especially important because thus far only the low energy photon experiments have been used in PDF fits to extract the gluon distribution.

The differences between many of the data sets and pQCD, seen in Fig. 1, may be due to the impact of the effective parton transverse momentum, $`k_T`$. In hadronic hard-scattering processes, there is generally a substantial amount of $`k_T`$ in the initial state resulting from gluon emission ktprd . The presence of $`k_T`$ impacts the final state and has been observed in measurements of Drell-Yan, diphoton, and heavy quark production; the amount of $`k_T`$ expected from NLO calculations is not sufficient to describe the data. The effective values of $`k_T`$/parton for these processes vary from $`1`$ GeV/$`c`$ at fixed target energies, as illustrated in Fig. 2 for diphoton distributions from E706 begel , to 3–4 GeV/$`c`$ at the Tevatron Collider — the growth is approximately logarithmic with center-of-mass energy ktprd . The size of the $`k_T`$ values, and their dependence on energy, argue against a purely “intrinsic” non-perturbative origin. Rather, the major part of this effect is generally attributed to soft-gluon emission. While the importance of including gluon emission through the resummation formalism has long been recognized and calculations have been available for some time for Drell-Yan altarelli , diphoton RESBOS ; fergani , and W/Z production RESBOS , they have only recently been developed for inclusive direct-photon production nason ; kidonakisowens ; laenen ; lilai ; li ; sterman .

In the absence of a rigorous theoretical treatment of the impact of gluon emission on high-$`p_T`$ inclusive production, a more intuitive phenomenological approach has proved successful ktprd . The soft gluon radiation was parametrized in terms of an effective $`k_T`$ that provided an additional transverse impulse to the outgoing partons. Because of the steeply falling cross section in $`p_T`$, such a $`k_T`$ can shift the production of final-state particles from lower to higher values of $`p_T`$, effectively enhancing the cross section.

As described in ktprd , a leading-order (LO) pQCD calculation owens  has been used to generate K-factors (ratios of calculations for any given $`k_T`$ to the result for $`k_T=0`$) for inclusive cross sections. These $`p_T`$-dependent factors have been then applied to the NLO pQCD calculations. The enhancements that would be expected for direct-photon production from parton-showering models pythia61 ; herwig61  have also been investigated osaka . These programs do not provide sufficient smearing at fixed-target energies because shower development is constrained by cut-off parameters that ensure the perturbative nature of the process. Consequently, these calculations allow additional input $`k_T`$ for Gaussian smearing, and are often used that way in comparisons to data. The respective corrections have been obtained using default settings for other program parameters and an input $`k_T`$ of 1.2 GeV/$`c`$ for the smearing, relative to these same settings with $`k_T=0`$, and then applied to NLO pQCD calculations. The resulting comparisons to data from E706 E706-kt  are displayed in Fig. 3 (similar results hold for pion production, not shown). The observed differences should be kept in mind when comparing these models to data for $`k_T`$-sensitive quantities.

Recently, there has been significant progress in more rigorous resummed pQCD calculations for single direct-photon production nason ; kidonakisowens ; laenen ; lilai ; li ; sterman . Substantial corrections to fixed-order QCD calculations are expected from soft-gluon emission, especially in regions of phase space where gluon emission is restricted kinematically. At large $`x`$, there is a suppression of gluon radiation due to the rapidly falling parton distributions and a complete description of the cross section in this region requires the resummation of “threshold” terms. Two recent threshold-resummed pQCD calculations for direct photons nason ; kidonakisowens  exhibit far less dependence on QCD scales than found in NLO theory. These calculations agree with the NLO prediction for the scale $`\mu p_T/2`$ at low $`p_T`$ (without inclusion of explicit $`k_T`$ or recoil effects), and show an enhancement in cross section at high $`p_T`$.

A method for simultaneous treatment of recoil and threshold corrections in inclusive single-photon cross sections is being developed sterman  within the formalism of collinear factorization. This approach accounts explicitly for the recoil from soft radiation in the hard-scattering subprocess, and conserves both energy and transverse momentum for the resummed radiation. The possibility of substantial enhancements from higher-order perturbative and power-law nonperturbative corrections relative to NLO are indicated at both moderate and high $`p_T`$ for fixed-target energies, similar to the enhancements obtained with the simple $`k_T`$-smearing model discussed above. Figure 3 (right) displays the results of an example calculation sterman  based on this approach compared with direct-photon measurements from E706.

While there is still no resummation calculation for inclusive pion production, the trend of recent developments in direct-photon processes has led to an increased appreciation of the importance of the effects of multiple gluon emission, and to the emergence of tools for incorporating these effects. These latest theoretical developments encourage optimism that the long-standing difficulties in developing an adequate description of these processes can eventually be resolved, making possible a global re-examination of parton distributions with an emphasis on the determination of the gluon distribution from the direct-photon data sterman-pdf .

###### Acknowledgements.

This work has been done in collaboration with L. Apanasevich, M. Begel, C. Bromberg, T. Ferbel, G. Ginther, J. Huston, S. Kuhlmann, P. Slattery, and V. Zutshi. We also wish to thank P. Aurenche, C. Balázs, S. Catani, M. Fontannaz, N. Kidonakis, J. F. Owens, E. Pilon, G. Sterman, W. Vogelsang and J. Womersley for helpful discussions.


format: text

text:
# An HST Lensing Survey of X-ray Luminous Galaxy Clusters: I. A 3831footnote 11footnote 1Based on observations with the NASA/ESA Hubble Space Telescope obtained at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy Inc., under NASA contract NAS 5-26555.

## 1. Introduction

Massive clusters of galaxies represent one extreme of the mass spectrum of collapsed structures at the present day. Their properties are expected to reflect predominantly gravitational processes and therefore provide unique insights into the nature and distribution of the dark matter which drives the formation of structure (e.g. Eke et al. 1996; Viana & Liddle 1996; Bahcall, Fan & Cen 1997; Kay & Bower 1999). In particular, accurate measurements of the mass distribution in clusters across a range of scales can be used to test the claim of a universal form for the profiles of dark matter halos (Navarro, Frenk & White 1997, NFW) and hence the nature of dark matter. At large scales the NFW profile falls off as $`\rho r^3`$, steeper than an isothermal model ($`\rho r^2`$). On smaller scales the NFW profile breaks to a shallower slope with $`\rho r^1`$, whilst retaining a central cusp. Analyses of the mass profiles of dwarf galaxies from dynamical studies have used the apparent lack of a central cusp in the density distribution to reject the form of the NFW profile and argue instead for the existence of self-interacting dark matter (e.g. Moore et al. 1998). However, on larger scales, (e.g. luminous elliptical galaxies and massive clusters of galaxies) the form of the mass profile has yet to be investigated systematically.

The central regions of massive, compact clusters at moderate redshifts ($`z<0.5`$) act as strong gravitational lenses forming multiple images of serendipitously placed background galaxies (e.g. Smail et al. 1996). The properties of the images of these background galaxies (position, relative surface brightness and parity) can be used to accurately model the distribution of total mass (both baryonic and non-baryonic) within the cluster core (Kneib et al. 1996, K96; Smail et al. 1996; Natarajan et al. 1998). Such mass maps provide the most direct and detailed view of the distribution and morphology of dark matter in the centers of galaxy clusters. Strong gravitational lensing by rich clusters is therefore an ideal observational tool with which to test the halo properties predicted by cosmological models.

The first mass maps constructed from lensing observations were produced in the early and mid-1990’s (Tyson et al. 1990; Mellier et al. 1993; Fahlman et al. 1994; Kneib et al. 1994; Smail et al. 1995). A significant improvement arrived with the refurbishment of HST which provided the resolution necessary to identify faint, lensed features in the crowded cores of rich clusters (e.g. K96), thereby allowing more detailed mass models to be constructed.

On small scales, unique constraints can be obtained from radial arcs (rare radially magnified images of background galaxies) which are found in the very centers of a few cluster lenses. Such features are very difficult to detect in ground-based observations, with only a single example known (MS 2137.3$``$23; Fort et al. 1992; Mellier et al. 1993). In contrast to this, HST has uncovered several radial arcs in previously well-studied cluster lenses (e.g. A 370, Smail et al. 1996; AC 114, Natarajan et al. 1998) which enabled the first tests of the form of cluster mass profiles on scales of less than $`100`$ kpc (Williams et al. 1999).

Previous HST studies of rich clusters covered a heterogeneous mix of clusters, either selected because they were previously well-studied (e.g. Smail et al. 1997) or because they were known to be strong lenses (e.g. K96; Smail et al. 1996). To study the form of the mass profile in rich clusters in an unbiased fashion we need HSTobservations of an objectively selected cluster sample. Ideally such a sample would be mass-selected; however, in the absence of samples compiled from large-scale weak lensing surveys, X-ray selected cluster samples (Gioia et al. 1990; Ebeling et al. 1998, 2000; DeGrandi et al. 1999; Ebeling, Edge & Henry 2000) are best suited to selecting well-defined samples of massive clusters.

Following this premise, we are conducting a survey of twelve of the most X-ray luminous clusters ($`L_\mathrm{X}8\times 10^{44}`$ erg s-1, 0.1–2.4 keV) in a narrow redshift slice at $`z0.2`$, selected from the XBACs sample (X-ray Brightest Abell-type Clusters; Ebeling et al. 1996). As XBACS is restricted to Abell clusters (Abell, Corwin & Olowin 1989), it is X-ray flux limited and not truly X-ray selected. However, a comparison with the X-ray selected ROSAT Brightest Cluster Sample (BCS, Ebeling et al. 1998, 2000) shows that $`75`$% of the BCS clusters in the redshift and X-ray luminosity range of our sample are in fact Abell clusters. Hence, our XBACs sample is, in all practical aspects, indistinguishable from an X-ray selected sample.

In this paper we describe HST observations of the core of one of the first clusters observed in our survey, A 383 ($`z=0.188`$), which reveal a multitude of strongly-lensed features. In combination with color information obtained from ground-based imaging observations in three passbands, the lensed features are used to produce an accurate model of the mass distribution within the central 500 kpc of the cluster. We complement these results with an analysis of the X-ray properties of A 383 using archival ROSAT HRI data. In §2 we describe the observational data, their reduction and analysis, followed in §3 by a discussion of the interpretation of these data and the construction of the lens model of the cluster. Finally in §4 we summarize the main results of our analysis and present our conclusions. We adopt $`H_0=50`$ kms-1 Mpc-1 and $`q_0=0.5`$ throughout. At the cluster redshift $`1^{\prime \prime }4.0`$ kpc in this cosmology.

## 2. Data Reduction and Analysis

A 383 is a massive cluster of galaxies at $`z=0.188`$ with a core dominated by a luminous cD galaxy (Figs. 1 & 3). The X-ray luminosity from a pointed observation with the ROSAT HRI is $`(9.8\pm 0.3)\times 10^{44}`$ erg s-1 in the 0.1–2.4 keV band (see §2.4).

### 2.1. HST Imaging

A 383 was observed with the HST WFPC2 camera on January 25, 2000. Three exposures totaling 7.5 ks were taken through the F702W filter. Each exposure was shifted relative to the others by 10 WFC pixels ($`1.0^{\prime \prime }`$) providing a partial overlap of the chip fields. After pipeline processing, standard iraf/stsdas routines were employed to shift and combine the frames to remove both cosmic rays and hot pixels. Corrections for under-sampling of the point spread function and geometric distortion of the optics were made using the dither package within iraf (Fruchter & Hook 1997). The final frame (Fig. 1) has a pixel scale of 0.05′′, an effective resolution of 0.15′′ and a 1$`\sigma `$ detection limit within the seeing disk of $`R_{702}31`$.

To produce a catalog of faint arclets from our data we first analysed the HST frame using the SExtractor package (Bertin & Arnouts 1996). All objects with isophotal areas in excess of 30 pixels ($`0.074`$ arcsec2) at the $`\mu _{702}=24.5`$ mag arcsec-2 isophote ($`3\sigma `$/pixel) were selected. A steep roll-over in the observed differential number counts arising from incompleteness occurs at $`R_{702}25.5`$. We estimate the completeness to be $`80`$% (approximately 5$`\sigma `$) at $`R_{702}25.5`$, which we adopt as our magnitude limit, giving a total of 457 sources within a 4.9 arcmin2 area (excluding the PC chip). Seventeen of these sources are classified as star-like on the basis of their profile shapes and are excluded from our analysis.

For the purposes of the photometric analysis presented in Table 1 we correct the F702W photometry ($`R_{702}`$) to the Cousins $`R`$-band in the following manner. The observed range in $`(BR)`$ colors of the arcs roughly translates to $`0.2<(VR)<1.0`$ giving a typical color of $`(VR)0.6`$. Adopting this color we obtain a correction of $`R_{702}R=0.2`$ (Holtzman et al. 1995), and estimate that this correction introduces a systematic uncertainty in the $`R`$-band photometry of $`\pm 0.06`$ magnitudes.

The final HST frame (Fig. 1) reveals many previously unknown strongly lensed features, including a giant arc and two radial arcs, making A 383 a striking new addition to the catalog of cluster lenses at intermediate redshifts.

### 2.2. Ground-based Imaging

With the aim of extending our lensing mass map to the turn-around radius of the cluster, we have used the 3.6m Canada-France-Hawaii Telescope (CFHT) with the CFH12K camera to obtain panoramic images of A 383 on the nights of November 14–16, 1999. Total exposure times of 7.2 ks, 6.0 ks and 3.6 ks, accumulated at 6–10 dither positions, were acquired in the $`B`$, $`R`$ and $`I`$ bands respectively. Data reduction was performed within iraf using the mscred package including standard bias subtraction and flat-fielding using twilight flats. The dithered exposures were aligned with the Digital Sky Survey frame of the same field to an rms accuracy of 0.15′′. More information on the reduction and analysis of these observations will be presented in a forthcoming paper (Czoske et al. 2000). Here we use the central regions of the $`B`$ and $`I`$ band frames to provide photometry of the lensed features in the cluster core (see Table 1). These two frames have seeing of 0.88′′ and 0.71′′ FWHM respectively.

In addition to the CFH12K imaging, we have obtained $`K`$-band images of the core of A 383 with the UFTI imager on the 3.8m United Kingdom Infrared Telescope (UKIRT), Mauna Kea, on October 14, 1999. The final frame was accumulated in 27 dithered sub-exposures of 90 s duration each to give a total on-source integration time of 2.4 ks, all in photometric conditions. Employing standard procedures these frames were reduced, combined and calibrated using observations of UKIRT faint standards bracketing the science exposures. The final frame has seeing of 0.42′′ FWHM with 0.09′′/pixel sampling and an effective 5$`\sigma `$ depth of $`K=20.3`$.

Optical and optical-infrared colors of the lensed features in A 383 were measured off the aligned and seeing-matched $`BRIK`$ frames using apertures designed for individual arcs. The arcs were first masked out of the science frame and the sky background estimated by median smoothing over the aperture defined by the absent arc. The sky subtraction was then performed and the individual arc apertures applied to the sky-subtracted frame to obtain the $`BRIK`$ photometry presented in Table 1. Error bars were estimated by varying the smoothing length used in estimating the sky background.

### 2.3. Spectroscopy

On January 26, 2000 we observed the cD galaxy in A 383 with the Wide Field Grism Spectrograph on the University of Hawaii’s 2.2-m telescope. The low dispersion spectrum obtained in a 1.8 ks exposure in extremely poor seeing (2.7′′) shows the full range of emission lines typically found in cD galaxies in cooling flow clusters (e.g. Crawford et al. 1999) and places the cD at a heliocentric redshift of $`z=0.1880\pm 0.0012`$.

On January 29, 2000 we obtained 28 additional spectra with the LRIS spectrograph (Oke et al. 1995) on the Keck-II 10-m telescope in MOS mode, with a total exposure time of 3.6 ks in average seeing of 0.9′′. Use of the 300/5000 grating centered at 7500 Å provided wide spectral coverage (5000–10,000 Å) at a spectral resolution of 2.55 Å/pixel. Redshifts were measured independently by three members of our team (OC, HE, JPK) and also using the rvsao cross-correlation package under iraf to accurately estimate the errors. Successive $`3\sigma `$ clipping around the peak of the radial velocity distribution yields a redshift of $`z=0.1880\pm 0.0002`$ for A 383 and a velocity dispersion of $`\sigma =1150_{200}^{+150}`$ kms-1 from 18 galaxy redshifts. Ten other galaxies were found to be background, including B16 at $`z=0.656`$, which lies close to the cluster core (see Fig. 1).

Fig. 2. – The spectrum of arc B0a in A 383 taken with the LRIS spectrograph on the Keck-II 10m telescope. The arc exhibits a blue continuum and we identify the strong emission line as \[Oii\]$`\lambda `$3727 and confirm this with several UV absorption features (marked) to give a redshift of $`z=1.0103\pm 0.0001`$. The upper spectrum is smoothed to the nominal resolution of the spectrograph and offset vertically for clarity. The shaded regions show areas of the spectrum strongly affected by night sky lines.

We also obtained a spectrum of the lensed feature B0a which we show in Fig. 2. We identify the strong emission line at 7493Å as \[Oii\] $`\lambda `$3727 which places the galaxy at a redshift of $`z=1.0103\pm 0.0001`$. This interpretation is confirmed by the identification of Feii$`\lambda `$2600 and Mgii$`\lambda `$2800 absorption features. The properties of this galaxy are discussed further in §3.

### 2.4. ROSAT/HRI Observation

A 383 was observed with the ROSAT HRI in February 1995 for a total integration time of 27.4 ks (ROR #800890). An adaptively smoothed map of the observed X-ray emission (created using the asmooth algorithm of Ebeling, White & Rangarajan 2000) is shown in Fig. 3. We used a modified version of Steve Snowden’s cast hri software to correct for exposure time variations across the field of view and subtract a particle background component.

From the fully processed ROSAT HRI image the total background-corrected HRI count rate from A 383 is measured to be $`(0.128\pm 0.005)`$ ct s-1 within 1.5 Mpc (6.2) of the cluster center, using the mean diffuse background at $`r>2`$ Mpc (8.2). Assuming a standard plasma emission spectrum, a Galactic hydrogen column density of $`4.1\times 10^{20}`$ cm-2 for this sight-line (Dickey & Lockman 1990), a metallicity of 0.3 and an ambient X-ray temperature of $`7.1`$ keV (estimated from the cooling-flow-corrected cluster $`L_\mathrm{X}`$$`kT`$ relation of Allen & Fabian 1998) we derive a total unabsorbed cluster flux of $`(6.55\pm 0.23)\times 10^{12}`$ erg cm-2 s-1 (0.1–2.4 keV). The corresponding X-ray luminosity is $`(9.8\pm 0.3)\times 10^{44}`$ erg s-1 (0.1–2.4 keV), in agreement with the XBACs value of $`(8.0\pm 2.4)\times 10^{44}`$ erg s-1 in the same energy band (Ebeling et al. 1996).

fig3.gif

Fig. 3. – X-ray flux contours of the adaptively smoothed emission from A 383 as seen with the ROSAT HRI, overlaid on an $`R`$-band image of the cluster taken with the CFH12K camera (see §2.2). All X-ray features are at least $`3\sigma `$ significant with respect to the local background. The lowest contour lies 30% above the background, adjacent contours differ by 30%. We adjusted the astrometry of the HRI image by 4′′ (well within the errors of the astrometry solution) to align the X-ray peak with the cD galaxy. The near circular symmetry of the X-ray flux contours reflects the symmetry of the lensing mass contours plotted in Fig. 1. The region shown in our Fig. 1 is marked by the dashed box.

## 3. Modelling and Results

In this section we first describe and interpret the multiply-imaged features observed in our HST frame of A 383 and explain how these were used to constrain the mass distribution in the central regions of the cluster. We then discuss the form of the mass profile we derive in the very center of the cluster and compare this to theoretical expectations from high resolution N-body simulations. Finally we discuss the X-ray properties of A 383 and compare the mass distribution determined from the lens model with the estimates based on the ROSAT data.

### 3.1. Multiply Imaged Features

As a result of the superb resolution of our HST images, a detailed analysis of the morphology of the observed arcs allows us to identify seven multiply-imaged systems in the field of A 383. This is significantly greater than the number typically seen in ground-based observations of clusters. We discuss each multiply-imaged system below, and summarize their photometric and spectroscopic properties in Table 1. The lensed features are identified in Fig. 1 and the multiple images are further illustrated below in Fig. 4–6.

#### 3.1.1 The Giant and Radial Arcs

We interpret the giant and radial arcs in A 383 as being the lensed images of two background galaxies, B0 and B1, at $`z=1.01`$ (§2.3) and $`z1.1`$ respectively. This interpretation arises from modelling (§3.2) of the radial arcs, which constrains the outer radial arc (B1d) to be a counter image of B1a/b/c at $`z1.1`$, and the inner radial arc (B0b) to be a counter image of B0a. B0a is also slightly closer to the center of the cluster than other components of the giant arc, supporting the suggestion that it is at a slightly lower redshift than B1.

Fig. 4. – The giant arc in A 383 (see text for interpretation). Each tick mark represents $`1^{\prime \prime }`$, and the orientation is as in Fig. 1. The dotted lines indicate the path of the $`z=1.1`$ critical line close to the observed arcs.

B0 (Figs. 4 & 5) — This system comprises the brightest component of the giant arc (B0a) and the inner of the two radial arcs (B0b). It is also plausible that B15 is a component of the B0 system given its proximity to B0a in the image plane, however it’s redder color (Table 1) suggests that it is either singly imaged, or that B0b has a significant color gradient. We are unable to measure the color of B0b to confirm this possibility due to the light from the central galaxy. We have also estimated the unlensed $`R`$-band magnitude of B0 (Table 1) which, along with it’s blue optical-infrared color, suggests that B0 is a low luminosity star-forming galaxy ($`M_V19.3`$).

Fig. 5. – A view of the radial arcs as seen on the HST frame (enhanced by subtracting a median-smoothed frame from the science frame). Each tick mark represents $`1^{\prime \prime }`$, and the orientation is as in Fig. 1.

B1 (Figs. 4 & 5) — B1a/b is a pair of merging images straddling the $`z1.1`$ critical line with counter-images B1c, at the extreme right hand end of the giant arc and B1d, the outer of the two radial arcs. The lens model (§3.2) and photometric analysis (Table 1) support this interpretation and appear to exclude B15 from being a further counter-image. The unlensed $`R`$-band magnitude of B1 has been estimated (Table 1), from which we determine that it is a very low luminosity galaxy ($`M_V18.3`$).

#### 3.1.2 Perturbations by Galaxy Halos

B2 and B3 (Fig. 6) — The lensing effect of neighboring cluster galaxy halos on this complex group of blue arcs demonstrates how galaxy-scale masses measurably perturb the cluster potential. We have modeled the effect of these perturbations and resolved B2 and B3 into two background galaxies, both at a redshift of $`z3.5`$. Although the photometric analysis is hampered by the nearby galaxy halos, the data presented in Table 1 lend additional support to this interpretation.

B2 comprises five images, of which B2a/b and B2d/e are pairs of images straddling the $`z3.5`$ critical line as it passes around the halos of cluster ellipticals #46 and #23 respectively. B2c lies outside this critical line and has the same parity as B2a. A sixth candidate counter image lies approximately 1′′ to the left of B3b under the halo of cluster elliptical #23. Multi-band HST imaging is necessary to achieve the accurate photometry required to confirm this hypothesis. There is also a low signal-to-noise feature lying in the saddle region between cluster ellipticals #46 and #76, which is probably an extension of B2c. Counter images of this portion of the B2 system are undetected, due to the intrinsic faintness of the galaxy and the relatively low magnification of the counter images relative to the influence of the saddle on the observed feature. Finally, we estimate that the B2 background galaxy has $`M_V22.0`$ and the blue colors typical of a star-forming galaxy.

Fig. 6. – The B2/B3 system of arcs (see text for more details). Cluster ellipticals #46, #76 and #23 are also labelled. Each tick mark represents $`1^{\prime \prime }`$, and the orientation is as in Fig. 1. The dotted lines indicate the path of the $`z=3.5`$ critical line close to the observed arcs.

B3 consists of three images; B3b/c is adjacent to cluster elliptical #23, and B3a, the counter image, lies outside the critical line, $`1^{\prime \prime }`$ below and to the right of B2c in Fig. 6. B3d appears to be a singly-imaged component of the background galaxy that produces B3a/b/c. This interpretation is supported by the photometric analysis in Table 1. We also estimate that B3 is a star forming galaxy with a similar luminosity ($`M_V22.4`$) to B2.

B4 (Fig. 1) — This is another complex and unusual arc, the precise nature of which is uncertain. One possibility is that B4 is part of a multiple image system containing B6; however, photometric analysis (Table 1) does not support this hypothesis. More complete color information will improve our ability to comment on this arc. If B4 is a singly imaged source it lies at a redshift of $`z<1.0`$.

#### 3.1.3 Faint Arcs

B5 (Fig. 1) — B5 is a very faint arc lying in the saddle region between the central galaxy and cluster elliptical #99. We interpret this arc as an image of a background galaxy which is magnified by this saddle and therefore constrain the redshift of the galaxy to be $`z>2.4`$.

B9a & B9b (Fig. 1) — This is a pair of extremely faint arcs detected only because of the superb resolution and low sky background of the HST observations. Both arcs appear to comprise three images resulting from the magnifying effect of the saddle potential between the cluster center and cluster elliptical #260. We estimate that they lie at a redshift of $`z>3`$.

B14 (Fig. 1) — B14 is barely detected above the sky background in the HST frame, however we detect it at $`K=19.70\pm 0.05`$ in the UKIRT frame, giving $`(RK)=6.02\pm 0.11`$. As B14 is singly-imaged, we obtain a redshift limit of $`z<3.9`$ from our lens model. A full analysis of galaxies exhibiting extreme colors and the general background population lying behind our cluster sample will be presented in a future paper.

#### 3.1.4 Arclets

We also detect numerous singly-imaged arclets in the field of A 383 (Fig. 1 and Table 1) and have obtained upper limits on the redshifts of three of these arclets (B6, B7 and B8) on the basis of their singly-imaged nature.

B6 (Fig. 1) — We obtain a redshift limit of $`z<0.9`$ for this arclet, however it’s blue $`(BR)`$ and $`(RI)`$ colors (Table 1) are consistent with it being a star forming galaxy at $`z1.5`$–2.0. The latter interpretation supports the possibility that B6 and B4 contain multiply imaged features of a single background galaxy. As noted in §3.1.2, improved photometry of B4 will help to resolve this uncertainty.

B7 (Fig. 1) — B7 lies at a redshift of $`z<1.3`$ which, on the basis of its $`(BR)`$ color (Table 1) is consistent with B7 being a moderately star forming galaxy.

B8 (Fig. 1) — This arclet has a redshift limit of $`z<2.0`$, while its colors are broadly consistent with this being a star forming galaxy at $`z1.0`$–1.4.

### 3.2. Lens Model

The mass distribution of A 383 was reconstructed using the parametric lens inversion method described in detail by K96. In addition to a cluster-size mass component, the model contains 29 individual galaxy masses (including the central galaxy) covering the entire HST field of view. The parameters describing the model were constrained were constrained within the cluster core using four of the multiple image systems identified in the HST frame (B0, B1, B2 and B3) and on larger scales using the weak shear field estimated from the shapes of the faint galaxies across the whole frame. The cluster scale mass that dominates the mass distribution has a core radius of ($`46\pm 3`$) kpc, velocity dispersion of ($`900\pm 20`$) kms-1 and a cut-off radius of ($`1200\pm 100`$) kpc (see K96 for details of how these quantities are defined).

The mass distribution obtained from the lens model is overlaid on the optical image in Fig. 1. We calculate the total projected mass within the radius of the giant arc (65 kpc) to be $`(3.5\pm 0.1)\times 10^{13}M_{}`$. Using the weak shear measurements we can extend our model to larger scales and constrain the projected mass within a 250 kpc radius of the cluster center to be: $`(1.8\pm 0.2)\times 10^{14}M_{}`$.

In Fig. 7 we plot the projected surface mass density derived from the lens model as a function of the distance from the cluster center ($`r`$). On scales $`r>10`$ kpc, the density profile is dominated by the cluster-scale potential, the core radius of the cluster-scale mass distribution ($`46`$ kpc) being consistent with previous estimates of core radii from lensing studies (Mellier et al. 1993; Smail et al. 1996; Natarajan et al. 1998). The shape of A 383’s profile is however only constrained out to $`r400`$ kpc, this being the extent of the WFPC2 field of view. The mass profile is most accurately measured within $`r100`$ kpc where information from the multiple images is available.

Fig. 7. – The projected surface mass density ($`\mathrm{\Sigma }`$) of A 383 as a function of radius, obtained from the lens model by measuring the mass density in concentric annuli. The dashed lines illustrate the difference between the profile slope at the positions of the two radial arcs (B0b and B1d), which are marked by arrows at 6.0 and 19.9 kpc respectively.

#### 3.2.1 Radial Arc Properties

B0b and B1d (see Figs. 1 & 5) lie at distances of ($`6.0\pm 0.4`$) kpc and ($`19.9\pm 0.4`$) kpc from the cluster center respectively. These two arcs therefore enable us to probe the shape of the mass distribution on scales of $`r<20`$ kpc. Specifically, the radial arcs constrain the slope of the density profile, as their radial positions depend on the local gradient of the projected mass, rather than on the mean enclosed surface density.

A detailed analysis of the shape of A 383’s density profile is currently limited by the pseudo-isothermal elliptical functional form adopted in the lens model (see K96). However, in order to estimate the slope of the density profile at the positions of the two radial arcs, we fitted various polynomial functions to the density profile data. De-projecting this analysis to obtain the slope of the three-dimensional density profile ($`\rho `$), we obtain the gradients: $`d(\mathrm{log}\rho )/d(\mathrm{log}r)=1.41`$ and $`1.29`$ at $`r=6.0`$ kpc (B0b) and $`r=19.9`$ kpc (B1d) respectively; the error on these slopes is of the order $`<0.03`$.

We interpret this difference in slope at the positions of the two radial arcs, specifically the steeper slope at $`r6`$ kpc, as being due to the presence of a distinct central galaxy halo within the overall cluster halo.

Williams et al. (1999) suggest that radial arcs are a useful test of the lensing role of such central galaxy halos within an overall cluster-scale NFW halo. They predict that $`\mathrm{log}(\theta _r/\theta _t)`$ (where $`\theta _r`$ and $`\theta _t`$ are the angular positions of the radial and tangential arcs respectively) should decrease with decreasing cluster mass as the lensing role of the central galaxy becomes more significant relative to the overall cluster mass distribution. From their Fig. 7 we estimate that for A 383: $`\mathrm{log}(\theta _r/\theta _t)=0.35`$ to $`0.65`$, taking into account the velocity dispersions obtained from spectroscopy (§2.3) and the lens model (§3.2).

We have computed angular position ratios of $`log(\theta _r/\theta _t)=1.04\pm 0.07`$ and $`0.51\pm 0.03`$ for A 383’s inner (B0b) and outer (B1d) radial arcs respectively. The outer radial arc appears to be consistent with Williams et al.’s predictions, however the inner radial arc falls well outside their prediction. This discrepancy may be due to simplifying assumptions made by Williams et al. regarding the nature and complexity of cluster substructure. We also note that, due to the (almost) circular symmetry of A 383, the radial arcs are not forced to lie at the radial critical line in the image plane. This is another potential source of discrepancy between the observations and the predictions.

#### 3.2.2 Comparison With Numerical Simulations

Ghigna et al. (2000) have performed numerical simulations of galaxy clusters at the highest resolution to date, resolving structure down to scales of $`r5`$ kpc. They find a de-projected slope of $`1.6\pm 0.1`$ on scales of $`r=5`$–100 kpc in a $`\sigma 800`$ km s-1 cluster. This is steeper than the density profile of A 383 on comparable scales, $`1.3`$ (see §3.2.1 and Fig. 7). We will present a detailed comparison of the lens model using NFW and Moore et al. (1998) density profiles in a subsequent paper.

Recent progress has also been made in the inclusion of gas in numerical simulations (e.g. Pearce et al. 1999). However, the resolution of these simulations precludes analysis of the form of the mass profile on the scales probed by the radial arcs in A 383 ($`r<20`$ pc). Higher resolution simulations including gas cooling have been performed by Lewis et al. (1999). However, these authors concede that the mass of the central galaxies in their simulations is unrealistically high, which hampers a reliable comparison with our observations of A 383.

The constraints on the density profile derived from the radial arcs in A 383 confirm the need for a massive central galaxy to explain the lensing properties of this X-ray luminous cluster. Our observations therefore support the case for increasing the resolution of N-body simulations of clusters including gas physics.

### 3.3. X-ray Properties of A 383

The radial X-ray surface brightness profile around the cluster X-ray centroid is shown in Fig. 8. We attempt to describe the X-ray emission profile by fitting a standard beta model (Cavaliere & Fusco-Femiano 1978)

$$I(s)=I_0(1+(\frac{r}{r_\mathrm{c}^{3\mathrm{d}}})^2)^{3\beta +0.5},$$

convolved with the HRI point-spread function (PSF), to the observed X-ray surface brightness distribution. The free fit parameters are the peak surface brightness ($`I_0`$), the 3-d core radius ($`r_\mathrm{c}^{3\mathrm{d}}`$), and the slope parameter ($`\beta `$). We note that the use of the beta model implies that the gas is isothermal and in hydrostatic equilibrium. While these assumptions are likely to be violated in most clusters, numerical simulations show that using the beta model introduces, on average, no bias and only a moderate scatter of typically 20% in the cluster mass estimates (Evrard, Metzler & Navarro 1996).

We fit two beta models, one to the full radial range of the data including the likely cooling flow region (model A), and a second one to the outer regions only (model B). While model A yields a poor description of the X-ray emission at large radii (see Fig. 8) it provides us with an acceptable fit to, and an analytic description of, the emission profile at $`r<100`$ kpc, which can be used to derive the gas density and, assuming hydrostatic equilibrium, the total binding mass within the core region (see below). Model B, obtained by fitting only the range from $`r=12`$ kpc to 850 kpc ($`3.5^{}`$), describes the non-cooling-flow component (Fig. 8), for which we find $`r_\mathrm{c}^{3\mathrm{d}}=(23\pm 3)^{\prime \prime }`$, corresponding to $`(93\pm 12)`$ kpc, and $`\beta =0.65\pm 0.03`$, consistent with the canonical value of 2/3. The quoted errors include systematic uncertainties which we explored by varying the radial range used for the fit.

Fig. 8. – The radial X-ray surface brightness profile of A 383 as observed with the ROSAT HRI. Overlaid are the best-fit beta models (convolved with the HRI PSF) and fitted to the full profile (model A, dashed line) or to the outer regions only (model B, solid line). The excess emission in the cluster core is caused by an additional cooling flow component, which is not included in the models.

Subtracting beta model B from the observed profile we find the cooling flow to contribute 13% to the total cluster flux within a radius of 500 kpc, corresponding to a luminosity of $`1\times 10^{44}`$ erg s-1 (0.1–2.4 keV). A rough (usually low) estimate of the mass deposition rate, $`\dot{M}`$, can be obtained by converting this excess luminosity to the bolometric band (by multiplying by a factor of 3.03 appropriate for the assumed cluster gas temperature of 7.1 keV) and equating it to the luminosity of the cooling flow

$$L_{\mathrm{cool}}=\frac{5}{2}\frac{\dot{M}}{\mu m_\mathrm{p}}kT$$

(Fabian 1994). We thus derive a lower limit to the mass deposition rate of $`\dot{M}>200`$$`M_{}`$ yr-1, a value typical of moderate cooling-flow clusters (e.g. Peres et al. 1998).

Under the assumption of hydrostatic equilibrium, and using a beta model to describe the radial gas density profile, the total gravitational mass within a projected radius $`r`$ can be derived from

$$M(<r)=\frac{3\beta kTr_\mathrm{c}^{3\mathrm{d}}}{\mu m_\mathrm{p}G}\frac{\pi }{2}\frac{(r/r_\mathrm{c}^{3\mathrm{d}})^2}{\sqrt{1+(r/r_\mathrm{c}^{3\mathrm{d}})^2}}$$

(Soucail et al. 2000).

Using the parameters of model A we estimate the projected mass within a radius of 65 kpc to be $`4.0_{1.7}^{+1.1}\times 10^{13};M_{}`$, in good agreement with our lensing mass measurement of $`(3.5\pm 0.1)\times 10^{13}M_{}`$ within the same projected radius (§3.2). The equivalent mass estimate within a radius of 250 kpc is more uncertain, but is approximately $`(1.2\pm 0.5)\times 10^{14}M_{}`$, again in reasonable agreement with the lensing measurement. We note that the bulk of the uncertainty in the X-ray mass estimates is caused by the adopted error of $`\pm 2`$ keV in the global gas temperature, $`kT`$. Forthcoming observations of our full cluster sample with Newton will allow us to obtain much more accurate X-ray mass measurements on scales of $`10`$–1000 kpc.

## 4. Discussion and Conclusions

We have discovered numerous new gravitationally lensed features, including a giant arc and two radial arcs in A 383, a massive, X-ray luminous cluster of galaxies at a redshift of $`z=0.188`$. This is the first cluster in which two radial arcs have been identified and these provide a detailed view of the cluster mass distribution on scales of $`r<20`$ kpc.

The morphologies and positions of the arcs have been used, in conjunction with ground-based multi-color photometry, to constrain the mass distribution and compute a precise mass density profile of this cluster within a radius of $`250`$ kpc.

Within the core radius ($`r=46\pm 3`$ kpc) derived from the lensing analysis, the density profile exhibits a shallower slope ($`d(\mathrm{log}\rho )/d(\mathrm{log}r)=1.29\pm 0.03`$) than a Moore profile (Moore et al. 1998; Ghigna et al. 2000) and a steeper slope than a single NFW profile.

We have also analysed the properties of the two radial arcs, interpreting the difference in angular position of these arcs as a signature of the lensing role of the massive central galaxy. This supports the proposal by Williams et al. (1999) that massive central galaxy halos are required to explain the lensing properties of the cores of massive clusters ($`\sigma 1000`$ km s-1). However, only the outer radial arc of A 383 falls within the range predicted by Williams et al. for the the angular position ratio of radial and tangential arcs. The inner radial arc lies well outside of their predicted range. This discrepancy is probably caused by simplifying assumptions in Williams et al.’s theoretical models.

The lensing analysis is complemented by an analysis of the cluster’s X-ray properties as obtained from archival ROSAT HRI data. We find strong evidence for a cooling flow in agreement with the optical properties of the cluster, and derive a lower limit on the mass deposition rate of 200 M yr-1, typical of moderate cooling flow clusters. The X-ray estimates for the total mass within the projected radius of the giant arc are in good agreement with the lensing measurements and are consistent with earlier claims (e.g. Allen 1998) that significant discrepancies between lensing and X-ray mass estimates are caused by substructure in unrelaxed, non-cooling-flow systems.

We note that the cluster sample used in our survey is effectively X-ray selected, thus allowing us, for the first time, to measure in detail the mass distribution in massive clusters of galaxies unaffected by selection biases. The theme of our future program will be to analyse the lensing signal and X-ray emission (as observed with Newton) from all twelve clusters in our sample and to trace the form of the mass profile from 50 kpc to 5 Mpc. We will investigate the dispersion in the profiles and the core properties of the clusters and determine how these correlate with the cluster’s dynamical state and the presence of a cooling flow. We shall also attempt to constrain the high end of the cluster mass function and calibrate directly the cluster mass-temperature and mass-luminosity relations at $`z0.2`$.

## Acknowledgements

We thank John Blakeslee, Alastair Edge, Ben Moore, Frazer Pearce and Liliya Williams for helpful conversations and assistance. Thanks also go to Laurence Jones for helping with the Keck/MOS observations.

GPS acknowledges support from PPARC. JPK acknowledges support from CNRS. HE gratefully acknowledges financial support by NASA LTSA grant NAG 5–8253. OC acknowledges support from the European Commission under contract no. ER–BFM–BI–CT97–2471. IRS acknowledges support from a Royal Society University Research Fellowship. We also acknowledge financial support from the UK-French ALLIANCE collaboration programme #00161XM.

UKIRT is operated by the Joint Astronomy Centre on behalf of the Particle Physics and Astronomy Research Council of the United Kingdom. The W. M. Keck Observatory is operated as a scientific partnership among the California Institute of Technology, the University of California and NASA.


text:
# Quantum nonlinear dynamics of continuously measured systems

## I Introduction

The Hamiltonian formulation of classical mechanics assigns simultaneous, arbitrarily accurate, values for the canonically conjugate position and momentum to distinguishable particles. Indeed in classical mechanics these simultaneous values are regarded as properties of the particles themselves; measurement simply reveals these values and need not, in principle, add uncertainty to their determination. In quantum mechanics the conjugate position and momentum are represented by noncommuting operators, $`\widehat{x}`$ and $`\widehat{p}`$, with $`[\widehat{x},\widehat{p}]=i\mathrm{}`$. It follows that there is no physical state for which position and momentum can take dispersion free values. This is the content of the standard formulation of the uncertainty principle $`(\mathrm{\Delta }\widehat{x})^2(\mathrm{\Delta }\widehat{p})^2\mathrm{}^2/4`$, where $`\mathrm{\Delta }\widehat{A}=\widehat{A}\widehat{A}`$. This does not mean, however, that simultaneous measurements of position and momentum are impossible, only that such simultaneous measurements cannot be made arbitrarily accurate Busch85 .

Our objective is to show that the observed phase space, reconstructed from the classical stochastic measurement record of continuous joint measurements of position and momentum, corresponds to the classical phase space description with added noise. By a careful specification of the measurement model, including the relevant states of the apparatus, together with a Markov assumption, we show that the observed classical stochastic measurement record corresponds to the conditional quantum averages $`\widehat{x}(t)`$ and $`\widehat{p}(t)`$. We obtain a stochastic Schrödinger equation for the conditional state of the system with complex noise, together with classical stochastic equations for the variables corresponding to the observed measurement records. This extends the results of Bhattacharya et al. habib  where only a position measurement was considered. Other authors gisin0 ; gisin ; gisin2 ; spiller ; schack ; brun  have considered a general type of stochastic Schrödinger equation (which may be interpreted as corresponding to a measurement) in the framework of quantum state diffusion to examine dissipative chaos in open systems.

A model for the simultaneous measurement of position and momentum was given long ago by Arthurs and Kelly arthurs . This model consists of two meters which are allowed to interact instantaneously with the system. The interaction couples one of the meters to position and the other to momentum, encoding the results of the measurement in the final states of the meters. Projective measurements are then made on each of the meter states separately. In more general terms the Arthurs-Kelly model is an explicit example of the Gelfand-Naimark-Segal theorem Conway  which enables us to write a positive joint probability density of two real variables (the actual measurement outcomes) in terms of an inner product on an extended Hilbert space (extended to include the two apparatus as well as the system to be measured). It can also be considered as an example of a positive operator valued measure POVM1 ; POVM2 . The Arthurs-Kelly model suggests a route to a phase space description of quantum mechanics based on the joint outcomes of simultaneous measurements of position and momentum. However, to achieve a phase space description of the dynamics it is necessary to consider a time sequence of such measurements. The Arthurs-Kelly model forces the conditional state of the system into a coherent state after a single measurement and is too strong for this purpose.

We generalize the Arthurs-Kelly model to allow for a weakening of the measurement. In this way repeated measurements can take place without invariably reducing the system to a coherent state. By simultaneously increasing the number of measurements and weakening the strength of each measurement, the continuous limit is achieved. We show that under continuous observation the evolution of the system state is described by an Itô stochastic Schrödinger equation. Unlike evolution via the ordinary Schrödinger equation, continuous observation forces the conditional state of the system to remain localized. Hence the quantum mean of the phase space variables can be thought of as a trajectory. It is in this light that we investigate quantum nonlinear dynamics under continuous observation.

The theory of continuous measurement in quantum mechanics was first developed by Barchielli and co-workers barchielli1 ; barchielli2 ; barchielli3  (see also diosi ; belavkin ) and some of the results in this paper are contained in their work. However our approach uses the Arthurs-Kelly measurement model as a starting point and follows the work on continuous position measurements by Caves and one of us in caves ; milburn . In this way our results have a clear physical interpretation.

In section II we generalize the Arthurs-Kelly measurement model. In section III we derive the stochastic Schrödinger equation for continuous measurement which is then investigated analytically in section IV and simulated numerically in section V for a chaotic system. Finally, in section VI we discuss our results.

## II Phase space measurements

Following Braunstein et al. braunstein  we rederive and generalize the result implicit in the work of Arthurs and Kelly arthurs . The simultaneous measurement of position and momentum of a one dimensional quantum system is achieved by constructing an interaction governed by the Hamiltonian

$$\widehat{H}_T=\widehat{H}(\widehat{x},\widehat{p})+\left(\frac{1}{s}\widehat{x}\widehat{p}_1+s\widehat{p}\widehat{p}_2\right)\delta (tt_r)$$
(1)
which couples the system Hamiltonian $`\widehat{H}`$, with two detectors $`d_1`$ and $`d_2`$. The measurement model also requires that the detectors have been prepared in the initial states

$`x_i|d_i`$ $`=`$ $`(\pi \mathrm{}\mathrm{\Delta }_i)^{\frac{1}{4}}\mathrm{exp}\left({\displaystyle \frac{x_i^2}{2\mathrm{}\mathrm{\Delta }_i}}\right)`$ (2)
$`p_i|d_i`$ $`=`$ $`\left({\displaystyle \frac{\pi \mathrm{}}{\mathrm{\Delta }_i}}\right)^{\frac{1}{4}}\mathrm{exp}\left({\displaystyle \frac{p_i^2\mathrm{\Delta }_i}{2\mathrm{}}}\right)`$ (3)
where $`\mathrm{\Delta }_1=\frac{\sigma ^2}{2s}`$ and $`\mathrm{\Delta }_2=\frac{s\sigma ^2}{2}`$. The parameter $`s`$ is called the squeezing parameter and biases the coupling so that one may obtain more information on either position or momentum. The parameter $`\sigma `$ will be used to weaken the measurement, decreasing the amount of information collected on both position and momentum.

Before the interaction the system state $`\widehat{\rho }`$ and the combined detector state

$$\widehat{\rho }_d=|d_1d_2d_1d_2||d_1d_1||d_2d_2|$$
(4)
are assumed to be uncorrelated with an initial combined density operator of $`\widehat{\rho }\widehat{\rho }_d`$. At $`t=t_r`$ the evolution operator for the interaction

$$\widehat{U}_I\mathrm{exp}\left[\frac{i}{\mathrm{}}\left(\frac{1}{s}\widehat{x}\widehat{p}_1+s\widehat{p}\widehat{p}_2\right)\right]$$
(5)
couples the system to the measurement apparatus, entangling the system and detector states. After the interaction the probability of finding the detector positions in the small area $`[x_1,x_1+dx_1]\times [x_2,x_2+dx_2]`$ is

$`\text{Prob}(x_1,x_2)dx_1dx_2`$ $`=`$ $`\text{tr}\text{tr}_d\left(\widehat{U}_I\widehat{\rho }\widehat{\rho }_d\widehat{U}_I^{}|x_1x_2x_1x_2|\right)dx_1dx_2`$ (6)
$`=`$ $`\text{tr}\left(\widehat{\mathrm{{\rm Y}}}_\sigma (x_1,x_2)\widehat{\rho }\widehat{\mathrm{{\rm Y}}}_\sigma (x_1,x_2)^{}\right)dx_1dx_2`$ (7)
where $`\text{tr}(\mathrm{})`$ and $`\text{tr}_d(\mathrm{})`$ are the traces over the system and detector states, respectively, $`\text{Prob}(x_1,x_2)`$ is the probability density, and the resolution operator

$`\widehat{\mathrm{{\rm Y}}}_\sigma (x_1,x_2)`$ $``$ $`x_1x_2|\widehat{U}_I|d_1d_2`$ (8)
$`=`$ $`(2\pi \mathrm{})^1{\displaystyle 𝑑p_1𝑑p_2\mathrm{exp}\left[\frac{i}{\mathrm{}}\left(p_1\left(\frac{1}{s}\widehat{x}x_1\right)+p_2\left(s\widehat{p}x_2\right)\right)\right]p_1p_2|d_1d_2}`$ (9)
$`=`$ $`\widehat{D}(sx_1,\frac{1}{s}x_2)\widehat{\mathrm{{\rm Y}}}_\sigma (0,0)\widehat{D}(sx_1,\frac{1}{s}x_2)^{}.`$ (10)
The displacement operator

$$\widehat{D}(\mu ,\nu )\mathrm{exp}\left[\frac{i}{\mathrm{}}\left(\mu \widehat{p}\nu \widehat{x}\right)\right]$$
(11)
and

$$\widehat{\mathrm{{\rm Y}}}_\sigma (0,0)=(2\pi \mathrm{})^1𝑑p_1𝑑p_2p_1,p_2|d_1d_2\widehat{D}(sp_2,\frac{1}{s}p_1)^{}.$$
(12)
We now define the annihilation operator

$$\widehat{a}\frac{1}{\sqrt{2\mathrm{}}}\left(\sqrt{s}\widehat{x}+i\frac{1}{\sqrt{s}}\widehat{p}\right)$$
(13)
which satisfies

$$[\widehat{a},\widehat{a}^{}]=1$$
(14)
and rewrite the displacement operator as

$$\widehat{D}(\zeta )\mathrm{exp}\left(\zeta \widehat{a}^{}\zeta ^{}\widehat{a}\right)$$
(15)
where $`\zeta =\frac{1}{\sqrt{2\mathrm{}}}\left(\sqrt{s}\mu +i\frac{1}{\sqrt{s}}\nu \right)`$. It can be easily shown that

$`\widehat{D}(\zeta )^{}=\widehat{D}(\zeta )=\widehat{D}(\zeta )^1`$ (16)
$`\widehat{D}(\zeta )\widehat{a}\widehat{D}(\zeta )^{}=\widehat{a}\zeta .`$ (17)
The coherent states perelomov  are now defined by applying the displacement operator onto the vacuum state

$$|\alpha \widehat{D}(\alpha )|0.$$
(18)
They satisfy the following relations

$`\widehat{a}|\alpha `$ $`=\alpha |\alpha `$ (19)
$`\widehat{D}(\beta )|\alpha `$ $`=\mathrm{exp}\left[\frac{1}{2}\left(\alpha ^{}\beta \alpha \beta ^{}\right)\right]|\alpha +\beta `$ (20)
$`\alpha |\beta `$ $`=\mathrm{exp}\left(\frac{1}{2}|\alpha |^2\frac{1}{2}|\beta |^2+\alpha ^{}\beta \right)`$ (21)
$`\pi ^1{\displaystyle d^2\alpha }`$ $`|\alpha \alpha |=1.`$ (22)
Now defining $`\varrho \frac{1}{\sqrt{2\mathrm{}}}\left(\sqrt{s}p_2+i\frac{1}{\sqrt{s}}p_1\right)`$ we can rewrite (12) as

$$\widehat{\mathrm{{\rm Y}}}_\sigma (0,0)=(2\pi \mathrm{})^{\frac{1}{2}}\sigma \pi ^1d^2\varrho \sigma \varrho |0\widehat{D}(\varrho )^{}.$$
(23)
Using the above relations and integrating we obtain

$`\alpha |\widehat{\mathrm{{\rm Y}}}_\sigma (0,0)|\beta `$ $`=`$ $`(2\pi \mathrm{})^{\frac{1}{2}}\frac{2\sigma }{\sigma ^2+1}\mathrm{exp}\left(\frac{1}{2}|\alpha |^2\frac{1}{2}|\beta |^2+\frac{\sigma ^21}{\sigma ^2+1}\alpha ^{}\beta \right)`$ (24)
$`=`$ $`(2\pi \mathrm{})^{\frac{1}{2}}\frac{2\sigma }{\sigma ^2+1}\alpha |\beta \mathrm{exp}\left(\frac{2}{\sigma ^2+1}\alpha ^{}\beta \right)`$ (25)
$`=`$ $`\alpha |(2\pi \mathrm{})^{\frac{1}{2}}\frac{2\sigma }{\sigma ^2+1}:\mathrm{exp}\left(\frac{2}{\sigma ^2+1}\widehat{a}^{}\widehat{a}\right):|\beta `$ (26)
where $`::`$ denotes normal ordering. Hence,

$$\widehat{\mathrm{{\rm Y}}}_\sigma (0,0)=(2\pi \mathrm{})^{\frac{1}{2}}\frac{2\sigma }{\sigma ^2+1}:\mathrm{exp}\left(\frac{2}{\sigma ^2+1}\widehat{a}^{}\widehat{a}\right):.$$
(27)
If we now define $`\chi =\chi _1+i\chi _2\frac{1}{\sqrt{2\mathrm{}}}\left(\sqrt{s}x_1+i\frac{1}{\sqrt{s}}x_2\right)`$ then the resolution operator becomes

$`\widehat{\mathrm{{\rm Y}}}_\sigma (x_1,x_2)`$ $`=`$ $`\widehat{D}(\chi )\widehat{\mathrm{{\rm Y}}}_\sigma (0,0)\widehat{D}(\chi )^{}`$ (28)
$`=`$ $`(2\pi \mathrm{})^{\frac{1}{2}}\frac{2\sigma }{\sigma ^2+1}:\mathrm{exp}\left[\frac{2}{\sigma ^2+1}(\widehat{a}^{}\chi ^{})(\widehat{a}\chi )\right]:`$ (29)
$``$ $`(2\mathrm{})^{\frac{1}{2}}\widehat{\mathrm{{\rm Y}}}_\sigma (\chi ).`$ (30)
The probability of finding the detector positions in the small area $`[\chi _1,\chi _1+d\chi _1]\times [\chi _2,\chi _2+d\chi _2]`$ is now

$`\text{Prob}(\chi )d^2\chi `$ $`=`$ $`2\mathrm{}\text{Prob}(x_1,x_2)d^2\chi `$ (31)
$`=`$ $`\text{tr}\left(\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )\widehat{\rho }\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )^{}\right)d^2\chi `$ (32)
$`=`$ $`\text{tr}\left(\widehat{F}_\sigma (\chi )\widehat{\rho }\right)d^2\chi `$ (33)
where $`d^2\chi =d\chi _1d\chi _2`$ and

$`\widehat{F}_\sigma (\chi )`$ $``$ $`\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )^{}\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )`$ (34)
$`=`$ $`\pi ^1\left(\frac{2\sigma }{\sigma ^2+1}\right)^2:\mathrm{exp}\left[\left(\frac{2\sigma }{\sigma ^2+1}\right)^2(\widehat{a}^{}\chi ^{})(\widehat{a}\chi )\right]:`$ (35)
is an effect density POVM1 . It can be easily shown that

$`{\displaystyle \chi ^n\widehat{F}_\sigma (\chi )d^2\chi }`$ $`=`$ $`\widehat{a}^n`$ (36)
$`{\displaystyle |\chi |^2\widehat{F}_\sigma (\chi )d^2\chi }`$ $`=`$ $`\widehat{a}\widehat{a}^{}+\left(\frac{\sigma ^21}{2\sigma }\right)^2.`$ (37)
Hence, defining the notion of a mean for this measurement process

$`f(\chi )_\sigma `$ $``$ $`{\displaystyle f(\chi )\text{Prob}(\chi )d^2\chi }`$ (38)
$`=`$ $`{\displaystyle f(\chi )\text{tr}\left(\widehat{F}_\sigma (\chi )\widehat{\rho }\right)d^2\chi }`$ (39)
we find that

$`\chi _\sigma `$ $`=`$ $`\widehat{a}`$ (40)
$`|\chi |^2_\sigma `$ $`=`$ $`\widehat{a}\widehat{a}^{}+\left(\frac{\sigma ^21}{2\sigma }\right)^2`$ (41)
or

$`x_1_\sigma =\widehat{x},`$ $`x_2_\sigma =\widehat{p}`$ (42)
$`x_{1}^{}{}_{}{}^{2}_\sigma =\widehat{x}^2`$ $`+\frac{\mathrm{}}{s}\left(\frac{1+\sigma ^4}{4\sigma ^2}\right)`$ (43)
$`x_{2}^{}{}_{}{}^{2}_\sigma =\widehat{p}^2`$ $`+\mathrm{}s\left(\frac{1+\sigma ^4}{4\sigma ^2}\right)`$ (44)
where $`\widehat{A}=\text{tr}(\widehat{A}\widehat{\rho })`$ is the quantum expectation. Thus the readout variables $`x_1`$ and $`x_2`$ give, respectively, the position and momentum of the system with additional noise dependent on $`\sigma `$. We obtain maximal information from the system when the variances are at a minimum. That is, when $`\sigma =1`$. In this case

$`\widehat{F}_1(\chi )`$ $`=`$ $`\pi ^1:\mathrm{exp}\left[(\widehat{a}^{}\chi ^{})(\widehat{a}\chi )\right]:`$ (45)
$`=`$ $`\pi ^1|\chi \chi |`$ (46)
and the probability density reduces to the Husimi density or $`Q`$ function

$`\text{Prob}(\chi )`$ $`=`$ $`\pi ^1\chi |\widehat{\rho }|\chi `$ (47)
$`=`$ $`\pi ^1Q(\chi ).`$ (48)
Suppose we take a measurement and obtain the outcome $`\chi ^{}`$. As a consequence of the strength of this measurement the system state collapses to

$`\widehat{\rho }^{}`$ $`=`$ $`{\displaystyle \frac{1}{\text{Prob}(\chi ^{})}}\widehat{\mathrm{{\rm Y}}}_1(\chi ^{})\widehat{\rho }\widehat{\mathrm{{\rm Y}}}_1(\chi ^{})^{}`$ (49)
$`=`$ $`{\displaystyle \frac{1}{\chi ^{}|\widehat{\rho }|\chi ^{}}}|\chi ^{}\chi ^{}|\widehat{\rho }|\chi ^{}\chi ^{}|`$ (50)
$`=`$ $`|\chi ^{}\chi ^{}|.`$ (51)
However when $`\sigma 1`$ the resolution operator has the asymptotic expansion

$`\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )`$ $`=`$ $`\pi ^{\frac{1}{2}}\frac{2}{\sigma }\left[1\frac{1}{\sigma ^2}\left(1+2(\widehat{a}^{}\chi ^{})(\widehat{a}\chi )\right)\right]+O(\frac{1}{\sigma ^5}).`$ (52)
Using this result one can show that the system state conditioned on the measurement outcome $`\chi ^{}`$ has the expansion

$$\widehat{\rho }^{}=\widehat{\rho }+\frac{2}{\sigma ^2}\{\text{tr}\left[(\widehat{a}^{}\chi ^{})(\widehat{a}\chi ^{})\widehat{\rho }\right](\widehat{a}^{}\chi ^{})(\widehat{a}\chi ^{}),\widehat{\rho }\}+O(\frac{1}{\sigma ^4})$$
(53)
and thus, if $`\sigma `$ is large enough, the process of measurement will have negligible effect on the system.

## III Continuous measurement

Consider a sequence of phase space measurements governed by the Hamiltonian

$$\widehat{H}_T=\widehat{H}(\widehat{x},\widehat{p})+\left(\frac{1}{s}\widehat{x}\widehat{p}_1+s\widehat{p}\widehat{p}_2\right)\underset{n=0}{\overset{\mathrm{}}{}}\delta (tn\delta t)$$
(54)
where after each measurement the detectors are reset into the initial states given by equation (2) or (3). The assumption that the detectors are reset is equivalent to making a Markov assumption for a single apparatus coupled to the system. By resetting the detector at each time step we ensure that no coherent memory of the system state survives in the states of the apparatus. Following caves  we will first derive the master equation for unconditional (or nonselective) evolution of the system density operator in the continuous limit $`\delta t0`$, $`\sigma \mathrm{}`$. By unconditional evolution we mean that no account is taken of the measured results. Thus after each measurement occurs we ignore the result and average over all possible measurement outcomes. If we denote the system density operator immediately before the $`n`$-th measurement by $`\widehat{\rho }(n\delta t)`$ then

$$\widehat{\rho }(n\delta t+\delta t)=\widehat{U}d^2\chi \widehat{\mathrm{{\rm Y}}}_\sigma (\chi )\widehat{\rho }(n\delta t)\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )^{}\widehat{U}^{}$$
(55)
where

$`\widehat{U}`$ $``$ $`\mathrm{exp}\left(\frac{i}{\mathrm{}}\widehat{H}\delta t\right)`$ (56)
$`=`$ $`1\frac{i}{\mathrm{}}\widehat{H}\delta t+O(\delta t^2).`$ (57)
For any operator $`\widehat{A}`$ it is possible to show that

$$d^2\chi \widehat{\mathrm{{\rm Y}}}_\sigma (\chi )\widehat{A}\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )^{}=\widehat{A}\frac{1}{\sigma ^2}[\widehat{a},[\widehat{a}^{},\widehat{A}]]+O(\frac{1}{\sigma ^4}).$$
(58)
Hence we obtain

$$\frac{\widehat{\rho }(n\delta t+\delta t)\widehat{\rho }(n\delta t)}{\delta t}=\frac{i}{\mathrm{}}[\widehat{H},\widehat{\rho }(n\delta t)]\frac{1}{\delta t\sigma ^2}[\widehat{a},[\widehat{a}^{},\widehat{\rho }(n\delta t)]]+O(\delta t)+O(\frac{1}{\sigma ^2})+O(\frac{1}{\delta t\sigma ^4}).$$
(59)
By setting $`t=n\delta t`$ and taking the continuous limit $`\delta t0`$, $`\sigma \mathrm{}`$, with $`\gamma =1/\delta t\sigma ^2`$ held constant, we obtain the master equation for unconditional evolution:

$`{\displaystyle \frac{d\widehat{\rho }}{dt}}`$ $`=`$ $`\frac{i}{\mathrm{}}[\widehat{H},\widehat{\rho }]\gamma [\widehat{a},[\widehat{a}^{},\widehat{\rho }]]`$ (60)
$`=`$ $`\frac{i}{\mathrm{}}[\widehat{H},\widehat{\rho }]\frac{1}{2\mathrm{}}\mathrm{\Gamma }_1[\widehat{x},[\widehat{x},\widehat{\rho }]]\frac{1}{2\mathrm{}}\mathrm{\Gamma }_2[\widehat{p},[\widehat{p},\widehat{\rho }]]`$ (61)
where $`\mathrm{\Gamma }_1\gamma s`$ and $`\mathrm{\Gamma }_2\frac{\gamma }{s}`$. This equation has already been derived by Barchielli et al. barchielli1 , but their approach is somewhat different. By setting $`\mathrm{\Gamma }_2=0`$ in (61) we obtain the unconditional master equation for continuous position measurements previously derived in caves .

We now wish to derive the conditional (or selective) master equation for the system density operator. In this case the evolution of the system is conditioned on a history of measurement readouts

$$\left\{\chi (n\delta t)\right\}=\{\chi (0),\chi (\delta t),\chi (2\delta t),\mathrm{}\}$$
(62)
where $`\chi (n\delta t)`$ is the detector position for the $`n`$-th measurement. Hence, if $`\widehat{\rho }(n\delta t)`$ is the system density operator immediately before the $`n`$-th measurement, then

$$\widehat{\rho }(n\delta t+\delta t)=\widehat{U}\left[\text{Prob}(\chi (n\delta t))\right]^1\widehat{\mathrm{{\rm Y}}}_\sigma (\chi (n\delta t))\widehat{\rho }(n\delta t)\widehat{\mathrm{{\rm Y}}}_\sigma (\chi (n\delta t))^{}\widehat{U}^{}.$$
(63)
To proceed we extend the definition of the readout variable by setting

$$\chi (t)=\chi (n\delta t)\text{ for }n\delta tt<(n+1)\delta t$$
(64)
and introduce the new variable

$$𝒳(t)_0^t\chi (t^{})𝑑t^{}.$$
(65)
Hence $`𝒳(0)=0`$ and

$$𝒳(n\delta t)=\delta t\underset{m=0}{\overset{n1}{}}\chi (m\delta t)\text{ for }n1.$$
(66)
Using (40) and (41) we obtain

$`E_c\left(\chi (n\delta t)\right)`$ $`=`$ $`\chi (n\delta t)_\sigma `$ (67)
$`=`$ $`\text{tr}\left(\widehat{a}\widehat{\rho }(n\delta t)\right)`$ (68)
$`V_c\left(\chi (n\delta t)\right)`$ $`=`$ $`\left|\chi (n\delta t)\chi (n\delta t)_\sigma \right|^2_\sigma `$ (69)
$`=`$ $`|\chi (n\delta t)|^2_\sigma \left|\chi (n\delta t)_\sigma \right|^2`$ (70)
$`=`$ $`\text{tr}\left(\widehat{a}\widehat{a}^{}\widehat{\rho }(n\delta t)\right)\left|\text{tr}\left(\widehat{a}\widehat{\rho }(n\delta t)\right)\right|^2+\left(\frac{\sigma ^21}{2\sigma }\right)^2.`$ (71)
where the subscript $`c`$ has been added to emphasize that the mean and variance are conditioned through $`\widehat{\rho }`$ on the entire history of measurement readouts. Now letting

$$\delta 𝒳(n\delta t)=𝒳(n\delta t+\delta t)𝒳(n\delta t)=\delta t\chi (n\delta t)$$
(72)
we find that

$`E_c\left(\delta 𝒳(n\delta t)\right)`$ $`=`$ $`\delta t\text{tr}\left(\widehat{a}\widehat{\rho }(n\delta t)\right)`$ (73)
$`V_c\left(\delta 𝒳(n\delta t)\right)`$ $`=`$ $`\delta t^2\left(\text{tr}\left(\widehat{a}\widehat{a}^{}\widehat{\rho }(n\delta t)\right)\left|\text{tr}\left(\widehat{a}\widehat{\rho }(n\delta t)\right)\right|^2+\left(\frac{\sigma ^21}{2\sigma }\right)^2\right)`$ (74)
$`=`$ $`\frac{1}{4}\gamma ^1\delta t+O(\delta t^2).`$ (75)
where we have set $`\gamma =1/\delta t\sigma ^2`$ to be constant in anticipation of the continuous limit. Hence for $`\delta t`$ small enough we can approximate increments in the variable $`𝒳`$ by

$$\delta 𝒳(n\delta t)=\text{tr}\left(\widehat{a}\widehat{\rho }(n\delta t)\right)\delta t+\frac{1}{2}\gamma ^{\frac{1}{2}}\delta \xi (n\delta t)$$
(76)
where it is understood that the complex Itô increment $`\delta \xi `$ is of order $`\delta t^{\frac{1}{2}}`$ and satisfies $`E\left(\delta \xi \right)=0`$, $`E\left(\delta \xi ^{}\delta \xi \right)=\delta t`$. In the continuous limit $`\delta t0`$ with $`t=n\delta t`$ constant we have

$$d𝒳(t)=\text{tr}\left(\widehat{a}\widehat{\rho }(t)\right)dt+\frac{1}{2}\gamma ^{\frac{1}{2}}d\xi (t)$$
(77)
where $`\xi (t)`$ is a complex Wiener process gardiner  and the Itô differential $`d\xi `$ satisfies the algebra

$`E\left(d\xi \right)`$ $`=`$ $`0`$ (78)
$`d\xi ^{}d\xi `$ $`=`$ $`dt`$ (79)
$`d\xi ^2`$ $`=`$ $`0`$ (80)
$`d\xi dt`$ $`=`$ $`0.`$ (81)
To simplify the following we will always replace $`\delta \xi ^{}\delta \xi `$ by $`\delta t`$ and set $`\delta \xi ^2=\delta \xi ^2=0`$ in anticipation of the above algebra in the continuous limit. Using (72) together with (76) we obtain the following asymptotic expansion for the resolution operator (30)

$$\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )=2\left(\frac{\gamma \delta t}{\pi e}\right)^{\frac{1}{2}}\left(1+\gamma ^{\frac{1}{2}}\left(\widehat{𝒜}^{}\delta \xi +\widehat{𝒜}\delta \xi ^{}\right)\gamma \delta t\left(\widehat{𝒜}^{}\widehat{𝒜}+\frac{1}{2}\right)\right)+O(\delta t^2)$$
(82)
where

$$\widehat{𝒜}\widehat{a}\text{tr}\left(\widehat{a}\widehat{\rho }\right)$$
(83)
and it is understood that $`\widehat{\rho }=\widehat{\rho }(n\delta t)`$, $`\chi =\chi (n\delta t)`$ and $`\delta \xi =\delta \xi (n\delta t)`$. Hence we find that

$$\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )\widehat{\rho }\widehat{\mathrm{{\rm Y}}}_\sigma (\chi )^{}=4\frac{\gamma \delta t}{\pi e}\left(\widehat{\rho }+\gamma ^{\frac{1}{2}}\{\widehat{𝒜}^{}\delta \xi +\widehat{𝒜}\delta \xi ^{},\widehat{\rho }\}\gamma \delta t[\widehat{𝒜},[\widehat{𝒜}^{},\widehat{\rho }]]\right)+O(\delta t^3)$$
(84)
and thus, using (57) we obtain

$`\widehat{\rho }(n\delta t+\delta t)`$ $`=`$ $`\widehat{\rho }\frac{i\delta t}{\mathrm{}}[\widehat{H},\widehat{\rho }]+\gamma ^{\frac{1}{2}}\{\widehat{𝒜}^{}\delta \xi +\widehat{𝒜}\delta \xi ^{},\widehat{\rho }\}\gamma \delta t[\widehat{𝒜},[\widehat{𝒜}^{},\widehat{\rho }]]+O(\delta t^{\frac{3}{2}})`$ (85)
$`=`$ $`\widehat{\rho }\frac{i\delta t}{\mathrm{}}[\widehat{H},\widehat{\rho }]\gamma \delta t[\widehat{a},[\widehat{a}^{},\widehat{\rho }]]+\gamma ^{\frac{1}{2}}[\widehat{a}^{}]\widehat{\rho }\delta \xi +\gamma ^{\frac{1}{2}}[\widehat{a}]\widehat{\rho }\delta \xi ^{}+O(\delta t^{\frac{3}{2}})`$ (86)
where we have defined the superoperator

$$[\widehat{A}]\widehat{\rho }\{\widehat{A}\text{tr}(\widehat{A}\widehat{\rho }),\widehat{\rho }\}.$$
(87)
In the limit $`\delta t0`$ we obtain the master equation for conditional evolution:

$`d\widehat{\rho }(t)`$ $`=`$ $`\widehat{\rho }(t+dt)\widehat{\rho }(t)`$ (88)
$`=`$ $`\frac{i}{\mathrm{}}[\widehat{H},\widehat{\rho }(t)]dt\gamma [\widehat{a},[\widehat{a}^{},\widehat{\rho }(t)]]dt+\gamma ^{\frac{1}{2}}[\widehat{a}^{}]\widehat{\rho }(t)d\xi (t)+\gamma ^{\frac{1}{2}}[\widehat{a}]\widehat{\rho }(t)d\xi (t)^{}.`$ (89)
It is easy to see that upon averaging this stochastic differential equation we reproduce our original master equation for unconditional evolution (60). However, note that unlike the unconditional equation, this equation preserves the pure state property of $`\widehat{\rho }`$. One can easily prove this by showing $`d\widehat{\rho }^2=d\widehat{\rho }`$ under the assumption that $`\widehat{\rho }=|\psi \psi |`$, where $`d\widehat{\rho }^2\{d\widehat{\rho },\widehat{\rho }\}+d\widehat{\rho }d\widehat{\rho }`$. As a consequence of this fact, the above master equation has an analogue for pure state evolution in terms of a stochastic Schrödinger equation:

$`d|\psi (t)=\frac{i}{\mathrm{}}\widehat{H}|\psi (t)dt\gamma `$ $`\left(\widehat{a}^{}\widehat{a}+\frac{1}{2}\widehat{a}^{}\widehat{a}\widehat{a}^{}\widehat{a}+|\widehat{a}|^2\right)|\psi (t)dt`$ (90)
$`+\gamma ^{\frac{1}{2}}\left(\widehat{a}^{}\widehat{a}^{}\right)|\psi (t)d\xi (t)+\gamma ^{\frac{1}{2}}\left(\widehat{a}\widehat{a}\right)|\psi (t)d\xi (t)^{}`$
where $`\widehat{A}=\psi (t)|\widehat{A}|\psi (t)`$. In terms of position and momentum variables the master equation for conditional evolution reads as

$$d\widehat{\rho }=\frac{i}{\mathrm{}}[\widehat{H},\widehat{\rho }]dt\frac{1}{2\mathrm{}}\mathrm{\Gamma }_1[\widehat{x},[\widehat{x},\widehat{\rho }]]dt\frac{1}{2\mathrm{}}\mathrm{\Gamma }_2[\widehat{p},[\widehat{p},\widehat{\rho }]]dt+\mathrm{}^{\frac{1}{2}}\mathrm{\Gamma }_{1}^{}{}_{}{}^{\frac{1}{2}}[\widehat{x}]\widehat{\rho }dW_1+\mathrm{}^{\frac{1}{2}}\mathrm{\Gamma }_{2}^{}{}_{}{}^{\frac{1}{2}}[\widehat{p}]\widehat{\rho }dW_2$$
(91)
where

$`E\left(dW_i\right)`$ $`=`$ $`0`$ (92)
$`dW_idW_j`$ $`=`$ $`\delta _{ij}dt`$ (93)
$`\frac{1}{\sqrt{2}}\left(dW_1+idW_2\right)`$ $`=`$ $`d\xi `$ (94)
and the readout variables $`x_1`$ and $`x_2`$ obey the following stochastic processes

$`dX_1(t)=\text{tr}\left(\widehat{x}\widehat{\rho }(t)\right)dt+`$ $`\frac{1}{2}\mathrm{}^{\frac{1}{2}}\mathrm{\Gamma }_{1}^{}{}_{}{}^{\frac{1}{2}}dW_1(t)`$ (95)
$`dX_2(t)=\text{tr}\left(\widehat{p}\widehat{\rho }(t)\right)dt+`$ $`\frac{1}{2}\mathrm{}^{\frac{1}{2}}\mathrm{\Gamma }_{2}^{}{}_{}{}^{\frac{1}{2}}dW_2(t)`$ (96)
with

$$X_1(t)_0^tx_1(t^{})𝑑t^{}X_2(t)_0^tx_2(t^{})𝑑t^{}.$$
(97)
By setting $`\mathrm{\Gamma }_2=0`$ in (91) we obtain the conditional master equation for continuous position measurements previously derived in milburn . Note that $`x_1`$ and $`x_2`$ are charged by stationary white noise

$`x_1=\widehat{x}+`$ $`\frac{1}{2}\mathrm{}^{\frac{1}{2}}\mathrm{\Gamma }_{1}^{}{}_{}{}^{\frac{1}{2}}\dot{W_1}`$ (98)
$`x_2=\widehat{p}+`$ $`\frac{1}{2}\mathrm{}^{\frac{1}{2}}\mathrm{\Gamma }_{2}^{}{}_{}{}^{\frac{1}{2}}\dot{W_2}`$ (99)
making their graph highly irregular. It is thus better to represent the measured trajectory by $`\widehat{x}`$ and $`\widehat{p}`$.

## IV Analytical investigations

We will now investigate the effect of measurement on the system state by setting $`\widehat{H}=0`$ and defining

$`V_x`$ $``$ $`E\left(\widehat{x}^2\widehat{x}^2\right)`$ (100)
$`V_p`$ $``$ $`E\left(\widehat{p}^2\widehat{p}^2\right)`$ (101)
$`C_{xp}`$ $``$ $`E\left(\frac{1}{2}\widehat{x}\widehat{p}+\frac{1}{2}\widehat{p}\widehat{x}\widehat{x}\widehat{p}\right)`$ (102)
where $`V_x`$ and $`V_p`$ are the expected variances in position and momentum, and $`C_{xp}`$ is the expected covariance between position and momentum. The average is taken over all possible measurement histories. One can then derive the following set of coupled differential equations:

$`{\displaystyle \frac{dV_x}{dt}}`$ $`=`$ $`{\displaystyle \frac{\gamma \mathrm{}}{s}}{\displaystyle \frac{4s\gamma }{\mathrm{}}}V_{x}^{}{}_{}{}^{2}{\displaystyle \frac{4\gamma }{s\mathrm{}}}C_{xp}^{}{}_{}{}^{2}`$ (103)
$`{\displaystyle \frac{dV_p}{dt}}`$ $`=`$ $`s\gamma \mathrm{}{\displaystyle \frac{4\gamma }{s\mathrm{}}}V_{p}^{}{}_{}{}^{2}{\displaystyle \frac{4s\gamma }{\mathrm{}}}C_{xp}^{}{}_{}{}^{2}`$ (104)
$`{\displaystyle \frac{dC_{xp}}{dt}}`$ $`=`$ $`{\displaystyle \frac{4\gamma }{\mathrm{}}}C_{xp}\left(sV_x+{\displaystyle \frac{1}{s}}V_p\right).`$ (105)
The solutions of which are

$`V_x(t)`$ $`=`$ $`{\displaystyle \frac{\mathrm{}}{2s}}{\displaystyle \frac{\left(2sV_x^0+\mathrm{}\mathrm{tanh}(2\gamma t)\right)\left(s\mathrm{}+2V_p^0\mathrm{tanh}(2\gamma t)\right)4sC_{xp}^{0}{}_{}{}^{2}\mathrm{tanh}(2\gamma t)}{\left(\mathrm{}+2sV_x^0\mathrm{tanh}(2\gamma t)\right)\left(s\mathrm{}+2V_p^0\mathrm{tanh}(2\gamma t)\right)4sC_{xp}^{0}{}_{}{}^{2}\mathrm{tanh}^2(2\gamma t)}}`$
$``$ $`{\displaystyle \frac{\mathrm{}}{2s}}\text{ as }t\mathrm{}`$
$`V_p(t)`$ $`=`$ $`{\displaystyle \frac{s\mathrm{}}{2}}{\displaystyle \frac{\left(2V_p^0+s\mathrm{}\mathrm{tanh}(2\gamma t)\right)\left(\mathrm{}+2sV_x^0\mathrm{tanh}(2\gamma t)\right)4sC_{xp}^{0}{}_{}{}^{2}\mathrm{tanh}(2\gamma t)}{\left(\mathrm{}+2sV_x^0\mathrm{tanh}(2\gamma t)\right)\left(s\mathrm{}+2V_p^0\mathrm{tanh}(2\gamma t)\right)4sC_{xp}^{0}{}_{}{}^{2}\mathrm{tanh}^2(2\gamma t)}}`$
$``$ $`{\displaystyle \frac{s\mathrm{}}{2}}\text{ as }t\mathrm{}`$
$`C_{xp}(t)`$ $`=`$ $`{\displaystyle \frac{s\mathrm{}^2C_{xp}^0\text{sech}^2(2\gamma t)}{\left(\mathrm{}+2sV_x^0\mathrm{tanh}(2\gamma t)\right)\left(s\mathrm{}+2V_p^0\mathrm{tanh}(2\gamma t)\right)4sC_{xp}^{0}{}_{}{}^{2}\mathrm{tanh}^2(2\gamma t)}}`$
$``$ $`0\text{ as }t\mathrm{}`$
where $`V_x(0)=V_x^0`$, $`V_p(0)=V_p^0`$ and $`C_{xp}(0)=C_{xp}^0`$. Hence the process of measurement induces the system state to collapse into a coherent state. If the measurement retrieves no information on momentum, i.e. $`\mathrm{\Gamma }_2=0`$, then

$`V_x(t)`$ $`=`$ $`{\displaystyle \frac{\mathrm{}V_x^0}{\mathrm{}+4V_x^0\mathrm{\Gamma }_1t}}`$ (109)
$`V_p(t)`$ $`=`$ $`V_p^0+\mathrm{}\mathrm{\Gamma }_1t{\displaystyle \frac{4C_{xp}^{0}{}_{}{}^{2}\mathrm{\Gamma }_1t}{\mathrm{}+4V_x^0\mathrm{\Gamma }_1t}}`$ (110)
$`C_{xp}(t)`$ $`=`$ $`{\displaystyle \frac{\mathrm{}C_{xp}^0}{\mathrm{}+4V_x^0\mathrm{\Gamma }_1t}}`$ (111)
and the growth in the momentum variance is unbounded. Similarly, if the measurement retrieves no information on position then the position variance grows unbounded. However, when both position and momentum are measured simultaneously the system state is forced into a coherent state. When $`\widehat{H}0`$ we expect that for a suitable choice of $`\gamma `$, any spreading of the quantum wavepacket caused by nonlinearities in the Hamiltonian will be counteracted by the measurement induced localization.

One might naively assume that if the measurement only retrieves information on position (or momentum) then the state will not localize. However this is not always the case. Note that when $`|C_{xp}^0|>\mathrm{}/2`$ in equation (110) the momentum variance will initially decrease. Thus if the system dynamics is such that it increases the covariance between position and momentum, then the continuous measurement of position may also localize momentum. For example, consider the Hamiltonian describing free particle motion, $`\widehat{H}=a\widehat{p}^2`$. When $`\mathrm{\Gamma }_2=0`$ the variances and covariance satisfy

$`{\displaystyle \frac{dV_x}{dt}}`$ $`=`$ $`{\displaystyle \frac{4\mathrm{\Gamma }_1}{\mathrm{}}}V_{x}^{}{}_{}{}^{2}+4aC_{xp}`$ (112)
$`{\displaystyle \frac{dV_p}{dt}}`$ $`=`$ $`\mathrm{\Gamma }_1\mathrm{}{\displaystyle \frac{4\mathrm{\Gamma }_1}{\mathrm{}}}C_{xp}^{}{}_{}{}^{2}`$ (113)
$`{\displaystyle \frac{dC_{xp}}{dt}}`$ $`=`$ $`{\displaystyle \frac{4\mathrm{\Gamma }_1}{\mathrm{}}}C_{xp}V_x+2aV_p.`$ (114)
Although we could not solve these equations analytically, it is easy to see that all physical solutions asymptotically attract to the stable fixed point

$$V_x=\sqrt{\frac{a}{2\mathrm{\Gamma }_1}}\mathrm{},V_p=\sqrt{\frac{\mathrm{\Gamma }_1}{2a}}\mathrm{},C_{xp}=\frac{\mathrm{}}{2}.$$
(115)
Hence, measurement of position does not introduce a diffusion in momentum, and the state localizes (this result has been derived previously in diosi2 ). However, for free particle motion, if we only measure momentum the state does not localize. The system dynamics accelerates the growth in the position variance. See gisin ; gisin2 ; schack ; brun ; percival ; strunz ; strunz2 ; doherty  for more on localization.

## V Numerical simulations

We will now numerically investigate the solution of the stochastic Schrödinger equation (90) for the driven system

$$H=ap^2+bx^2+cx^4+dx\mathrm{cos}(\omega t).$$
(116)
The numerical method to solve this equation is simple. To take advantage of the measurement induced localization we use a local moving number basis

$$|n\frac{1}{\sqrt{n!}}\widehat{a}^n|0n=0,1,\mathrm{}N(s=1).$$
(117)
truncated at some finite value $`N`$. The stochastic terms are integrated using the first-order Euler method while other terms are integrated by diagonalizing the position and momentum operators and using the split-operator formula.

We will first consider an integrable case when $`a=b=5`$, $`c=1`$, $`d=0`$ and $`\mathrm{}=0.05`$. The initial state was chosen to be a coherent state ($`s=1`$) centered at $`(x,p)=(2,1)`$ when $`t=0`$. The Husimi density of the initial state together with the contours of the Hamiltonian are plotted in Fig. 1(a). The Husimi density of the evolved state ($`t=4`$) together with the trajectories $`(\widehat{x},\widehat{p})`$ for different measurement schemes is plotted in figures 1(b)-(d). The evolved state in Fig. 1(b) is the result when no measurements occur ($`\gamma =0`$). In this case, nonlinearities in the Hamiltonian cause the state to shear as it evolves, spreading it along the contours. The trajectory has little meaning when $`\gamma =0`$. In Fig. 1(c) the evolved state is the result of continuous simultaneous measurement of position and momentum with $`\gamma =1/\sqrt{2}`$ and $`s=1`$. In this case the state has remained localized as it follows the contours. The continuous measurement of position only ($`\mathrm{\Gamma }_2=0`$) when $`\mathrm{\Gamma }_1=1`$ has also kept the state localized. This is shown in Fig. 1(d). The combined variance of position and momentum is plotted in Fig. 1(e). We must emphasize that only when both position and momentum are measured together does the trajectory correspond via (98,99) to the outcome of an actual measurement. If only position is measured, only $`\widehat{x}`$ observed while $`\widehat{p}`$ is simply the result of a mathematical calculation.

Now consider the chaotic case when $`a=5`$, $`b=8`$, $`c=1`$, $`d=15`$, $`\omega =2\pi `$ and $`\mathrm{}=0.05`$. A Poincaré stroboscopic map with unit strobing frequency is plotted in Fig. 2(a). For an initial state the same as above, the evolved state ($`t=5`$) together with the trajectories for different measurement schemes is plotted in figures 2(b)-(d). When no measurement occurs (Fig. 2(b)) the chaotic action of stretching and folding spreads the state across the phase space. However when we continuously measure position and momentum (Fig. 2(c)) the state remains localized and the trajectory resembles classical motion with noise. This noise will vanish as we approach the classical limit $`\mathrm{}0`$. In Fig. 3(a) we have plotted the quantum trajectory for the same parameter values as above except with $`\mathrm{}=10^6`$. In this case the total variance remains below $`10^5`$ (Fig. 3(b)) and the noise is not visible. The corresponding classical trajectory is plotted in grey and is only visible when it deviates from the quantum at $`t5`$. The evolution of this system under the continuous measurement of position has already been studied by Bhattacharya et al. habib . They also find that the measurement keeps the system state localized. It is not surprising that this is also the case when only momentum is measured (Fig. 2(d)).

## VI Conclusion

We have derived an Itô stochastic Schrödinger equation (90) describing the evolution of a quantum system under the continuous simultaneous measurement of position and momentum. The outcome of this measurement is a classical stochastic record obeying (77). As a consequence of the measurement, the system state is forced to remain localized allowing a classical interpretation of the quantum mean of the phase space variables as the trajectory of the system. This trajectory corresponds to the actual measured trajectory minus noise (98,99). Furthermore, the localization property allows a well-defined classical limit via Ehrenfest’s theorem. Indeed, for small $`\mathrm{}`$, numerical results show that the quantum system approximately follows classical trajectories. However a more complete theoretical understanding of the classical limit under continuous measurement is needed.

## VII Figure captions

FIG. 1. (a) Contours of $`H=5p^2+5x^2+x^4`$ and the initial state ($`\mathrm{}=0.05`$). (b-d) The trajectory and final state at $`t=4`$ when (b) $`\gamma =0`$, (c) $`\gamma =1/\sqrt{2}`$ and $`s=1`$, and, (d) $`\mathrm{\Gamma }_1=1`$ and $`\mathrm{\Gamma }_2=0`$. (e) The combined variance of position and momentum. All quantities are dimensionless.

FIG. 2. (a) Poincaré map for $`H=5p^28x^2+x^4+15x\mathrm{cos}(2\pi t)`$ and the initial state ($`\mathrm{}=0.05`$). (b-d) The trajectory and final state at $`t=5`$ when (b) $`\gamma =0`$, (c) $`\gamma =1/\sqrt{2}`$ and $`s=1`$, and, (d) $`\mathrm{\Gamma }_1=0`$ and $`\mathrm{\Gamma }_2=1`$. (e) The combined variance of position and momentum.

FIG. 3. (a) The quantum (black) and classical (grey) trajectories for the same Hamiltonian as in Fig. 2 except with $`\mathrm{}=10^6`$. (b)The combined variance of position and momentum.