MILLIONS OF MULTIPLES: DETECTING AND CHARACTERIZING CLOSE-SEPARATION BINARY SYSTEMS IN SYNOPTIC SKY SURVEYS

Hoekstra et al. (2005, hereafter HWU) describes a method to measure the ellipticity of a single object with corrections for PSF anisotropy, which we summarize here. HWU quantify the ellipticity of a PSF with ellipticity parameters e₁ and e₂, which relate the second moments of flux as measured from the centroid of the object:

$$\begin{equation} e_{1} =\frac{I_{11} - I_{22}}{I_{11} + I_{22}},\quad e_{2} =\frac{2I_{12} }{I_{11} + I_{22}}, \end{equation} \tag{ 1 }$$

where

$$\begin{equation} I_{11} = \sum \limits _{x,y}{f(x,y)x^{2}W}, \end{equation} \tag{ 2 }$$

$$\begin{equation} I_{22} = \sum \limits _{x,y}{f(x,y)y^{2}W}, \end{equation} \tag{ 3 }$$

$$\begin{equation} I_{12} = \sum \limits _{x,y}{f(x,y)xyW}. \end{equation} \tag{ 4 }$$

Here, x and y are the x and y pixel numbers as measured in a coordinate system whose origin is at the centroid of the object, f(x, y) is the source flux falling on the pixel located at the point (x, y) in that same coordinate system, and W is a weight function which accounts for the increase in Poisson noise with decreasing photon count away from the centroid.

The ellipticity of the PSF of an unresolved binary or an extended object is a combined measure of the true ellipticity caused by the astrophysical configuration of the object, and the ellipticity induced by directional smearing of starlight by the atmosphere and telescope optics. Thus, any raw ellipticity parameters need to be corrected for the induced ellipticity in order to obtain astrophysically significant measures of the true ellipticity. HWU do this using the smear polarizability tensor P, which quantifies the changes in PSF ellipticity resulting from perturbations in the PSF anisotropy:

$$\begin{equation} P_{ij}=X_{ij}-e_{i}e_{j}^{sm}, \end{equation} \tag{ 5 }$$

where

$$\begin{equation} X_{11}\,{=}\,\frac{1}{I_{11} \,{+}\, I_{22}} \sum \limits _{x,y}{f(x,y)(W\,{+}\,2W^{\prime }(x^{2}\,{+}\,y^{2})+W^{\prime \prime }(x^{2}\,{-}\,y^{2})),} \end{equation} \tag{ 6 }$$

$$\begin{equation} X_{22}=\frac{1}{I_{11} + I_{22}} \sum \limits _{x,y}{f(x,y)(W+2W^{\prime }(x^{2}\,{+}\,y^{2})+4W^{\prime \prime }x^{2}y^{2}),} \end{equation} \tag{ 7 }$$

$$\begin{equation} e_{1}^{sm}=\frac{1}{I_{11} + I_{22}} \sum \limits _{x,y}{f(x,y)(x^{2}-y^{2})(2W^{\prime }+W^{\prime \prime }(x^{2}+y^{2})),} \end{equation} \tag{ 8 }$$

$$\begin{equation} e_{2}^{sm}=\frac{1}{I_{11} + I_{22}} \sum \limits _{x,y}{f(x,y)(2xy)(2W^{\prime }+W^{\prime \prime }(x^{2}+y^{2})).} \end{equation} \tag{ 9 }$$

Here, the differentiation of the weight function W is taken to be with respect to x² + y².

The smear polarizability tensor of a source allows HWU to measure the PSF anisotropy parameters p₁ and p₂, which quantify the anisotropy in the image:

$$\begin{equation} p_{i}=\frac{e_{i}}{P_{ii}}. \end{equation} \tag{ 10 }$$

To correct the ellipticity of a star in an image, it is necessary to examine how the anisotropy parameters vary across the CCD in the vicinity of the star. HWU measure the anisotropy parameters of point sources near the star of interest, and create a polynomial fit of p₁ and p₂ as a function of position on the CCD. They then correct the observed ellipticity of the star via the relationship

$$\begin{equation} e_{i}^{{\rm cor}}=e_{i}^{{\rm obs}}-P_{ii}p_{i}-\alpha \left(p_{1}^{2}+p_{2}^{2}\right)^{0.5}, \end{equation} \tag{ 11 }$$

where α is a constant proportional to the magnitude of the anisotropy at the point on the CCD where the star is located. This provides a final, anisotropy-corrected measure of the PSF ellipticity.

HWU measure the ellipticity of the object in many separate observations. The observations are taken at varying seeings, and thus they determine the ellipticity at a reference FWHM by deriving a fit of ellipticity as a function of FWHM and evaluating the fit at the chosen reference radius.

In order to measure the necessary parameters of all the stars in our images (using Equations (1)–(10)), we require centroid coordinates and FWHM measurements of the PSF of each source. We obtain a catalog of objects (these include stars, blended binary systems, foreground/background star blends, and galaxies) for each image using SExtractor (Bertin & Arnouts 1996). SExtractor provides the sky coordinates (R.A. and decl.) and pixel coordinates (x and y) of the centroid of each object, the FWHM assuming a Gaussian core, the flux, and the SExtractor error flag. We remove all objects with non-zero SExtractor flags to avoid saturation, bad pixels, and the wider-separation binaries which SExtractor can directly detect and deblend.

Using the catalog of objects obtained from SExtractor, we measure the raw e₁ and e₂ ellipticity parameters and the PSF anisotropy parameters p₁ and p₂ of every object in an image. The weight function we employ to account for changing signal-to-noise ratio (S/N) across the PSF is a Gaussian with dispersion equal to the FWHM calculated by SExtractor, centered at the SExtractor-derived centroid coordinates. We define an object aperture with a radius equal to twice the FWHM, and a sky-subtraction aperture with radii between 2.5 times and 3.5 times the FWHM.

The anisotropy profile of images can vary significantly in wide-field images due to geometric distortions. In PTF data, we found it necessary to partition the images into subsections, and derive separate fits for each partition (see Figure 1). It is important to make the regions as close to square in shape as possible, as an accurate fit requires a comparable span in each direction. Additionally, the subset of sources used to create the fits should always extend further than the boundaries of the region being corrected (as shown by the blue square in Figure 1), to ensure that the ellipticities of all sources are corrected based on a fit obtained from sources surrounding them on all sides.

Zoom In Zoom Out Reset image size

Following the method established by HWU, we derive a polynomial fit of the anisotropy parameters for each group to determine how the PSF anisotropy varies across the CCD. We group objects with similar flux levels together, and we carry out a separate PSF anisotropy correction for each group. This makes using S/N weighted fits unnecessary as the similar flux levels within each group imply similar S/Ns and thus similar measurement errors. We perform a 3σ clip on the values of e₁, e₂, p₁, and p₁ to remove severe outliers. The points that remain may not all be point sources (some are blends, galaxies, etc.). In contrast to instrumentally induced ellipticity, this will not cause any directional bias in the fit, as we expect astrophysically extended objects to appear with equal probabilities in all orientations.

Our anisotropy fit is of the form

$$\begin{equation} p_{1}^{{\rm fit}}=c_{1}x\,{+}\,c_{2}x^{2}\,{+}\,c_{3}y\,{+}\,c_{4}y^{2}\,{+}\,c_{5}xy\,{+}\,c_{6}\,{+}\,c_{7}f_{{\rm tot}}\,{+}\,c_{8}f_{{\rm tot}}^{2}, \end{equation} \tag{ 12 }$$

where c_i are the constant coefficients to be found which minimize the residuals of the measured p₁ parameters of the objects being fitted and x and y are the CCD pixel coordinates where the objects are located.

We extended the Hoekstra et al. (2005) model to include terms which are proportional to the flux of the objects (f_tot). While applying the algorithm to the PTF images, we noticed that the ellipticities of high-flux objects were not corrected accurately to zero using a fit dependent only on x and y object coordinates. Upon closer examination of the images, we realized that for high flux (yet unsaturated) objects, slightly less than 100% charge transfer efficiency (CTE) leads to a smearing in the direction of readout. This causes an additional elongation dependent on the flux of the object, requiring flux terms in the anisotropy fit.

A second-order polynomial fit appears to be sufficient for PTF images. We found that, when using a third degree instead of second degree fit of the p₁ and p₂ parameters, the resulting change in the corrected ellipticities of the objects in each image was approximately nine times smaller than if we used a second degree fit but varied the raw ellipticities and PSF anisotropy parameters by random values consistent with their measurement error. Thus, for the PTF data set, the measurement error of the raw ellipticities and anisotropy parameters was the dominant source of error in the corrected ellipticities, and increasing the fit order does not significantly increase the accuracy of the corrected ellipticities.

A similar fit is also derived for the p₂ parameter, although the flux-dependent terms are not necessary as readout-induced elongation exactly parallel to the pixel grid (whether it is the x or y axis) does not affect the e₂ parameter, by definition (see Equations (1)–(4)).

Once we have derived these fits we evaluate them at the coordinates of the objects to be corrected and obtain the corrected ellipticities according to Equation (12). For the purposes of evaluating the ellipticity of an object compared to others in the field, we use e, the $\sqrt{\vphantom{A^A}\smash{{{e_1^2 + e_2^2}}}}$ combination of the ellipticities.

We find that the typical instrument-induced ellipticity has a magnitude between 0.01 and 0.1, roughly the same amount as our final measured ellipticities for detected binary systems.

When we measure the corrected ellipticities of an object across multiple images, we measure many values for e₁ and e₂ at a variety of FWHMs. However, ellipticity is not constant with seeing; close objects induce smaller ellipticities when the PSFs are large. To make ellipticity measurements comparable between fields with different seeing distributions, we create a fit of object ellipticities as a function of FWHM, and evaluate this fit at a reference seeing value (Figure 2).

Zoom In Zoom Out Reset image size

Undersampled images with FWHMs close to the pixel scale of the CCD cannot accurately measure ellipticities because the ellipticity measurement is then very sensitive to the location of the object relative to the pixel grid. We found that pixel effects become negligible for roughly Nyquist-sampled images with seeing FWHM of two pixels of greater, and we therefore use only those well-sampled PTF images. We require at least 30 separate observations of a target in order for us to obtain an accurate fit of ellipticity with seeing, and we also require several hundred sufficiently bright sources within our image for PSF measurement (the last requirement is satisfied by essentially all PTF images). For R-band PTF images, our simulations and on-sky measurements showed that the faintest useful sources have m_R ∼ 18.

Figure 3 provides an example of the distribution of e₁, e₂, and e (all measured at a common reference FWHM) of the objects appearing in a single chip of a PTF field (1/12 of a PTF field). The distributions have two components—the bulk of the objects which are closely centered around zero ellipticity, and the extended wings. The first component represents true point sources or blends of objects too close together for our algorithm to differentiate from point sources. The thickness of this component of the distribution is due to measurement errors in the ellipticities. The wings of the distribution represent either blends of multiple point sources, or other extended objects such as galaxies. Figure 4 shows a sample of PTF objects from the same chip, ordered by increasing ellipticity.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Robo-AO is a visible and near-infrared laser guide star adaptive optics system specifically engineered for 1–3 m class telescopes (Baranec et al. 2012). The Robo-AO system located at the Palomar 60 inch (150 cm) telescope comprises an ultraviolet Rayleigh laser guide star, an integrated adaptive optics and science camera system, and a robotic control system. The system currently incorporates both an electron-multiplying CCD and an InGaAs infrared array camera for imaging. The robotic system is designed for high-efficiency observing, allowing large numbers of targets to be imaged rapidly. The diffraction-limited visible-light capability of Robo-AO allowed us to make a similar wavelength comparison between the PTF-measured ellipticities and the observed binary contrast ratios.

We obtained Robo-AO images of the 44 targets on the nights of 2012 July 16–18 and August 4–7 (UT). The July targets, with R.A.s around 13 hr, were observed with 60 s total exposure times in the i-band filter; the August targets (R.A.s of 22–23 hr) used a long-pass filter with a 600 nm cut-on to obtain increased signal compared to a bandpass filter. We operated Robo-AO without tip-tilt correction, instead relying on post-facto shift-and-add processing of the individual frames and used the pipeline described by Law et al. (2009b, 2012) to perform the image alignment and co-addition. The Robo-AO system provided diffraction-limited resolution (∼0.1 arcsec at these visible wavelengths) for all observed targets, in 1–2 arcsec seeing conditions; the entire target list was observed in a total of only ∼2.5 hr.

The high-resolution images of our test sample (Figure 8) confirm the expected variation of binarity with measured ellipticity. Targets with small ellipticities (e  < 0.01) are all confirmed to be single stars. We found that five of the six targets with 0.01 < e < 0.02 are binaries; all the systems in that range are high-contrast and/or <1 arcsec separation. Targets with e > 0.02 are all confirmed to be binaries, with separations increasing and contrast ratios decreasing as the ellipticities increase. A single target, PTF13.341, is detected as a binary with a faint, close companion by Robo-AO, but with a different position angle and closer separation than that measured from PTF data; this is likely to be a close companion not measured by our PTF data and we do not include it in further analysis.

Zoom In Zoom Out Reset image size

The crossover point between single stars and high-confidence binaries occurs at the ellipticities predicted by our artificial point source simulations (Section 3), although the >80% binarity fraction of the targets with ellipticities between 0.01 and 0.02, compared to the zero binarity fraction at lower ellipticities, suggests that a less conservative limit could be set for many science programs.

The detected Robo-AO binaries show that our ellipticity measurements can predict the orientation of the binary from the relative magnitudes of the e₁ and e₂ parameters (see the ellipses in Figure 8). The e₁ parameter varies as cos(2θ) and the e₂ parameter varies as sin(2θ), where θ is the angle that the vector connecting the binaries makes with the horizontal. This makes it possible to predict θ via

$$\begin{equation} \theta = \frac{\tan ^{-1}\left(\frac{e_{2}}{e_{1}}\right)}{2}. \end{equation} \tag{ 13 }$$

Since the bulk of the central distribution of e₁ and e₂ consists of point sources with zero ellipticity, the scatter of the central distribution is a good estimate for the measurement error of the ellipticity parameters. We fitted a two-dimensional Gaussian around the central distribution and found its variance to determine the measurement error of the parameters. We propagated the error to the angle θ according to Equation (13) to find the uncertainty in position angle. For the PTF set, at the minimum definite binary detection ellipticity, the predicted binary orientation has a 1σ uncertainty of 2° for our brightest objects (14–15 mag) and 4° for the least bright objects (17–18 mag). As the binary separation and contrast ratios increase, the uncertainty decreases, down to as low as 1° for the larger separations we have measured. In our test data set, where the objects all had high S/N in the PTF data set, the median difference between the PTF-predicted orientation and the observed Robo-AO orientation is 21, which closely matches the derived measurement error in the position angle.

The ellipticity of a close binary increases when the flux ratio of the stars (luminosity of the fainter star/luminosity of the brighter star) in the binary system increases, or alternatively when the separation between the binaries increases (Figure 5). Because of this degeneracy between separation and contrast, it is not possible to estimate the binary parameters by examining solely the magnitude of the ellipticity. However, we can differentiate between the two cases based on how the ellipticity varies with seeing.

When the members of a binary are very close together, the image of their combined light will only be noticeably elliptical for smaller FWHM values. As the FWHM increases, the image will quickly start to resemble a circular PSF. On the other hand, when the two members of a binary are separated by a large angular distance, their combined image will remain elliptical for much larger FWHMs than the close binaries. Thus, if the magnitude of the ellipticity of a binary is relatively constant with seeing, the binary has a wide separation, and if the ellipticity decreases quickly with seeing, the binary members must be close together.

For example, two of the binaries from Figure 8, PTF22.573 and PTF23.513, have nearly identical ellipticities at the chosen reference FWHM of 2.5 arcsec (as measured from the PTF images). The high-resolution Robo-AO images reveal that these two binaries differ in contrast and separation—PTF22.573 has greater angular separation and lower contrast ratio than PTF23.513. Figure 9 shows the seeing dependence of the ellipticity of both of these binaries. It is apparent than the closer, higher contrast binary has an ellipticity which is much more sensitive to the FWHM, as expected. This difference in response to changes in seeing allows us to resolve the degeneracy in the ellipticity measurement between separation and contrast.

Zoom In Zoom Out Reset image size

For the binaries we found with Robo-AO, we examined how the measured ellipticity in the PTF images varied with seeing (Figure 10). We derived a simple linear fit of the ellipticity as a function of seeing for each of the binaries. We evaluated these fits at a reference FWHM of 2.5 arcsec, and we parameterized the response of the ellipticity to seeing with the fractional decrease in seeing per arcsecond (Δe/e, where Δe is the magnitude of the slope of the linear fit of ellipticity with respect to seeing, and e is the ellipticity measured at the reference radius of 2.5 arcsec). As can be seen in Figure 10, the binaries with widest angular separations have close to zero variance of ellipticity with seeing. Using the measured relation to predict the Robo-AO measured separation on the basis of PTF data alone confirms that BinaryFinder can measure the binary separations to ∼25% precision.

Zoom In Zoom Out Reset image size

The measured ellipticity is a function of flux ratio and separation only, and so with a measured separation the flux ratio of the binary can also be estimated. Also in Figure 10, we use the estimated separations to constrain the flux ratios of the systems, assuming a simple linear relationship between separation and ellipticity independent of contrast ratio (which is approximately correct; see Section 3). The resulting distribution is (weakly) correlated with the Robo-AO measured flux ratios; although there are outlier points, a low value of ellipticity over separation strongly suggests a high flux ratio between the components of the binary systems.

False positives can be generated in our algorithm by instrumental effects that mimic PSF ellipticity, or extended astrophysical objects that are not physically associated binaries.

The careful attention to input data quality and PSF anisotropy removal results in a low level of instrumental false positives. Of the 44 targets observed with Robo-AO, 25 had ellipticities greater than our blend detection criterion of 0.02 ellipticity. All 25 of these targets were confirmed to be binary systems, suggesting a false-positive rate of less than 11% with 95% confidence.

Proper motion sample cuts can ensure that the BinaryFinder targets are relatively nearby stars rather than galaxies and other extended extragalactic objects. We note, however, that if the proper motion cuts are relaxed, the technique detailed in this paper can detect and measure the characteristics of any objects that lead to an elliptically shaped image. Beyond binary detection, the algorithm could thus be used to search for and statistically characterize faint barely extended galaxies, strong gravitational lenses and other extragalactic objects.

We have demonstrated that BinaryFinder can detect binaries with separations down to ≈1/5 of the seeing limit in wide-field survey data, and directly measure their position angles, separations and contrast ratios. When applied to the synoptic sky survey data currently being collected in great quantities by synoptic sky surveys like PTF and Pan-STARRS, we will be able to perform a search for binary systems over very large sky areas. This can be performed in targeted fields—such as surveying the entire Kepler field (Borucki et al. 2010) for blended binaries which could lead to false-positive exoplanet detections (e.g., Torres et al. 2004)—or BinaryFinder can be applied to the entire sky to build very large samples of binary stars.

At the time of writing PTF has imaged 1200 fields with sufficient numbers of epochs for full BinaryFinder binary detection and measurement of position angles and separations. Averaging across the sky, each PTF field contains ≈10,000 objects which are bright enough (m_R < 18) for the detection of companions at separations down to ∼1/5 of the seeing limit. In typical fields, we find that approximately 35% of those objects have high enough proper motions to confirm that they are nearby stars, and ≈11% of those are detected as <2 arcsec separation binaries by our algorithm (either binaries or blended unassociated foreground/background stars). Estimating the background blend probability from the stellar density in our images, we estimate that only ≈1% of those targets will be background blends.

In total, a full PTF-data set BinaryFinder search could target ∼12,000,000 stars; we would expect to detect ∼450,000 binary systems with detectable proper motions and a high confidence of physical association. A binary sample of that size will allow unprecedentedly detailed statistical analysis of the binarity fraction as a function of stellar mass, age, metallicity, and galactic population.

PTF and similar surveys operate in multiple filters. Multi-color imaging ellipticity measurements, if available, would scale according to the relative colors of the binary components, and could be used to reduce the false-positive probability or further constrain the nature of detected systems.

Many current and planned synoptic sky surveys image the same areas of sky repeatedly over many years. With the ability to directly measure position angle, separation, and contrast ratio of individual binaries, in many cases sufficient data will be collected to follow the motion of binary systems with time. Since this can be done over large areas and very large samples, BinaryFinder could be used to obtain orbital measurements (and thus mass constraints) for a very large sample of binaries.

We have shown that BinaryFinder is capable of detecting and characterizing close binary systems in wide-field synoptic survey data. Using Robo-AO, we have confirmed the PTF direct detection and separation measurement of binaries with separations as small as 0.4 arcsec, or less than one-fifth of the PTF FWHM. Images from the PTF camera have a pixel sampling of 1.01 arcsec and a median FWHM of ≈2.0 arcsec. Because our algorithm requires well-sampled images, we restricted its operation to PTF images with 2–3 arcsec FWHMs. New synoptic sky surveys with smaller pixels and smaller PSFs will enable the detection of much closer binary systems; for well-sampled images the detection limit scales with the PSF FWHM (assuming similar systematic noise performance). Using BinaryFinder, a survey with improved sampling such as LSST is likely to be capable of the detection, characterization and orbital motion measurement of millions of multiple systems with separations as close as 0.1 arcsec.