text:
# A Submillimetre Selected Quasar in the Field of Abell478 Based on observations made with ESO Telescopes at the Paranal Observatory under programme IDs 63.O-0087 and 68.A-0111. Also based on observations with ISO, an ESA project with instruments funded by ESA member states (especially the PI countries: the Netherlands, the United Kingdom, Germany, and France) with the participation of ISAS and NASA.

## 1 Introduction

Major advances in submillimetre (submm) continuum observations came with the Submillimetre Common-User Bolometer Array (SCUBA, Holland et al. Holland (1999)), which is mounted at the 15 m James Clerk Maxwell Telescope (JCMT) at Mauna Kea, Hawaii. This dual-channel instrument for the first time allowed sensitive mapping, making it possible to survey larger areas of the sky to greater depths than previously possible at submm wavelengths. This development led to the discovery of a new class of objects of high infrared (IR) luminosity, located at cosmological distances (e.g., Smail et al. Smail (1997)). Even though these objects are less common than Lyman-break galaxies at similar redshifts, they would dominate the cosmic star formation rate density at these redshifts, if star formation is indeed the source of their high luminosity (e.g., Blain et al. Blainetal99 (1999)). Since the importance of these objects was realized, a number of submm surveys have been performed or are in progress (e.g. Smail et al. SIBK (2002) ; Eales et al. Eales (1999); Scott et al. scott (2002); Chapman et al. chapman\_lens (2002); Knudsen et al. in prep.). One of the biggest challenges for those surveys has turned out to be the follow-up observations and the identification of the counterparts causing the submm emission. As a result, the determination of the nature and redshift of these objects has been significantly hampered. While currently more than a hundred submm sources have been detected, less than 20 of these have reliably been identified with sources at other wavelengths, in spite of painstaking attempts. The number of published spectroscopic redshifts is even significantly smaller. The majority of the reliably identified counterparts are very or extremely red objects (e.g., Frayer et al. smmj00266 (2000); Smail et al. smail\_smmero (1999)); several are also exhibiting active galactic nuclei (AGN) features in their spectra (Ivison et al. smmj02399 (1999)). Because of the small number of secure identifications, any new identification adds important information to our understanding of the submm population. A survey with different selection criteria, radio-preselected and submm-detected, has produced 10 spectroscopic redshifts in a sample of 34 sources (Chapman et al. chapman\_nat (2003)).

We are carrying out an extensive SCUBA survey of a number of galaxy clusters fields, aimed at detecting gravitationally amplified background galaxies: the Leiden-SCUBA Lensed Survey (Knudsen et al. in prep.). In the course of doing the optical identifications and follow-up of this survey we discovered one of our submm sources to be a previously unknown type-1 quasar (previously reported in Knudsen et al. knudsen00 (2001)). While submm surveys of optically selected quasars have been quite succesful (Isaak et al. Isaak (2002)), this object is the first type-1 quasar first discovered by its submm emission. In contrast, type-2 quasars have been detected in small numbers in other submillimetre surveys (e.g. SMM J02399-0136, Vernet & Cimatti vernet (2002)), and IRAS-radio-optical quasars have been selected before at a wide range of redshifts (e.g. APM 08279+5255 in Irwin et al. irwin (1998)). In this paper we present the observations of the quasar. We discuss unusual properties of the object, its optical spectrum, and its IR spectral energy distribution, and compare the results to optically selected quasars. We adopt an $`\mathrm{\Omega }_0=0.3`$ and $`\mathrm{\Lambda }=0.7`$ cosmology with $`H_0=70\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$.

## 2 Observations and results

### 2.1 Submillimetre data

The SCUBA data of the $`z=0.088`$ galaxy cluster Abell 478 have been obtained during five nights in September and December 1997, March 1998 and December 1999. The first data were obtained in a program to study the cooling flow in the cluster itself. In these data a bright point source was detected. Consequently, extra data was obtained to study this object better. The total integration time was 6.6 hours (excluding overheads), recording data at both 850 $`\mu \mathrm{m}`$ and 450 $`\mu \mathrm{m}`$ simultaneously in jiggle-map mode. The data were obtained mostly under good conditions with 850 $`\mu \mathrm{m}`$ zenith atmospheric opacity typically around 0.2. The pointing was checked regularly and was found to be stable. Calibration maps of CRL618 were also obtained. The data were reduced using the SURF (SCUBA User Reduction Facility) and KAPPA software packages (Jenness & Lightfoot surf (1998)). The resulting images have an angular resolution of $`15^{\prime \prime }`$ at $`850\mu `$m and $`8^{\prime \prime }`$ at $`450\mu `$m.

Source extraction and estimation of the uncertainties were carried out using a method based on Mexican Hat wavelets (Cayón et al. mexhat (2000); Barnard et al. in prep; Knudsen et al. in prep), which was adopted for the entire Leiden-SCUBA Lensed Survey, and which will be described in a forthcoming publication (Knudsen et al. in prep.), where the full survey will be presented. This method was adopted because it is mathematically rigorous and its performance on SCUBA jiggle maps can be fully characterized. Monte Carlo simulations have been performed to determine the noise and uncertainties of the derived parameters. The area-weighted noiselevels of the maps are 2 mJy at $`850\mu \mathrm{m}`$ and 14 mJy at $`450\mu \mathrm{m}`$.

In the $`850\mu \mathrm{m}`$ map four sources were detected of which the brightest has a flux of $`S_{850}=25\pm 2.8\mathrm{mJy}`$. This is the only source in the map with detected $`450\mu \mathrm{m}`$ emission, $`S_{450}=55\pm 17\mathrm{mJy}`$. It was detected with a signal-to-noise of 15, for which the formal positional uncertainty including the pointing uncertainty of the JCMT is $`3.2^{\prime \prime }`$. This is the object SMM J04135+10277 for which we are here presenting the follow-up observations. Fluxes and positions are presented in Table 1.

### 2.2 Optical identification

For identification and redshift determination of the SCUBA source(s) optical imaging and spectroscopy was obtained with FORS1 at VLT-UT1 (Antu) in Chile, in September 1999. Four 15 min exposures in $`I`$band were aquired in photometric conditions. The frames were bias-subtracted, flatfielded and stacked. The resulting image is shown in Fig. 1 with the SCUBA $`850\mu \mathrm{m}`$ contours overlayed. The seeing measured in the final image is $`0\stackrel{}{.}9`$. The standard star field PG0231+051 (Landolt landolt (1992)) was used for the calibration. The source detection and photometry was performed using SExtractor (Bertin & Arnouts SEx (1996)). The center of SMM J04135+10277 is coincident with an $`I=20.5\pm 0.1\mathrm{mag}`$ point source at $`\alpha =04^\mathrm{h}13^\mathrm{m}27\stackrel{s}{.}28`$, $`\delta =10\mathrm{°}27\mathrm{}40\stackrel{}{.}4`$ (J2000). There are no other apparent candidate counterparts. The optical position is within the error circle of the submm observation. One of the other SCUBA sources (SMM J04134+10270) coincides with a galaxy, which, given its size and magnitude, is a probable cluster member. There are no obvious candidate counterparts for the two other SCUBA sources. Using the DIRBE/FIRAS maps (Schlegel et al. dustmaps (1998)), a Galactic reddening $`E(BV)0.52\mathrm{mag}`$ is derived — a substantial reddening. The corrected $`I`$ magnitude of the optical counterpart for SMM J04135+10277 is thus $`19.4\pm 0.1\mathrm{mag}`$.

### 2.3 Optical spectroscopy

FORS1 spectroscopy of SCUBA sources in the A478 field was also obtained in September 1999. We used FORS1 in Multi-Object Spectroscopy (MOS) mode to obtain spectra of a number of targets, using grism 150I+17, without order sorting filter. This setup gives a spectral resolution of 260 at 720 nm with the $`1^{\prime \prime }`$ slit which we employed. Overlap of the second spectral order may affect the wavelength region longwards of 650 nm, but was in the present case found not to affect the spectra. Two exposures of 1800 sec were obtained in a seeing of $`1.3^{\prime \prime }`$. The spectra were bias-subtracted and flatfielded. Wavelength calibration was achieved using exposures of He and Ar lamps. Correction for telluric absorption and flux calibration was carried out using observations of the white dwarf EG274 ($`V=11.03`$), which we corrected for photospheric absorption features. The multislit mask included both the bright SCUBA source SMM J04135+10277 and the fainter SCUBA source SMM J04134+10270. The extracted spectra were corrected for galactic foreground absorption using the DIRBE/FIRAS maps.

The optical spectrum of SMM J04135+10277 (Fig. 2) shows broad emission lines, of which the four most prominent can be identified with Ly$`\alpha `$+N v, Si iv+O iv\], C iv and C iii. In addition the spectrum shows a power-law continuum. All of these features are characteristic of quasars. Bluewards of the Ly $`\alpha `$ emission line, Ly $`\alpha `$ forest absorption is seen. We use the C iii, C iv and the Si iv+O iv\] lines, with the largest weight on the symmetric C iii line, to determine the redshift. We find the value $`z=2.837\pm 0.003`$, consistently for the peak value of all three profiles.

The spectrum of SMM J04134+10270 confirms its membership of the A478 cluster. It shows the characteristic spectrum of a quiescent elliptical galaxy with no evidence for nuclear activity. This source will be discussed together with the rest of the survey in Knudsen et al. (in prep.).

### 2.4 Near-infrared spectroscopy

Since restframe ultraviolet emission lines of quasars can be significantly blueshifted with respect to the systemic velocity (e.g., Carswell et al. Carswell (1991)), we also attempted to obtain additional redshift information using restframe optical lines. Unfortunately, at $`z2.84`$, the brightest lines (H$`\alpha `$, H$`\beta `$, \[O iii\] 5007 Å) are all in wavelength regions where the earth atmosphere is opaque. We therefore concentrated on the \[O ii\] 3727 Å line which is redshifted to the blue edge of the H-band window, a region strongly affected by atmospheric absorption lines.

We used ISAAC on VLT-UT1 (Antu) in February 2002, to take H-band spectra of SMM J04135+10277. We used the medium resolution grating with a $`1^{\prime \prime }`$ slit to obtain an $`R=3000`$ spectrum between 1.41 and $`1.49\mu `$m, which should contain the \[O ii\] line for redshifts between 2.78 and 3.00. In addition, we obtained a low-resolution spectrum ($`R=500`$) of the entire H-band, in an attempt to detect H$`\gamma `$, which although intrinsically faint, should at least lie in a clear part of the spectrum. Both spectra were obtained in photometric conditions and in an optical seeing of $`0.8^{\prime \prime }`$, by nodding the object along the slit. Total integration time was 44 minutes in both spectra. The individual frames were pairwise subtracted in order to remove the bright OH nightsky lines, flatfielded and coadded. Wavelength calibration was derived from the OH nightsky lines. Correction for telluric absorption and flux calibration were achieved using the B5V star Hip25499 (H=5.62) and the B2V star Hip28142 (H=7.497), corrected for photospheric absorption. Unfortunately, while the continuum of the quasar is clearly detected in both spectra, no emission features are seen. Undoubtedly, this is due to atmospheric absorptions in the region of the redshifted \[O ii\] line, and the faintness of the relevant features in the rest of the H-band spectrum.

### 2.5 ISO data

We also inspected the archive of the Infrared Space Observatory (ISO, Kessler et al. ISO (1996)) and extracted observations of A478 using both the mid-infrared camera (CAM, Cesarsky et al. CAM (1996)) and the mid/far-infrared photometer (PHOT, Lemke et al. PHOT (1996)).

The PHOT data were obtained using the P22 raster mode at $`180\mu `$m, with $`92^{\prime \prime }`$ pixels, on February 21, 1998. The data were reduced using the PHOT Interactive Analysis (Gabriel et al. PIA (1997)). Initial data reduction steps included discarding of corrupted data, non-linearity correction, and deglitching of individual ramps. After fitting all integration ramps with a first order polynomial, further deglitching and data editing, and dark current subtraction, the data were corrected for detector drifts and for vignetting, and calibrated using the internal Fine Calibration Sources. The resulting image shows a $`0.62\pm 0.19`$Jy source, the centroid of which is however displaced by $`46^{\prime \prime }`$ from SMM J04135+10277. This displacement is less than the ISO angular resolution at $`180\mu `$m, but much more than the nominal ISO pointing uncertainty of $`2^{\prime \prime }`$. While the centroid of a faint source can be displaced somewhat when placed on top of a highly structured background, in the present case the offset is so large that the $`180\mu `$m detection cannot reliably be associated with the quasar. Hence in the following we label this as an upper limit.

The CAM data were obtained using the LW3 filter (effective wavelength $`14.3\mu `$m) on February 21, 1998, and using the LW2 filter (effective wavelength $`6.7\mu `$m) on March 21, 1998 using the CAM01 raster observing mode with $`6^{\prime \prime }`$ pixels. The data were reduced using the CAM Interactive Analysis package (Ott et al. CIA (1997)). Processing steps consisted of dark current subtraction, deglitching and correction of transients using the PRETI method (Aussel et al. Aussel (1999)), which is particularly suited for the detection of faint sources, flatfielding using a flatfield derived from the stacked dataframes, and mosaic contruction taking into account the image distortion. This resulted in clear detections of SMM J04135+10277 at $`14.3\mu `$m with a flux density of $`0.47\pm 0.08`$Jy, and at $`6.7\mu `$m with a flux density of $`0.20\pm 0.03`$Jy. Positions of these sources are listed in Table 1.

### 2.6 CO $`J=32`$ emission

In a recent commisioning project of the new COBRA spectrometer on the Owens Vally Radio Observatory, the CO $`J=32`$ emission line has been detected from the quasar. This detection confirms the nature of SMM J04135+10277 as a hyperluminous IR quasar. The redshift implied by the CO line is $`z=2.84`$, which is indeed somewhat higher than the optically determined redshift. This result will be discussed in detail in Hainline et al. (in prep).

## 3 Discussion

### 3.1 SMM J04135+10277 and the optical quasar population

Only little is known about the importance of AGNs in the submm population. Most studies comparing X-ray and submm observations conclude that the submm population is powered by star formation rather than AGNs and especially quasars (e.g., Almaini et al. almaini (2003)), based on the lack of overlap of X-ray and submm sources in deep studies. This has been confirmed in a study combining very deep Chandra observations with SCUBA observations of the HDF-N, where Alexander et al. (alexander (2002)) found that a significant fraction of bright submm sources ($`f_{850\mu \mathrm{m}}>5\mathrm{mJy}`$) harbour an AGN, however, the AGN is not powerful enough to power the submm emission. This makes SMM J04135+10277 a particularly interesting object, since here we have a bright submm source that is unequivocally identified with a type-1 quasar. Given what is known about the abundance of type-1 quasars, is this an ordinary object that we should have expected to find in our survey, or are we dealing with an exceptional case? We here estimate the probability of finding a high redshift submm emitting quasar in our survey. The total area of our survey is $`65\mathrm{arcmin}^2`$ (Knudsen et al., in prep.). Using the optical spectrum we estimate that the quasar has a $`B`$ magnitude $`B21.021.5\mathrm{mag}`$. Based on the counts of Kennefick et al. (Kennefick (1997)), we find that there is only a 20% probability of finding a quasar with $`z>2.3`$ and $`16.5\mathrm{mag}<B<22\mathrm{mag}`$ in our survey. Furthermore, the probability that such a quasar is a bright submm source is also less than unity, as shown by Priddey et al. (Priddey (2002)), who did a submm study of optically selected quasars at $`1.5<z<3`$. For the sub-sample of quasars with $`z>2.3`$, only 30% of these had detectable submm emission down to 6.8 mJy and all of those are fainter than SMM J04135+10277. Combining the numbers we estimate only a 6% chance of detecting a submm bright quasar at $`z>2.3`$ in our survey, if that quasar was drawn from the population of optically selected quasars.

We also estimate the expected number of bright submm sources in the surveyed area, regardless of their physical nature. According to the number counts from Smail et al. (SIBK (2002)) we should expect to find two sources with $`850\mu `$m fluxes between 20 and 25 mJy. Our observations (Knudsen et al., in prep.) are in agreement with that number. Comparing this to the small chance of finding a high redshift submm emitting quasar in our survey, this result suggests that the bright part of the submm population does not originate from dusty quasars, and that SMM J04135+10277 is an unusual object.

### 3.2 Optical spectrum

Since SMM J04135+10277 is the first quasar selected based on its submm emission, it is of interest to compare its properties to those of optically selected quasars. Turning first to the optical spectrum, the shapes of the C iii, C iv and the Si iv+O iv\] lines appear as expected. However, the Ly $`\alpha `$+N v emission line has a more unexpected shape. The peak and blue wing appear to be absorbed. Furthermore, the strength of the line relative to the other emission lines is unusually low for a quasar. Since dust is present in this quasar, it is natural to assume that atomic hydrogen will also be present, so that associated absorption may play a role in suppressing the Ly $`\alpha `$ emission. However, for a more detailed assessment of this effect, a higher resolution spectrum is needed. Comparing the optical spectrum of SMM J04135+10277 to that of quasars selected at other wavelengths (see e.g. Francis et al. (Fran92 (1992)) for a composite spectrum), there are no significant differences except for the suppressed Ly $`\alpha `$ emission.

We note that SMM J04135+10277 belongs to the optically fainter part of the quasar population. Can this be the effect of strong absorption by dust, which would then simultaneously account for the luminous dust emission from SMM J04135+10277? We obtain a measure of the isotropic luminosity of the quasar nucleus using the C iv emission line; obviously, the Ly $`\alpha `$ line cannot be used since it appears to be absorbed, and the continuum cannot be used because of the effects of relativistic beaming, which cannot reliably be quantified. The observed flux of the C iv line is $`1.010^{14}`$erg s-1 cm-2. A comparison sample can be constructed from the optically selected submm emitting quasars studied by Priddey et al. (Priddey (2002)), using the spectra from Hagen et al. (Hagen (1999)). This comparison sample covers redshifts from 2.60 to 2.79 and can therefore be compared directly to SMM J04135+10277. The comparison sample has C iv fluxes from 1.7 to $`4.110^{14}`$erg s-1 cm-2, roughly a factor of 3 higher than SMM J04135+10277. The rest-frame equivalent width of C iv on the other hand shows the opposite trend: while SMM J04135+10277 has a C iv restframe equivalent width of $`170`$Å, values in the comparison sample are approximately a factor of 10 lower, ranging from 13 to $`25`$Å. In other words, the quasar continuum is fainter by about a factor of 30 than would be expected for its C iv flux. It is highly unlikely that extinction could account for this, since the quasar continuum and the broad line region should be viewed through approximately the same obscuring column. Furthermore, as Fig. 2 shows, the quasar continuum is characterized by a blue power law. The slope of this continuum does not indicate the presence of abnormally large absorption. Therefore a more likely explanation of the optical faintness of this quasar is a large viewing angle away from the direction of relativistic beaming. The beamed flux density is proportional to $`\delta ^p`$ with $`p4`$, where the Doppler factor $`\delta =[\gamma (1\beta \mathrm{cos}\theta )]^1`$, where $`\beta `$ is the bulk velocity in units of the speed of light, and $`\gamma =(1\beta )^{1/2}`$ is the corresponding Lorentz factor, and $`\theta `$ is the angle away from the beam (Urry & Padovani UrryPadovani (1995)). Therefore a decrease in $`\delta `$ of a factor 2.3 would be sufficient to produce a factor 30 decrease in the beamed continuum with respect to the lines. The required angle away from the beam cannot be calculated since $`\beta `$ is not known. However, as shown by Urry & Padovani (UrryPadovani (1995)), variations in $`\delta `$ of this magnitude are entirely reasonable for angles $`\theta <20\mathrm{°}`$, provided $`\gamma >2`$. This estimate confirms the viability of our suggestion that the optical faintness of the quasar is due to a large viewing angle away from the direction of relativistic beaming, and not to abnormally large extinction. If in fact the optical spectrum is still dominated by the doppler boosted jet then our detection of this one object suggests that a much larger number of yet unidentified sources are similar AGNs viewed from a larger angle to the jet axis. We finally note that it would be interesting to make the same comparison with low-$`z`$ far-IR detected quasars, addressing also the properties of the dust emssion spectrum. This comparison would require spectrophotometry of quasars in the vacuum ultraviolet.

Going further, we can investigate whether the observed submm emission from SMM J04135+10277 is likely powered by the AGN or whether the presence of an additional power source is indicated. In the comparison sample, the observed $`850\mu `$m fluxes range from 6.8 to 10.0 mJy, increasing monotonically with C iv flux. The three times fainter C iv flux of SMM 04135+10277 thus would suggest an AGN-powered $`850\mu `$m flux of approximately 3 mJy. The observed flux is almost a factor of 10 higher. This result suggests that the bulk of the submm emission from SMM 04135+10277 is not powered by the AGN but by an additional source of energy, most likely vigorous star formation in the host galaxy. If this interpretation is correct, high resolution imaging of CO lines and dust emission with ALMA should reveal an extended source.

### 3.3 Gravitational magnification

The low redshift ($`z=0.088`$) of A478 is far from the optimal lensing redshift ($`z0.2`$), and no arcs are detected in the vicinity of the quasar. This suggests that the gravitational magnification is small and that the quasar is not subject to differential lensing, which otherwise would influence the shape of the spectral energy distribution (SED). We calculate the magnification of the quasar using LENSTOOL (Kneib et al. lenstool (1993)). The cluster is modelled using two components: the overall cluster potential with a core radius $`r_\mathrm{c}=250\mathrm{kpc}`$ and a velocity dispersion $`\sigma =905\mathrm{km}\mathrm{s}^1`$ and the potential of the cD galaxy with $`r_\mathrm{c}=2\mathrm{kpc}`$ and $`\sigma =350\mathrm{km}\mathrm{s}^1`$ (Allen et al. A478numbers1 (1993); Zabludoff et al. A478numbers2 (1990)). We find that the quasar is magnified by a factor of 1.3. Hence, all fluxes should be corrected for this value. In all calculations in the following sections of this paper the fluxes have been corrected for the gravitational lensing.

### 3.4 Spectral energy distribution

From archival data from the NRAO Very Large Array a radio source near the position of the SMM J04135+10277 was found. The fluxes measured are $`S(4.86\mathrm{GHz})=220\pm 35\mu \mathrm{Jy}`$ and $`S(1.4\mathrm{GHz})<750\mu \mathrm{Jy}`$ ($`3\sigma `$) (M. Yun, private comm.). Based on this low radio flux, it is concluded that the quasar is radio-quiet (according to the radio-power criterion given by Stocke et al. stocke (1992) to divide quasars into radio-loud and radio-quiet types). Hence, the non-thermal contribution to the submm flux is expected to be small and is here neglected.

Combined, we have then six points on the SED and two upper limits, ranging from the radio to the optical regime (see Table 1). The SCUBA points and the two ISO points at $`14.3\mu \mathrm{m}`$ and $`170\mu \mathrm{m}`$ are in the rest frame all at wavelengths typical for thermal emission by dust. The ISO $`6.7\mu \mathrm{m}`$ point in rest frame is $`1.7\mu \mathrm{m}`$, which together with the optical point is expected to originate from stellar light, possibly contaminated with non-thermal emission from the AGN. The radio emission is attributed to non-thermal synchrotron emission.

In the following we will focus on the thermal dust emission. Of the thermal emission, the two SCUBA points are most likely due to the cool dust typically described by a modified blackbody, whereas the $`14.3\mu \mathrm{m}`$ point arises from a hot component. The shape of the SED between the cold and hot component seems to be different for different quasars (see e.g. Haas et al. haas\_pg (2000)). We have no measurements between $`450\mu \mathrm{m}`$ and $`14.3\mu \mathrm{m}`$, which makes an analysis of the IR SED difficult. We do, however, attempt to make a tentative analysis in which we compare with known objects and also estimate parameters like temperature and luminosity.

We first compare the quasar with other known quasars. Comparing to high-$`z`$ quasars is not trivial, since the high-$`z`$ quasars which have well-sampled IR SEDs, are often strongly lensed and their observed SEDs may have suffered differential lensing. We therefore first focus on low-$`z`$ quasars. Haas et al. (haas\_pg (2000)) have made a detailed study of the IR SED of Palomar-Green (PG) quasars. The majority of these quasars are at fairly low redshift. We compare to three low-$`z`$ PG quasars with well-sampled SEDs, PG 0050+124, PG 1206+439 and PG1613+638 (all shown in Fig. 3). All three SEDs are redshifted to $`z=2.837`$. If the SEDs are scaled to the quasar $`850\mu \mathrm{m}`$ point, the comparison gives the impression of a deficit in the mid/near-IR emission of SMM J04135+10277. Alternatively, inspired by the findings of Archibald et al. (2001) and Page et al. (2001), that the star formation rate observed in AGNs is higher at higher redshift, leading to enhanced long-wavelength emission at higher redshift, we may choose instead to scale the low-$`z`$ SEDs to the observed $`14.3\mu \mathrm{m}`$ point, i.e., the hot dust emission associated with the AGN. This, not unexpectedly, then suggests an excess in the far-IR-submm emission from SMM J04135+10277. This result corroborates our earlier conclusion that a significant portion of the observed $`850\mu `$m emission of SMM J04135+10277 results from extended star formation, and is not powered directly by the AGN. The SED of the strongly lensed $`z=3.87`$ quasar APM 08279+5255 (Lewis et al. lewis (1998) and references therein) is also shown in Fig. 3. It has also been appropriately shifted and scaled to the $`14.3\mu \mathrm{m}`$ point. In this case the submm/FIR deficit relative to the SMM J04135+10277 is even more pronounced, corroborating the discussion above.

Given this result, it is also of interest to compare the SED of SMM J04135+10277 to the SEDs of well-studied starburst galaxies. We use the SEDs of the starburst galaxy NGC 253 and the ultraluminous infrared galaxy (ULIRG) NGC 6240 (extracted from the NASA Extragalactic Database), redshifted to $`z=2.837`$, for comparison. Scaled to the $`850\mu \mathrm{m}`$ point, the far-IR/submm range matches quite well, whereas the mid/near-IR emission is much brighter for SMM J04135+10277. This result is expected, as starbursts are known not to have the hot dust component that is characteristic of AGNs, especially quasars (Sanders et al. Sanders (1989); Barvainis Barvainis (1987)). Considering that the SED of the quasar is not well-sampled, and that we are looking at only one quasar, no definite conclusions can be drawn about the precise shape of the IR SED and the power source of the dust emission of SMM 04135+10277 in particular, or of the submm-selected quasar population in general. Tentatively, though, these observations support the suggestion that quasars at higher redshift have a high submm and far-IR flux, suggesting a higher star formation rate. To study this in detail, however, observations filling the big gaps in the IR SED are needed. Such data can possibly be obtained with SIRTF.

Finally, we calculate the dust temperature, dust mass, and total luminosity in dust emission of SMM J04135+10277. For the dust emission, we use blackbody emission, modified by the frequency-dependent mass absorption coefficient

$$k_\mathrm{d}(\nu _{\mathrm{rest}})=1.5\mathrm{cm}^2\mathrm{g}^1\left(\frac{\nu _{\mathrm{rest}}}{375\mathrm{GHz}}\right)^\beta ,$$
(1)
using the average value from the literature at $`800\mu `$m (see Hughes et al. 1997 for a discussion of the assumptions and uncertainties in this parameter) and assume $`\beta =1.5`$. Since the IR SED cannot be fit by a single modified blackbody, we first fit the cold dust component sampled by the SCUBA points. Using only the SCUBA points, we find a temperature $`T_{\mathrm{cold}}=29\pm 2\mathrm{K}`$ and a dust mass of $`M_\mathrm{d}=(1.8\pm 0.3)\times 10^9\mathrm{M}_{}`$ for the cold dust component. The total luminosity of a modified blackbody spectrum can be calculated using the following analytical expression:

$$L_{\mathrm{tot}}=4\pi D_\mathrm{L}^2\frac{S(\nu _{\mathrm{obs}})}{\nu _0^\beta B(\nu _0,T)}\frac{2h}{c^2}\left(\frac{kT}{h}\right)^{4+\beta }\mathrm{\Gamma }(4+\beta )\zeta (4+\beta ),$$
(2)
where $`h`$ and $`k`$ are Planck’s respectively Boltzmann’s constants, $`D_\mathrm{L}`$ is the luminosity distance, and the two last factors are the Gamma function and Riemann’s zeta function. $`S(\nu _{\mathrm{obs}})`$ is the flux density at the observing frequency, and $`\nu _0`$ is the corresponding rest frequency. The total luminosity of the cold component is found to be $`(2.4\pm 0.5)10^{13}\mathrm{L}_{}`$. To estimate the total IR luminosity, we assume a powerlaw between the peak of the modified blackbody curve and the observed $`14.3\mu \mathrm{m}`$ point. The powerlaw is integrated from $`200\mu \mathrm{m}`$ (observed frame), where the modified blackbody and the powerlaw balance eachother, to $`14.3\mu \mathrm{m}`$ with the result of $`(5.4\pm 1)10^{12}\mathrm{L}_{}`$. In total the IR luminosity (corrected for gravitational amplification) is then $`L_{\mathrm{IR}}=(2.9\pm 0.5)10^{13}\mathrm{L}_{}`$, dominated by the cold dust component. Using a different method based on the analysis in Blain et al. (blain03 (2003)), where the whole IR SED is fit with a single temperature modified blackbody with a powerlaw on the Wien side ranging all the way into the mid-IR, a temperature of 38 K is found and a total IR luminosity of $`1.8\times 10^{13}\mathrm{L}_{}`$. This gives a higher temperature, though a slightly lower luminosity, compared to the fit above where a cold component was fitted to the two SCUBA points.

The temperature as we find is lower than that found in other high-$`z`$ quasars such as APM 08279+5255, which has a temperature of $`120220`$K determined for a pure blackbody (Lewis et al. lewis (1998)), or BR 1202$``$0725, which has a dust temperature of $`5068`$K (Leech et al. leech (2001)). Both quasars have luminosities in order of $`10^{1415}\mathrm{L}_{}`$, thus brighter than SMM J04135+10277, so that higher dust temperatures might be expected. On the other hand, these two quasars are strongly lensed and it is possible that differential lensing distorts the integrated SED and overemphasizes warm dust components.

We finally attempt to compare the radio-submm flux density ratio with the relevant simulations performed by Blain (blain (1999)), which are based on the IR-radio correlation observed at low redshift. As we do not have a 1.4 GHz flux density measurement, we estimate it by assuming that the radio SED is a power law, $`f\nu ^\alpha `$, with slope $`\alpha =0.8`$ and scale it to the observed flux density at 4.86 GHz. We find $`f_{1.4\mathrm{GHz}}=595\mu \mathrm{Jy}`$ (not corrected for the gravitational lensing). Still assuming $`\beta =1.5`$, we use Fig. 4 in Blain (blain (1999)) by interpolating between his two models with $`T=20\mathrm{K}`$ respectively $`T=40\mathrm{K}`$. For $`z=2.837`$ this gives a flux density ratio of between $`1.4`$GHz and $`850\mu `$m of $`3.510^3`$. Therefore, the observed $`850\mu `$m flux would imply $`f_{1.4\mathrm{GHz}}88\mu \mathrm{Jy}`$ if the quasar strictly followed the local IR-radio correlation. This number is however a factor $`67`$ lower than what we had just estimated based above. This indicates that SMM J04135+10277 has more radio emission (for its IR emission) than e.g., the ULIRG Arp220, which was used for the template SED in Blain (1999). This result is not surprising, as quasar radio emission is powered by both the synschrotron emission from stellar remnants and the synchrotron emission from the central black hole.

## 4 Conclusions

We have discovered a type-1 quasar behind the cluster of galaxies A478. The quasar, SMM J04135+10277, was discovered by it submm emission. The quasar has a redshift of $`z=2.837`$ and is radio quiet. The quasar is optically faint, but has a large submm flux. Using number counts of quasars and of submm sources, we argue that SMM 04135+10277 is an unusual object. It is in any case a remarkable object since there is little overlap between deep submm and X-ray samples, suggesting that the incidence of powerful AGNs among submm galaxies is low. The slope of the rest-frame UV continuum is similar to that of optically selected quasars (Francis et al. Fran92 (1992)) and does not exhibit any signs of extraordinary dust extinction. This leads us to suggest that the line of sight to the quasar nucleus is not abnormally obscured. Analysis of optical continuum, spectral lines and submm emission leads us to conclude that the optical faintness of the quasar results from a large viewing angle from the direction of relativistic beaming, and that a significant amount of the submm flux is not powered by the active nucleus. More likely the cold dust is heated by a high rate of star formation in the environment surrounding the quasar. Comparison of the sparsely sampled IR SED to that of other objects tentatively supports this conclusion. The total IR luminosity is found to be $`(2.9\pm 0.5)10^{13}\mathrm{L}_{}`$ and is dominated by the emission from cool dust.

###### Acknowledgements.

We thank Remo Tilanus for taking most of the SCUBA data presented in this paper, Min Yun for providing us with the VLA archive data, and Jean-Paul Kneib for making his LENSTOOL program available for us. We also thank the referee, Andrew Blain, for useful comments. KKK is supported by the Netherlands Organization for Scientific Research (NWO). The JCMT is operated by the Joint Astronomy Centre on behalf of the United Kingdom Particle Physics and Astronomy Research Council (PPARC), the Netherlands Organization for Scientific Research, and the National Research Council of Canada. The National Radio Astronomy Observatory is a facility of the National Science Foundation, operated under cooperative agreement by Associated Universities, Inc. This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. The ISOCAM data presented in this paper was analysed using “CIA”, a joint development by the ESA Astrophysics Division and the ISOCAM Consortium. The ISOCAM Consortium is led by the ISOCAM PI, C. Cesarsky. The ISOPHOT data presented in this paper was reduced using PIA, which is a joint development by the ESA Astrophysics Division and the ISOPHOT consortium, with the collaboration of the Infrared Analysis and Processing Center (IPAC) and the Instituto de Astrof sica de Canarias (IAC).


format: text

text:
# Differential Astrometry of Sub-arcsecond Scale Binaries at the Palomar Testbed Interferometer

## 1 Introduction

Long-baseline optical interferometry promises high precision astrometry using modest ground-based instruments. In particular the Mark III Stellar Interferometer (Shao et al., 1988) and Navy Prototype Optical Interferometer (NPOI, Armstrong et al. 1998) have achieved global astrometric precision at the 10 mas (1 mas = $`10^3`$ arcseconds) level (Hummel et al., 1994), while the Palomar Testbed Interferometer (PTI, Colavita et al. 1999 ) has demonstrated an astrometric precision of 100 $`\mu `$as ($`1\mu `$as = $`10^6`$ arcseconds) between moderately close (30 arcsecond) pairs of bright stars (Shao & Colavita, 1992; Colavita, 1994; Lane et al., 2000). While interferometric and astrometric methods have proven very useful in studying binary stars, and have long been argued to be well-suited to studying extra-solar planets (Colavita & Shao, 1994; Eisner & Kulkarni, 2001), to date results using these techniques have been limited (Benedict et al., 2002b).

There are several reasons why it is desirable to develop viable astrometric planet-detection methods. Most importantly, the parameter space explored by astrometry is complementary to that of radial velocity (astrometry is more sensitive to larger separations). Second, unlike current radial velocity detections, astrometric techniques can be used to determine the orbital inclination of a planet. Finally, astrometry is particularly well-suited to studying binary stellar systems; such systems challenge other planet-finding techniques. For example, radial velocimetry can suffer from systematic velocity errors caused by spectral contamination from the light of the second star (Vogt et al., 2000). Similar problems are faced by coronographic techniques, where the light from the second star is not usually blocked by the occulting mask.

In this paper we describe recent efforts to obtain very high precision narrow-angle astrometry using PTI to observe binary stars with separations less than one arcsecond, i.e. systems that are typically observed using speckle interferometry (Saha, 2002) or adaptive optics. Such small separations allow us to achieve astrometric precision on the order of 10 $`\mu `$as, which for a typical binary system in our target sample (binary separation of 20 AU), should allow us to detect planets with masses down to 0.5 Jupiter masses in orbits in the 2 AU range. This approach has been suggested (Traub, Carleton & Porro, 1996) and tried (Dyck, Benson & Schloerb, 1995; Bagnuolo et al., 2003) before, though with limited precision. However, this work is unique in that it makes use of a phase-tracking interferometer; the use of phase-referencing (Lane & Colavita, 2003) removes much of the effect of atmospheric turbulence, improving the astrometric precision by a factor of order 100.

The Palomar Testbed Interferometer (PTI) is located on Palomar Mountain near San Diego, CA (Colavita et al., 1999). It was developed by the Jet Propulsion Laboratory, California Institute of Technology for NASA, as a testbed for interferometric techniques applicable to the Keck Interferometer and other missions such as the Space Interferometry Mission, SIM. It operates in the J ($`1.2\mu \mathrm{m}`$),H ($`1.6\mu \mathrm{m}`$) and K ($`2.2\mu \mathrm{m}`$) bands, and combines starlight from two out of three available 40-cm apertures. The apertures form a triangle with 86 and 110 meter baselines.

The paper is organized as follows: in Section 2 we describe the experiment and derive expected performance levels. In Section 3 we describe initial observations as well as the extensive data analysis processing required to achieve the desired astrometric precision. In Section 4 we discuss our preliminary results, and in Section 5 we discuss the prospects of a larger search.

## 2 Interferometric Astrometry

In an optical interferometer light is collected at two or more apertures and brought to a central location where the beams are combined and a fringe pattern produced. For a broadband source of wavelength $`\lambda `$ the fringe pattern is limited in extent and appears only when the optical paths through the arms of the interferometer are equalized to within a coherence length ($`\mathrm{\Lambda }=\lambda ^2/\mathrm{\Delta }\lambda `$). For a two-aperture interferometer, neglecting dispersion, the intensity measured at one of the combined beams is given by

$$I(\delta )=I_0\left(1+V\frac{\mathrm{sin}\left(\pi \delta /\mathrm{\Lambda }\right)}{\pi \delta /\mathrm{\Lambda }}\mathrm{sin}\left(2\pi \delta /\lambda \right)\right)$$
(1)
where $`V`$ is the fringe contrast or “visibility”, which can be related to the morphology of the source, and $`\mathrm{\Delta }\lambda `$ is the optical bandwidth of the interferometer assuming a flat optical bandpass (for PTI $`\mathrm{\Delta }\lambda =0.4\mu `$m). The differential (between arms of the interferometer) optical path $`\delta `$ is found from geometric considerations to be

$$\delta =\stackrel{}{B}\stackrel{}{s}d$$
(2)
$`\stackrel{}{B}`$ is the baseline – the vector connecting the two apertures. $`\stackrel{}{s}`$ is the unit vector in the source direction, and $`d`$ is any differential optical path introduced by the instrument. For a 100-m baseline interferometer an astrometric precision of 10 $`\mu `$as corresponds to knowing $`d`$ to 5 nm, a difficult but not impossible proposition.

The dominant source of measurement error is atmospheric turbulence above the interferometer, which adds varying amounts of optical path and hence makes the fringes appear to move about rapidly (Roddier, 1981). This atmospheric turbulence, which changes over distances of tens of centimeters and millisecond time-scales, forces one to use very short exposures to maintain fringe contrast, and hence limits the sensitivity of the instrument. It also severely limits the astrometric accuracy of a simple interferometer, at least over large sky angles. However, this atmospheric turbulence is correlated over small angles, and hence it is still possible to obtain high precision “narrow-angle” astrometry.

### 2.1 Narrow-Angle Astrometry

Dual-star interferometric narrow-angle astrometry (Shao & Colavita, 1992; Colavita, 1994) promises astrometric performance at the 10-100 $`\mu `$as level for pairs of stars separated by 10-60 arcseconds, and has been demonstrated at PTI. Achieving such performance requires simultaneous measurements of the fringe positions of both stars via the use of a “dual-star” optical beam-train.

For more closely spaced stars, it is possible to operate in a simpler mode. We have recently used PTI to observe pairs of stars separated by no more than one arcsecond. In this mode, the small separation of the binary results in both binary components being in the field of view of a single interferometric beam combiner. The fringe positions are measured by modulating the instrumental delay with an amplitude large enough to record both fringe packets. This eliminates the need for a complex internal metrology system to measure the entire optical path of the interferometer, and dramatically reduces the effect of systematic error sources such as uncertainty in the baseline vector (error sources which scale with the binary separation).

However, since the fringe position measurement of the two stars is no longer truly simultaneous it is possible for the atmosphere to introduce pathlength changes (and hence positional error) in the time between measurements of the separate fringes. To reduce this effect we split off a fraction of the incoming starlight and direct it to a second beam-combiner. This beam-combiner is used in a “fringe-tracking” mode (Shao & Staelin, 1980; Colavita et al., 1999) where it rapidly (10 ms) measures the phase of one of the starlight fringes, and adjusts the internal delay to keep that phase constant. This technique – known as phase referencing – has the effect of stabilizing the fringe measured by the astrometric beam-combiner. In addition, the residual phase error measured by the fringe tracker is a high time-resolution trace of the phase error introduced by the atmosphere and not fully corrected by the fringe tracker; it can be applied to the measured fringe position in post-processing.

### 2.2 Expected Performance

In making an astrometric measurement we modulate the optical delay applied internally in a triangle-wave pattern around the stabilized fringe position, while measuring the intensity of the combined starlight beams. Typically we obtain one such “scan” every 1-3 seconds, consisting of up to 3000 intensity samples. The range of the delay sweep is set to include both fringe packets; typically this requires a scan amplitude on the order of $`100\mu `$m. We then fit a double fringe packet based on Eqn. 1 to the data, and solve for the differential optical path between fringe packets.

In calculating the expected astrometric performance we take into account three major sources of error: errors caused by fringe motion during the sweep between fringes (loss of coherence with time), errors caused by differential atmospheric turbulence (loss of coherence with sky angle, i.e. anisoplanatism), and measurement noise in the fringe position. We quantify each in turn below, and the expected measurement precision is the root-sum-squared of the terms (Figure 2).

#### 2.2.1 Loss of Temporal Coherence

The power spectral density of the fringe phase of a source observed through the atmosphere has a power-law dependence on frequency (Figure 1); at high frequencies typically

$$A(f)f^\alpha $$
(3)
where $`\alpha `$ is usually in the range 2.5–2.7. The effect of phase-referencing is to high-pass filter this atmospheric phase noise. In our case, the servo is an integrating servo with finite processing delays and integration times, with the residual phase error “fed forward” to the second beam combiner (Lane & Colavita, 2003). The response of this system to an input atmospheric noise can be written in terms of frequency (see Appendix A in Lane 2003) as

$`H(f)`$ $`=`$ $`{\displaystyle \frac{12\mathrm{s}\mathrm{i}\mathrm{n}\mathrm{c}(\pi fT_s)\mathrm{cos}(2\pi fT_d)+\mathrm{sinc}^2(\pi fT_s)}{12\frac{f_c}{f}\mathrm{sinc}(\pi fT_s)\mathrm{sin}(2\pi fT_d)+\left(\frac{f_c}{f}\right)^2\mathrm{sinc}^2(\pi fT_s)}}`$ (4)
where $`\mathrm{sinc}(x)=\mathrm{sin}(x)/x`$, $`f_c`$ is the closed-loop bandwidth of the fringe-tracker servo (for this experiment $`f_c=10`$ Hz), $`T_s`$ is the integration time of the phase sample (6.75 ms), and $`T_d`$ is the delay between measurement and correction (done in post-processing, effectively 5 ms). The phase noise superimposed on the double fringe measured by the astrometric beam combiner has a spectrum given by $`A(f)H(f)`$.

The sampling of the double fringe packet takes a finite amount of time, first sampling one fringe, then the other. In the time domain the sampling function can be represented as a “top-hat” function convolved with a pair of delta functions (one positive, one negative). The width of the top-hat is equal to the time taken to sweep through a single fringe, while the separation between the delta functions is equal to the time to sweep between fringes. In the frequency domain this sampling function becomes

$$S(f)=\mathrm{sin}^2(2\pi f\tau _p)\mathrm{sinc}^2(\pi f\tau _{})$$
(5)
where $`\tau _p`$ is the time taken to move the delay between stars ,$`\mathrm{\Delta }d/v_s`$, and $`\tau _{}`$ is the time to sweep through a single stellar fringe, $`\mathrm{\Lambda }/v_s`$. $`v_s`$ is the delay sweep rate.

The resulting error in the astrometric measurement, given in radians by $`\sigma _{tc}`$, can be found from

$$\sigma _{tc}^2=\left(\frac{\lambda }{2\pi \mathrm{B}}\right)^2\frac{1}{N}_0^{\mathrm{}}A(f)H(f)S(f)𝑑f$$
(6)
where $`N`$ is the number of measurements. It is worth noting that if phase-referencing is not used to stabilize the fringe, i.e. $`H(f)=1`$, the atmospheric noise contribution increases by a factor of $`10^2`$$`10^3`$.

#### 2.2.2 Anisoplanatism

The performance of a simultaneous narrow-angle astrometric measurement has been thoroughly analyzed in Shao & Colavita (1992). Here we restate the primary result for the case of typical seeing at a site such as Palomar Mountain, where the astrometric error in arcseconds due to anisoplanatism ($`\sigma _a`$) is given by

$$\sigma _a=540\mathrm{B}^{2/3}\theta t^{1/2}$$
(7)
where B is the baseline, $`\theta `$ is the angular separation of the stars, and $`t`$ the integration time in seconds. This assumes a standard (Lindegren, 1980) atmospheric turbulence profile; it is likely that particularly good sites will have somewhat (factor of two) better performance.

#### 2.2.3 Photon Noise

The astrometric error due to photon-noise ($`\sigma _p`$) is given in radians as

$$\sigma _p=\frac{\lambda }{2\pi \mathrm{B}}\frac{1}{\sqrt{N}}\frac{1}{\mathrm{SNR}}$$
(8)
where N is the number of fringe scans, and SNR is the signal-to-noise ratio of an individual fringe.

## 3 Observations & Data Processing

We used the Palomar Testbed Interferometer to observe the binary star HD 171779 (HR 6983, K0III+G9III, $`m_K=2.78`$, $`m_V=5.37`$) on the nights of 11-14 August 2003. Conditions were calm with good seeing although with high-altitude clouds limiting observations, particularly after the first night. Based on a visual orbit obtained from speckle interferometry (Hartkopf et al., 2001), the predicted separation was 0.248 arcsec at a position-angle of 120.9 degrees. The orbit quality was listed as 3 (on a scale of 1–5), indicating a “reliable” orbit, though no uncertainties are given. The published $`\mathrm{\Delta }V=0.21`$ (Worley, Mason & Wycoff, 2001), and the orbital period is 191.49 years.

We used the longest available baseline (110 m). The observations were done in the K band ($`2.2\mu `$m), with $``$70% of the K-band light being used for fringe tracking, and the remaining light going to the astrometric measurement. The fringe tracker operated in the standard configuration (Colavita et al., 1999), with a sample time of 10 ms and a closed-loop servo bandwidth of $`10`$ Hz. The delay modulation used for the astrometric measurement typically had an amplitude of 150 $`\mu `$m and period of 3 seconds. The modulation was done using a PZT-actuated mirror, the position of which was measured using a laser metrology system.

The intensity vs. delay position measurements produced by the interferometer were processed into astrometric measurements as follows: (1) Detector calibrations (gain, bias and background) were applied to the intensity measurements. (2) The residual phase errors from the primary fringe tracker were converted to delay and applied to the data. Note that while the intensity measurements were spaced regularly in time, and the delay scanned linearly in time, the variable amount of delay correction applied from the fringe tracker resulted in the intensity measurements being unevenly spaced in delay. This somewhat complicated the downstream processing, in that FFT-based algorithms could not be used. (3) The data were broken up into “scans” either when the delay sweep changed direction or when the fringe tracker lost lock. (4) For each scan, a power spectrum was calculated using a Lomb-Scargle (Scargle, 1982; Press et al., 1992) algorithm. This spectrum provided an SNR estimate based on the ratio of the power in and out of the instrument bandpass. Only the scans with an SNR greater than unity were kept. (5) The intensity measurements were optionally bandpass-filtered to remove the effects of atmospheric turbulence changing the amount of light being coupled into the detector. We note however that the final results of the fit (see below) did not depend on whether or not the filter was applied; the results differed by $`<<1\sigma `$ (the results shown were based on unfiltered data). See Figure 3 for an example of filtered fringes. (6) The differential delay of each scan was found by evaluating the least-squares goodness-of-fit ($`\chi ^2`$) parameter for a fit of a double fringe packet to each scan, for the range of plausible delay separations and mean delay values, and selecting the separation corresponding to the $`\chi ^2`$-minimum. See the discussion below for how plausible ranges were determined. (7) Given a set of delay separations as a function of time, a simple astrometric model was fit to the data. The free parameters were separation and position angle, or equivalently, differential right ascension and declination (Figure 4).

In order to characterize any systematic dependence on assumptions made in the data reduction, such as the shape of the fringe envelope, the effect of dispersion, the effect of filtering and the level of the fringe SNR-cutoff, we exhaustively varied these parameters and and re-fit the data each time. In all cases the effect was small compared to the claimed uncertainties. This suggests that our results are robust against systematic errors in the data processing algorithm; many potential systematic error sources are eliminated due to the differential nature of the measurement.

### 3.1 Resolving Fringe Ambiguity

The oscillatory nature of the fringe means that there will be many local minima in the $`\chi ^2`$ surface, separated from the global minium (at the center of the fringe packet) by integer multiples of the wavelength $`\lambda `$. In the presence of noise it is not always obvious which is the global minimum. In previous work (Dyck, Benson & Schloerb, 1995) this ambiguity was avoided by forming an “envelope” of the fringe packets, then fitting these envelopes. However, this effectively removed the high-resolution phase information, which is what provides the high precision astrometry.

We find that it is possible to co-add the $`\chi ^2(\mathrm{\Delta }d)`$ functions from many scans. This increases the signal-to-noise ratio and so reduces the likelihood of picking the “wrong” fringe. In order to perform the co-adding, the two-dimensional (separation and mean delay) $`\chi ^2`$ surface was projected into a one-dimensional space of delay separation by selecting the best (lowest $`\chi ^2`$) mean delay position for each possible delay separation. The minimum of this co-added $`\chi ^2`$ function was used to determine the fringe separation to within $`\lambda `$; the range of plausible separations (used in step 6 above) was limited to $`\pm \lambda /2`$ around this value.

The number of scans that could be co-added in this manner depended on the rate of change of delay separation, with too long of an integration time “smearing” the fringes. We found that in the case of our observations, co-adding all the $`\chi ^2`$ functions from a 75-second period gave adequate SNR to reliably (85–95% of the time) determine the central fringe. The remaining 5–15% of scans produced delays that were shifted off the median fringe position by a wavelength. However, rather than artificially adjusting the delays by a multiple of the wavelength, these scans were simply discarded. Note that no a priori information was used to constrain the location of the central fringe, and each 75-second group was treated independently.

We also point out that it should be possible to co-add $`\chi ^2`$ functions not as a function of separation but as a function of differential right ascension and declination; such an approach should allow one to use the co-adding approach with very long integration times (e.g. for all the scans of a given star in a night). We are currently developing such an approach.

## 4 Results

We show the results of an astrometric fit to 45 minutes of data (taken over the course of 70 minutes of observation) in Figures 4 and 6. With 1769 scans used, we find the residual delay errors to be well modeled by a Gaussian distribution with a full-width at half-maximum of 0.160 $`\mu `$m. In order to characterize the residuals, and in particular determine if they could be considered to be independent, we plot the Allan variance of the residuals in Figure 5. The Allan variance (Thompson et al., 2001) at lag $`l`$ is given by

$$\sigma _A^2(l)=\frac{1}{2(M^{}+12l)}\underset{n=0}{\overset{M^{}2l}{}}\left(\frac{1}{l}\underset{m=0}{\overset{l1}{}}x_{n+m}x_{n+l+m}\right)^2$$
(9)
where $`M^{}`$ is the total number of data points. As can be seen in the figure, the residuals are white out to lags of more than 500 scans, implying a final astrometric precision of 10 $`\mu `$as. We list the results from 4 nights of observation in Table 1.

The 1-sigma error region (found by plotting the $`\chi ^2=\chi _{min}^2+\mathrm{\Delta }\chi ^2`$ contour, Press et al. 1992) is highly elliptical with the major axis oriented roughly parallel to the R.A. axis. Such error ellipses are to be expected in single-baseline interferometric data, which has limited sensitivity in the direction perpendicular to the baseline. It should be noted however that for sufficiently long observations, Earth-rotation will provide an orthogonal baseline. The major and minor axes of the uncertainty ellipse are easily found by diagonalizing the covariance matrix: the magnitude of the uncertainty in the direction of the minor axis was 8.3 $`\mu `$as for the 10 August data and 12 $`\mu `$as for the 11 August data; consistent with the delay residuals. The uncertainty in the major axis direction was 144 and 143 $`\mu `$as respectively. Figure 6 shows the four error ellipses superimposed. We fit the measured differential declinations to a linear trend of $`18\pm 7\mu `$as/day; the $`r.m.s.`$ of the residuals is $`16\mu `$as, and the reduced $`\chi ^2=1.1`$. The predicted change in differential declination due to orbital motion is approximately $`11\mu `$as/day.

## 5 Discussion

We have achieved an astrometric measurement precision of $`9`$ $`\mu `$as in 70 minutes of observation on the 0.25 arcsecond separation binary star HD 171779. The measurements were repeated on four nights, yielding answers consistent at the level of $`16\mu `$as. These early results show the promise of very narrow-angle astrometry for achieving very high precision measurements of stellar binaries with a separations on the order of 0.05–1 arcseconds. This is particularly useful in searching for planets around such binaries, although the resulting high-precision visual orbits for these stars should also be useful in providing stellar masses and distances at the level of $`10^3`$ . It should be noted that such high precision mass determinations may require an extended period of observation as these systems can have very long orbital periods. It is also necessary to obtain radial velocity data accurate to 10 m/s or better, which can be challenging.

This newly demonstrated approach to narrow-angle astrometry compares quite favorably with other astrometric methods used to date, such as speckle interferometry (Horsch et al., 2002; Woitas et al., 2003), single-telescope astrometry with (Lane et al., 2001) and without (Pravdo & Shaklan, 1996, 2000) adaptive-optics, and space-based astrometry using HST (Benedict et al., 2002a). These have all produced measurement precisions on the order of 0.1-1 mas. The improved performance of the interferometric approach is primarily due to the longer available baseline. However, the combination of phase-referencing and fringe-scanning also compares favorably to other interferometric methods, such as using fringe visibility measurements to find binary star orbits ( i.e. a fractional precision of $`10^4`$ as compared to $`10^2`$, Boden et al., 2000). In this case the improvement comes from the unbiased nature of the phase estimator.

Recent years have seen a veritable explosion in the number of known extra-solar planets (Schneider, 2003), starting with the famous example of 51 Peg (Mayor & Queloz, 1995). To date, most of these systems have been found using high-precision radial velocity techniques, although recent photometric transit searches have observed several planets (Charbonneau et al., 2000; Konacki et al., 2003). However, these methods have to a large extent avoided searching close (separation less than a few arcseconds) binary stellar systems, primarily because the techniques used are not well suited to such systems. Therefore it is particularly interesting to note that despite the deliberate avoidance of binaries, of the more than 100 known planets, 14 of them are in wide binary stellar systems. In addition, given the high frequency of binary stellar systems (57% among systems older than 1 Gyr, Duquennoy & Mayor (1991)), it is clear that any comprehensive planetary census must address the question of how frequently planets occur in such systems. This is all the more relevant given that several theoretical investigations have indicated that there exist regions in binary parameter space where planets can form and exist in stable orbits over long periods of time, though this does remain controversial (Whitmire et al., 1998; Boss, 1998; Marzari & Scholl, 2000; Nelson, 2000; Barbieri, Marzari & Scholl, 2002).

An instrument such as PTI, capable of 10 $`\mu `$as very-narrow-angle astrometry, could be used to search many of the brightest speckle binary systems for planets. We have compiled a list of approximately 50 suitable systems, with $`m_K<4.5`$, separations less than 1 arcsecond, and within the field of regard of PTI. The median orbital separation between the binary components in these systems is 19 AU, and hence there should be regions where planets can remain stable for long periods. In particular, adopting the result from Holman & Wiegert (1999) we calculate the largest stable orbit in each system. We find that the median detectable planetary mass in such an orbit is 0.5 Jupiter masses (assuming $`3\sigma `$ confidence detections). The corresponding median orbital period is 2.2 years. A limited survey could quickly begin to provide useful constraints on the frequency of planets in binary stellar systems. As a new generation of long-baseline optical interferometers become operational in the next few years, this type of survey could be easily extended to sample sizes of several hundred stars, hence providing strong constraints on planetary formation in the binary environment.

We would like to thank E. Bertschinger, A. Boden, B. F. Burke, M. Colavita, M. Konacki, S. R. Kulkarni & N. Safizadeh for their contributions to this effort. We also particularly acknowledge the extraordinary efforts of K. Rykoski, whose work in operating and maintaining PTI is invaluable and goes far beyond the call of duty. Observations with PTI are made possible through the efforts of the PTI Collaboration, which we acknowledge. Part of the work described in this paper was performed at the Jet Propulsion Laboratory under contract with the National Aeronautics and Space Administration. Interferometer data was obtained at the Palomar Observatory using the NASA Palomar Testbed Interferometer, supported by NASA contracts to the Jet Propulsion Laboratory. This research has made use of the Simbad database, operated at CDS, Strasbourg, France. MWM acknowledges the support of the Michelson Graduate Fellowship program. BFL acknowledges support from a Pappalardo Fellowship in Physics.


text:
# Vortices in the Co-orbital Region of an Embedded Protoplanet

## 1 INTRODUCTION

The discovery of a large number of short-period extrasolar planets has generated renewed interests in the study of tidally induced migration of protoplanets in disks. The basis of migration is the gravitational interaction between a gaseous disk and an orbiting perturber, and has been analyzed well before extrasolar planets were actually discovered (Goldreich & Tremaine, 1979, 1980; Lin & Papaloizou, 1986a, b). It has been suggested that this tidal interaction induces the embedded protoplanet to migrate under two different circumstances. In Type I migration, the protoplanet mass is small and the disk structure is weakly perturbed. In this limit, the disk response can be analyzed with a quasi-linear approximation. Linear calculations have revealed an intrinsic imbalance of inner and outer disk torques (Ward, 1997), so the planet could migrate quite quickly inward. In Type II migration, the protoplanet mass is sufficiently large for it to open up a gap in the disk (Papaloizou & Lin, 1984). The protoplanet’s orbital evolution is locked to the evolution of the disk, which is likely to be on the viscous timescale (Nelson et al., 2000). Both types of migration have been basically confirmed by the non-linear numerical simulations.

The high mobility of small mass protoplanets or embryo cores of giant planets from Type I migration has presented a challenge to the formation of giant planets. The absence of deep oceans is contrary to the inference that any terrestrial planet in the solar system may be formed well outside the snow-line and migrated to their present location. It is thus imperative to explore the cause which may have inhibited the Type I migration. One potential contributor to this process is the disk response to low mass planets, especially in the co-orbital region where only sketchy analytical arguments have been developed. There are, however, two notable recent exceptions. Masset (2001, 2002) has evaluated the co-orbital corotation torque and found that they can play a crucial role for the migration rate. Extensive numerical simulations of viscous disks have been performed to understand the torque dependence on the planet mass, sound speed, disk viscosity, and the torque saturation mechanism. One generic behavior he found is that the total torque on the planet first shows an oscillation consistent with the librating fluid motion in the co-orbital region then it levels off to a constant due to the disk viscosity.

Another study is by Balmforth & Korycansky (2001), who used matched asymptotic expansions to explore the non-linear dynamics of the flows in the co-orbital region. They investigated the role of re-arranging the potential vorticity profile as a response to the planet’s driving and its role in the corotation torque saturation. For inviscid flows, they also noticed the formation of vortices which they attributed to the secondary instabilities in a forced critical layer problem (the critical layer in this case is the corotation resonance).

In this Letter, we present results on the non-linear dynamics of the co-orbital region, based on direct, high resolution numerical simulations (§2). Our approach is quite similar to that of Masset’s except that we have focused on the inviscid limit simulations. As a consequence, the total torque on the planet exhibits some interesting new behavior, instead of leveling off to a constant as in the case of viscous disks. Some of the findings are similar to what Balmforth & Korycansky have found in the excitation of vortices (§3). We provide an interpretation of why these vortices are produced in §4 and discuss the implications of our results in §5.

## 2 Initial Set-up

We assume that the protoplanetary disk is thin and can be described by the inviscid 2-D Euler equations in a cylindrical $`(r,\varphi )`$ plane with vertically integrated quantities. Two gravitating bodies, a central star and a protoplanet, are located at $`r=0`$ and $`r=1`$ respectively and the 2-D disk is modeled between $`0.4r2`$. The planet is on a fixed circular orbit at $`r=1`$. The self-gravity of the disk is not included.

The mass ratio between the planet and the central star is taken to be $`\mu =M_p/M_{}=10^4`$. Its Hill (Roche) radius is $`r_H=(\mu /3)^{1/3}=0.032`$. The disk is assumed to be isothermal with a constant temperature throughout the simulation region (i.e., it is attached to a thermal bath). The isothermal sound speed, scaled by the Keplerian rotation speed $`v_\varphi `$ at $`r=1`$ is $`c_s/v_{\varphi ,r=1}=H/r`$ where $`H`$ is the disk scale height. We used a set of three sound speeds $`c_s=0.04,0.05,0.06`$. In the simulations presented here, we choose an initial surface density profile with $`\mathrm{\Sigma }(r)r^{3/2}`$, in order that the ratio of vorticity to surface density (potential vorticity) has a flat radial profile $`\zeta =(\times 𝐯)_z/\mathrm{\Sigma }const.`$ (Note that $`\zeta `$ has a small deviation from being a precise constant due to the finite pressure gradient which slightly modifies $`v_\varphi `$ from its Keplerian value.) With this choice of initial condition, we avoid the generation of inflection points due to the rearrangement of $`\zeta `$ distribution as it is carried by the stream lines. We have also made other runs which have an initial potential vorticity gradient and we find that the results presented here also hold for them. The detailed analysis of these results will be presented elsewhere.

Since $`r_H<H`$ for all three different sound speeds, the planet is embedded. The planet’s gravitational potential is softened by an approximate 3-D treatment. The key approach is to take the surface density at each cell spreading it out vertically based on a Gaussian density profile with a scale height determined by the local temperature. The gravitional force between the planet and the disk matter is calculated using this 3-D configuration but all the forces with the same $`(r,\varphi )`$ positions are summed together. This gives a correction factor $`C(r,\varphi )`$ to the gravitational forces, which effectively softens the gravitational potential. The torque overestimation in a 2-D geometry was studied in D’Angelo, Kley, & Henning (2003) which arises only from locations near the planet. As a test if our key finding is affected by the 3-D treatment, we excluded several Roche radii from the torque calculation and show that it is not the case.

Although there are some accumulation of gas near the planet, we do not allow any disk gas to be accreted onto it. This boundary condition is equivalent to that of sufficiently low-mass protoplanets with quasi hydrostatic envelopes (Pollack et al. 1996). The planet’s potential is switched on over 10 orbits in order for the disk to make an adiabatic adjustment. (We have used other, longer turn-on times but the results do not depend on them.)

Simulations are carried out using a 2-D code whose basic algorithm was based on FARGO by Masset (2000) and Li et al. (2001). The co-rotating frame is used so that the positions of the central star and the planet are fixed at $`(r,\varphi )=(0,0)`$ and $`(1,0)`$ (acceleration due to frame rotation is also included, see Kley (1998)). Runs are made using several radial and azimuthal grids to study the influence of resolution, and those presented here have $`(n_r\times n_\varphi )=300\times 1200`$. Outflow boundary conditions are used at both the inner and outer radial boundaries. The simulations typically last several hundred orbits at $`r=1`$.

One key difference between our simulations and previous studies is that our simulations are performed in the inviscid limit, i.e., we do not explicitly include a viscosity term, though numerical viscosity is inevitable and is needed to handle for example shocks. In addition, our simulations tend to be of higher resolution. With a radial $`n_r=300`$ (some runs go up to $`n_r=600,n_\varphi =2400`$) and $`0.4r2`$, the diameter of the Hill sphere of a planet with $`\mu =10^4`$ is resolved by at least $`12`$ cells in each direction. Using a nested-grid technique, D’Angelo et al. (2002) have also made very high resolution simulations, though their simulations include an explicit viscosity.

## 3 Torque Evolution

Our simulation results are generally in agreement with previous simulations for the planet on a fixed circular orbit (see Masset, 2002) in terms of the generic structures produced in the disk by the tidal interactions between the planet and its surrounding flows. Since the planet mass is relatively small, only density dips are produced without gaps (at least during a few hundred orbits simulated here). We define the co-orbital region here as $`|\mathrm{\Delta }r|=|r1|6r_H`$. The co-orbital region can be separated into several parts depending on their streamline behavior (see Masset, 2002), including the horseshoe (or librating) region $`|\mathrm{\Delta }r|<r_H`$, the separatrix region $`r_H<|\mathrm{\Delta }r|<\sqrt{12}r_H`$, and the streaming region $`|\mathrm{\Delta }r|>\sqrt{12}r_H`$. The horseshoe region can be further divided into real horseshoe streamlines enclosing all three libration points ($`L_{3,4,5}`$) and tadpole like streamlines enclosing either $`L_4`$ or $`L_5`$.

Our key finding is shown in Fig. 1 which displays the total torque on the planet from the disk. It has two very distinct phases: Phase I: the total torque is negative and smooth, modulated by a period that is similar to the libration period near the L4 and L5 points ($`\tau _{lib}\sqrt{4/27\mu }40P`$ where $`P`$ is one orbit at $`r=1`$). Phase II: the total torque shows very large amplitude and fast oscillations with a quasi-period of a few orbits. We see such a behavior in all the runs with three different sound speeds, though the times when the “phase” transition occurs are different. They are at $`t40P,90P,`$ and $`210P`$ for $`c_s=0.04,0.05`$ and $`0.06`$, respectively.

## 4 Shocks, Potential Vorticity Evolution and Instabilities

The obvious question is what causes such large amplitude and fast changes in torques. Closer inspection of the co-orbital region at times around $`90P`$ (using $`c_s=0.05`$ as an example) reveals that vortices are generated in the separatrix region. This development is shown in Fig. 2, which depicts maps of potential vorticity $`\zeta `$ through the times when the system undergoes its “phase” transition. One can see that these vortices are anti-cyclones which are also blobs of higher densities (not shown here). Being in the separatrix region, they experience very close encounters with the planet repeatedly, thus exerting strong “impulses” to the planet. This results in the large amplitude and fast variations in the torque evolution.

To understand this phenomenon further, in Fig. 3, we plot the azimuthally averaged radial potential vorticity profile $`\zeta `$ for $`c_s=0.05`$ at different times $`t=0,20P,40P,80P,120P`$. One sees that this profile deviates progressively away from the initially flat radial profile, developing two prominent structures at $`|\mathrm{\Delta }r|23r_H`$, with a minimum at $`\mathrm{\Delta }r2.4r_H`$ and a maximum at $`\mathrm{\Delta }r=3.4r_H`$. Comparing with Fig. 2, vortices apparently emerge from the minima.

The existence of the extrema (or inflection points) in the potential vorticity profile is usually regarded as a necessary condition for non-axisymmetric instabilities in rotating shear flows (Drazin & Reid 1981). The progressive steepening of this profile, say, at $`\mathrm{\Delta }r2.4r_H`$ until $`t<90P`$ and the subsequent flattening (see the curve at $`t=120P`$) strongly suggest that there is a secondary instability which is excited around $`\mathrm{\Delta }r2.4r_H`$ and vortices are the non-linear outcome of such an instability.

How could the potential vorticity profile develop such inflection points/structures from an initially flat and stable profile? As can be seen from the evolution equation for potential vorticity

$$D(\zeta \widehat{z})/Dt=(\zeta \widehat{z})𝐯=0\mathrm{in}2\mathrm{D},$$
(1)
where $`D/Dt`$ is the Lagrangian derivative, that $`\zeta `$ is conserved along the streamlines in 2D (Korycansky & Papaloizou 1996). Thus, the initial uniform $`\zeta `$ distribution should be preserved in an inviscid disk. This conclusion, however, is only true provided that streamlines can be defined everywhere in a dissipationless flow.

Disk flows in the co-orbital region around the planet, however, do not satisfy this condition (except perhaps for flows within $`|\mathrm{\Delta }r|<r_H`$) because of the two sets of spiral shocks produced by the planet. These two sets of shocks emanate from close to the planet and cut through the whole disk. Flow lines which encounter these shocks are broken and the dissipation in shocks breaks the potential vorticity conservation. The consequence of the shock is equivalent to that of a strong viscous stress in the momentum equation, the curl of which, does not vanish.

To establish the role of shocks in breaking the potential vorticity conservation, we have performed three runs with different sound speeds $`c_s=0.04,0.05,0.06`$. Even though they only span a small range, but the effects are quite clear. In Fig. 4, we plot the azimuthally averaged potential vorticity profiles for these three runs at a time just before the emergence of vortices. Vortices first appear at the minimum of the “valley” in the profiles which moves outwards as $`c_s`$ increases. This growth pattern is consistent with the fact that the starting location of the shocks is also moving away from the planet when the sound speed increases.

Note that the cause for the emergence and growth of inflection points in the azimuthally averaged $`\zeta `$ distribution presented here is different from that due to its re-arrangment along the horseshoe orbits (Balmforth and Korycansky 2001). In our calculation with a more general initial $`\mathrm{\Sigma }`$ (and $`\zeta `$) distribution, we confirm that onset of the shearing instability occurs at an earlier epoch. The appearance of the subsequently self-excited vortices is quite similar to the development of the Kelvin-Helmholtz (KH) instability observed at shear interfaces. A detailed study on the properties of this secondary instability, such as its threshold and growth rate, will be presented in a forthcoming paper.

## 5 Conclusions and Discussions

In this paper, we consider the tidal interaction between a protostellar disk and an embedded planet with a modest mass. Linear torque calculations suggest that such interaction leads to angular momentum transfer and inward migration of the planets. For planets with mass comparable to that of the Earth, the inferred migration time scale (Ward 1997) is much shorter than the observationally determined depletion time scale of the disk (Haisch, Lada, & Lada, 2001). This analytic result raises a problem for the formation of giant planets through the core accretion process (Pollack et al. 1996). We reinvestigate this process in an attempt to resolve this paradox.

The embedded protoplanet induces a circulation flow pattern near its corotation region in the disk. We show here that even in disks with initially uniform potential vorticities, inflection points may emerge and grow spontaneously as potential vorticity is generated near the spiral shock in the vicinity of the protoplanet. These inflection points lead to the onset and growth of secondary instabilities.

The existence of secondary instabilities and their non-linear outcome as vortices indicate that the flows in the co-orbital region are perhaps more complicated than linear analyses have suggested. These vortices should exert strong torques on the planet.

One limitation of our current study is that the planet is being artificially held on a fixed circular orbit. Such strong and rapidly varying torques from the vortices on the planet might cause the planet to “oscillate” and lose its phase coherence with the surrounding flow. This feedback process should have important implications for the Type I migration problem. Allowing the planet to migrate radially must be included in order to understand the planet’s response to such torques. These issues will be addressed in a forthcoming publication.

We wish to thank Neil Balmforth, Darryl Holm and Edison Liang for useful discussions. This research was performed under the auspices of the Department of Energy. It was supported by the Laboratory Directed Research and Development Program at Los Alamos and by LANL/IGPP. We are also supported in part by NASA through NAG5-11779, NAG5-9223, and NSF through AST-9987417.