Time Delay and Accretion Disk Size Measurements in the Lensed Quasar SBS 0909+532 from Multiwavelength Microlensing Analysis

We present three complete seasons and two half-seasons of SDSS r-band photometry of the gravitationally lensed quasar SBS 0909+532 from the U.S. Naval Observatory, as well as two seasons each of SDSS g-band and r-band monitoring from the Liverpool Robotic Telescope. Using Monte Carlo simulations to simultaneously measure the system's time delay and model the r-band microlensing variability, we confirm and significantly refine the precision of the system's time delay to \Delta t_{AB} = 50^{+2}_{-4} days, where the stated uncertainties represent the bounds of the formal 1\sigma\ confidence interval. There may be a conflict between the time delay measurement and a lens consisting of a single galaxy. While models based on the Hubble Space Telescope astrometry and a relatively compact stellar distribution can reproduce the observed delay, the models have somewhat less dark matter than we would typically expect. We also carry out a joint analysis of the microlensing variability in the r- and g-bands to constrain the size of the quasar's continuum source at these wavelengths, obtaining log[(r_{s,r}/cm) [cos{i}/0.5]^{1/2}] = 15.3 \pm 0.3 and log[(r_{s,g}/cm) [cos{i}/0.5]^{1/2}] = 14.8 \pm 0.9, respectively. Our current results do not formally constrain the temperature profile of the accretion disk but are consistent with the expectations of standard thin disk theory.


INTRODUCTION
Much of the standard picture for the detailed structure of accretion disks surrounding supermassive black holes in active galactic nuclei is based on theoretical models rather than observational measurements, because these very compact regions cannot be resolved with existing telescopes. Yet, for gravitationally lensed quasars, the relative motions of the observer, the background source, the foreground lens galaxy, and its stars cause uncorrelated variations in the source magnification as a function of time and wavelength which depend on the projected area of the continuum source. By analyzing these microlensing brightness fluctuations with numerical simulations, one can measure the continuum source size, permitting invaluable observational tests of theoretical models of accretion disk structure (e.g., Eigenbrod et al. 2008;Anguita et al. 2008;Poindexter et al. 2008) and orientation (Poindexter & Kochanek 2010). Such tests in the literature typically fall into two categories. One type of study employs single-epoch multi-band photom-etry of lensed quasars in which flux ratios of the images exhibit deviations from the predictions of macroscopic lens models or significant wavelength dependence (e.g., Pooley et al. 2006;Bate et al. 2008;Blackburne et al. 2011;Motta et al. 2012). The alternative method, described in detail in Kochanek (2004), analyzes the time variability of the quasar's flux ratio and requires monitoring a quasar over a significant period of time. Although such investigations can be observationally and computationally challenging, they have the advantage that it is not necessary to assume a value for the mass of the microlenses or an extinction law for the lens galaxy, nor are the results highly sensitive to assumed priors. Studies employing both methods of analysis have revealed accretion disk temperature profiles in agreement with the simple thin-disk model of Shakura & Sunyaev (1973), although the individual observed disk sizes tend to be larger than those predicted by the theoretical model (e.g. Pooley et al. 2007; Morgan et al. 2010;Blackburne et al. 2011).
SBS 0909+532 (hereafter SBS 0909; α J2000 = 09 h 13 m 01. s 05, δ J2000 = +52 d 59 m 28. s 83) is a doublyimaged quasar lens system in which the background quasar has redshift z s = 1.377 and the foreground earlytype lens galaxy has redshift z l = 0.830 (Kochanek et al. 1997;Oscoz et al. 1997;Lubin et al. 2000). SBS 0909 is a somewhat challenging system to study because discrepant results in the literature expose significant uncertainties about some of its fundamental properties and the nature of the quasar's variability. Most notably, Lehár et al. (2000) encountered difficulty when attempting to use the imfitfits routine to measure the lens galaxy photometry in H-band NICMOS images from the Hubble Space Telescope (HST ), eventually settling and B (bottom) including measurements from MDM Observatory (stars), the WIYN 3.5-m telescope (diagonal crosses), the Liverpool Telescope (triangles), and USNO (squares). The measurements for image A have been offset by +0.2 mag to minimize empty space in the plot area. The light curve of image B exhibits a substantially steeper slope over the time period 4500 HJD − 2450000 5200 days. on a low surface brightness de Vaucouleurs model with effective radius r eff = 1. ′′ 58 ± 0. ′′ 9 and magnitude H = 16.75. In stark contrast, Sluse et al. (2012) used an iterative deconvolution technique (e.g., Magain et al. 1998;Chantry & Magain 2007) on the same data to find a significantly smaller (r eff = 0. ′′ 54 ± 0. ′′ 02) and less luminous (H = 19.44 ± 0.01) lens galaxy, formally inconsistent with the Lehár et al. (2000) result. Additionally, past optical monitoring observations have not shown evidence of significant microlensing variability (Ullán et al. 2006;Goicoechea et al. 2008), although Mediavilla et al. (2011) find evidence of chromatic microlensing (wavelength-dependent microlensing magnification) through an analysis of the differences between the continuum and emission line flux ratios from the quasar's ultraviolet (UV) -near-infrared (near-IR) spectra, which can separate microlensing effects from differential extinction present in the lens galaxy.
We have compiled a new data set consisting of monitoring observations of SBS 0909 in two optical bands from two different observatories, the analysis of which provides some resolution to these discrepant results from the literature. Here we will show that in the four years spanning 2008 -2012, SBS 0909 has exhibited significant uncorrelated time variability in the rest-frame near-UV. We will analyze the uncorrelated variability, which we attribute to microlensing, to determine a size for the accretion disk in the two different bands, and then compare our results to that derived from the observed chromatic microlensing by Mediavilla et al. (2011). Since our multi-band data set allows us to constrain the size of the accretion disk at two different wavelengths, we also gain a glimpse at The g-band light curves exhibit similar intrinsic variability to the r-band light curves over the same period of time, although with increased scatter due to the lower quasar flux in g and poorer observing conditions on some occasions. the temperature profile of the accretion disk, with less reliance on priors than single-epoch microlensing analyses.
In order to analyze the uncorrelated variability of a lensed quasar using the methods of Kochanek (2004), we must first have accurate knowledge of the time delay between the lensed images. In general, imprecise time delays can result in significant uncertainties in microlensing analyses, since residual variability from an improperly corrected time delay must be modeled as uncorrelated microlensing variability in simulated light curves. In the case of SBS 0909, the lens mass model provides very little help with reducing the time delay uncertainty given the very discrepant lens galaxy photometric fits of Lehár et al. (2000) and Sluse et al. (2012). Unfortunately, the time delay of SBS 0909 has been rather difficult to determine because there have been relatively few large-amplitude ( 0.1 mag), short-duration extrema in the light curves in the years since its discovery. The most recent time delay measurement for this quasar, ∆t AB = 49 ± 6 days (Goicoechea et al. 2008), still had relatively significant uncertainty. In our new r-band data set for SBS 0909, the quasar images exhibit substantial intrinsic flux variation, but this flux variation is modulated by the uncorrelated microlensing variability so that a simple polynomial-based cross-correlation analysis (e.g, Kochanek et al. 2006;Poindexter et al. 2007) fails, a problem explored by Eigenbrod et al. (2005). Morgan et al. (2008a) introduced a method of making this problem tractable, by analyzing the microlensing variability with the techniques from Kochanek (2004) while simultaneously solving for the time delay. We will apply this methodology to our r-band data set to make a new, independent measurement of the time delay for SBS 0909 which provides strong evidence in favor of a more compact lens galaxy photometric model (e.g. Sluse et al. 2012) but at the astrometric position favored by Lehár et al. (2000).
In Section 2 we describe our g-and r-band monitoring observations from the United States Naval Observatory -Flagstaff and the Liverpool Telescope and how we com- Note. -HJD is the Heliocentric Julian Day. The magnitudes listed in the QSO A and B columns are measured relative to the comparison stars. The magnitudes in the Stars column are the mean magnitudes of the comparison stars for that epoch relative to their mean over all epochs. The light curve points listed in parentheses have not been included in the analysis. Table 1 will be published in its entirety in the electronic edition of The Astrophysical Journal. A portion is shown here for guidance regarding its form and content.
bined the two data sets. In Section 3 we outline our Bayesian Monte Carlo method for simultaneously estimating time delays and quasar structure. In Section 4 we present the results of our analysis and compare our findings to those of Mediavilla et al. (2011). Throughout our discussion we assume a flat cosmology with Ω M = 0.3, Ω Λ = 0.7, and H 0 = 70 km s −1 Mpc −1 (Hinshaw et al. 2009).

USNO Monitoring
We observed SBS 0909 regularly as part of the United States Naval Academy/United States Naval Observatory (USNA/USNO) Lensed Quasar Monitoring Program. Using the 1.55-m Kaj Strand Astrometric Reflector at the USNO -Flagstaff Station, we take three five-minute exposures of the quasar in the Sloan Digital Sky Survey (SDSS) r-band per epoch at a cadence of two to three nights per month, using either the 2048 × 2048 Tek2K CCD camera (0. ′′ 33 pixel −1 ) or the 2048 × 4096 EEV CCD camera (0. ′′ 18 pixel −1 ). The details of our photometric analysis techniques are discussed in detail in Kochanek et al. (2006). In summary, we measure the quasar image fluxes relative to three reference stars, located at (−12. ′′ 8, 71. ′′ 9), (67. ′′ 2, 11. ′′ 4), and (−13. ′′ 8, 7. ′′ 1) with respect to image A of SBS 0909, using a three-component elliptical Gaussian point-spread function (PSF) model. We hold the relative positions of the quasar images fixed to those derived from the HST H-band images of SBS 0909 for the PSF fitting process. The photometric model of the very red z l = 0.83 lens galaxy is a Gaussian approximation to a de Vaucouleurs profile of fixed effective radius and flux. We use the effective radius derived from the HST images by Lehár et al. (2000) and for the flux we use the value which minimizes the total χ 2 in the residuals over all epochs. We also attempted to measure the quasar photometry using the more compact and dimmer lens galaxy photometric fit from Sluse et al. (2012), but the changes in the quasar image fluxes were negligible. Unlike the case of Q 0957+561 in Hainline et al. (2012), no color offset is required between the two different detectors used for our observing program. In Table 1, we present the rband measurements of SBS 0909 components A and B Note. -HJD is the Heliocentric Julian Day. The magnitudes listed have been calibrated to the SDSS photometric system using the flux of the "b" star in the field of SBS 0909 (see Kochanek et al. 1997). Table  2 will be published in its entirety in the electronic edition of The Astrophysical Journal. A portion is shown here for guidance regarding its form and content.
from 61 nights between 2008 March and 2012 February. The images from which our measurements are derived are characterized by a median stellar FWHM (seeing) of 1. ′′ 3. Because the quasar images are closely spaced (1. ′′ 17), they are blended in our USNO images, causing our photometric analysis to break down for seeing conditions somewhat larger than the image separation, so we keep only epochs for which the seeing is better than 1. ′′ 6 in our analysis. This removes 10 epochs from our USNO data set, which are identified in Table 1 by the parentheses surrounding the measurements. We discarded an additional 6 epochs not listed in Table 1 due to partial cloud cover and bright sky conditions.

Liverpool Telescope Monitoring
We also monitored SBS 0909 in the r-band with the 2.0-m Liverpool Robotic Telescope independently of the observations at USNO. Our r-band monitoring program used the RATCam CCD camera, providing a 4. ′ 6 field of view with pixel scale 0. ′′ 28 pixel −1 , and was carried out over two different periods:  Note. -HJD is the Heliocentric Julian Day. The magnitudes listed have been calibrated to the SDSS photometric system using the flux of the "b" star in the field of SBS 0909 (see Kochanek et al. 1997). Table  3 will be published in its entirety in the electronic edition of The Astrophysical Journal. A portion is shown here for guidance regarding its form and content. In addition, we present here data from the g-band monitoring program of SBS 0909 at the Liverpool Telescope, contemporaneous with the r-band data set published in Goicoechea et al. (2008) and thus spanning two years (2005 January to 2007 January). We used the RAT-Cam instrument at the Liverpool Telescope and obtained 167 individual exposures (frames) of 100 or 200 s each. After bias subtraction, overscan region trimming, and flat fielding of the images, a crowded-field PSF photometry pipeline measures instrumental fluxes for bright field stars and the quasar images. We then transform the photometry to the SDSS photometric system, correcting the instrumental fluxes for the systematic effects of color and inhomogeneity (see Goicoechea et al. 2010). The transformation pipeline is only applied to the frames in which the signal-to-noise ratio (SNR) of the "c" field star, measured through an aperture of radius equal to twice the FWHM, is greater than 100, and for which the seeing is less than 2 ′′ . We also discard frames requiring anomalous color coefficients, frames which produce photometry outliers, and frames in which the quality of the PSF fits to the quasar region are poor. After averaging together the photometry from individual frames obtained on the same night, the final Liverpool Telescope g-band data set consists of 43 epochs of SDSS magnitudes with average uncertainties of 0.016 mag (image A) and 0.036 mag (image B).
The new r-band photometry from the Liverpool Telescope is listed in Table 2 and the complete g-band data set is provided in Table 3. We list the Liverpool Telescope light curves separately from the USNO light curves because the two data sets have different photometric calibrations: the Liverpool data frames have been calibrated onto an absolute system using the absolute flux of a reference star, while the USNO measurements are not calibrated to a standard system. By presenting the two data sets separately, we preserve the original photometric system of each and provide transparency of origin for future users of the data sets.

Construction of Light Curves and Difference Light
Curves In order to construct light curves spanning the longest possible time baseline, we combined all the r-band USNO and Liverpool data for SBS 0909 along with one R-band epoch obtained at the MDM Observatory's Hiltner 2.4m telescope, using the 1024 × 1024 "Templeton" CCD camera, and two SDSS r-band epochs obtained with the WIYN Tip-Tilt Module (WTTM) at the Wisconsin-Indiana-Yale-NOAO (WIYN) 3.5-m telescope. Since these epochs were not contemporaneous with each other or our USNO/Liverpool light curves, we were unable to make an empirical measurement of any magnitude offsets arising from differences between the detectors and filters. We accounted for this unknown offset by applying an additional 0.02 mag of uncertainty to the time-delay corrected flux ratio for these observations. We determined the magnitude offset between our USNO and Liverpool data sets (−14.108±0.018) by making a weighted average of the offset ∆(USNO − LRT) found for the 7 individual nights with contemporaneous observations.
In Figure 1 we show the combined r-band light curves for images A and B from all of our data sources. The dominant feature in the light curves is the intrinsic variability, as much as 0.3 -0.4 mag, the analysis of which we will present in a future paper. Closer inspection, though, of the time range of 4500 HJD − 2450000 5600 reveals uncorrelated variability, as the slope of the increase in the brightness of image B over this time period is noticeably steeper than the slope of image A's light curve. We attribute this component of the variability to microlensing by stars in the lens galaxy. A careful examination of the lightcurves and Tables 1 and 2 will reveal that our data cadence is somewhat variable. In the first few seasons there are many periods with several observations per week, but much of the USNO data is at a cadence of 1-2 observations per month. The strength of our analysis technique is that while the sparsely sampled periods do not constrain the time delay directly, their long time baseline provides very strong constraints on the microlensing model, thereby severely limiting the number of trial light curves from our Monte Carlo simulation with adequate fits to the data during the densely sampled intervals (see Section 3.1).
We display the shorter observed g-band light curves for SBS 0909, which are composed entirely of Liverpool Telescope data, in Figure 2. The intrinsic variability in the g-band is similar to that in the r-band during the same time period. The g-band curves exhibit considerably more scatter than the r-band curves, particularly for image B. The increased scatter is due primarily to the lower quasar flux in g-band and the increased difficulty in extracting bluer g-band fluxes under poor observing conditions. Although the increased scatter impedes visual identification of microlensing variability in the gband light curves, the measurements are still valuable as a complementary data set for a simultaneous two-band microlensing analysis.
The microlensing variability in the g-and r-band light curves is displayed in Figure 3. These time-delay-shifted difference light curves are generated by shifting the light curve of the less variable image (image A, in both bands) Fig. 3.-Difference light curves in g (top panel) and r-band (bottom panel) for SBS 0909, shown with an example of a simulated light curve from our Monte Carlo simulations that is a good fit to the observations. To construct the light curves, image A's data has been shifted by ∆t AB = t A − t B = 50 days. Significant uncorrelated variability is apparent in the r-band.
by the system's 50-day time delay ∆t AB and then performing a linear interpolation of image A's shifted light curve to generate a set of photometric measurements at the same epochs of observation as those in image B's (unshifted) light curve. We discard any data points that were interpolated in the inter-season gaps. Finally, we subtract the light curve of image B from the shifted and interpolated light curve of image A. We do not apply any corrections for lens galaxy reddening to the light curves or difference light curve. Rather, we allow for a 0.5 mag systematic uncertainty in the intrinsic value flux ratios in our simulations, which accounts for the uncertainty in the intrinsic flux ratio from both microlensing and dust extinction in the lens galaxy ). The resulting difference light curves are shown in Figure 3 for both g and r. As expected, the microlensing signal is difficult to identify in the shorter and noisier g-band difference light curve. However, the uncorrelated variability in the r-band is quite clear. We observe a steep fall and rise in ∆m A − ∆m B in the time interval 4000 HJD − 2450000 5200, followed by a second, smaller-amplitude oscillation in the time interval 5500 HJD − 2450000 5900.

ANALYSIS
As discussed in Section 1, we wish to use our combined r-and g-band light curves to derive a more precise time delay for SBS 0909 as well as to place new constraints on the quasar continuum source size and structure. We do so here utilizing the alternate implementation of the Kochanek (2004) microlensing analysis techniques presented in Morgan et al. (2008aMorgan et al. ( , 2012 in which the time delay and source size of a lensed quasar are simultaneously determined through a Bayesian analysis of Monte Carlo microlensing simulations for a range of trial time delays. Our approach here differs only in that we analyze two optical light curves (r-and g-band) for each quasar image rather than an optical light curve and an X-ray light curve for each image. We carry out the microlensing and time delay analysis in two major steps: first, we analyze Monte Carlo simulations of the longer, better-sampled light curves (in this case, r-band) for a range of time delays and source sizes. Then, we simultaneously analyze Monte Carlo simulations of the g and r light curves incorporating the best time delay value and the continuum source size distribution produced in the first stage. The final results for the r-and g-band continuum source sizes are derived in this second step.

Time Delay Analysis of r-band Light Curves
We analyze only the r-band light curves for the time delay analysis, as the r-band monitoring spans a longer length of time and includes the largest-amplitude flux variation. To begin the analysis, we generate r-band light curve pairs in which image A's (the less variable image) light curve is shifted by a set of trial time delays, ∆t AB = t A − t B , according to the procedure outlined in Morgan et al. (2008a). Our set of trial time delays spans the range −70 days ≤ ∆t AB ≤ 70 days in time steps of one day. The shifting of the light curves by the trial time delays must take into account two issues. First, as explained in Morgan et al. (2008aMorgan et al. ( , 2012, each A/B light curve pair must have measurements for the same dates, necessitating interpolation of the shifted image A measurements when they fall in the middle of an observing season and extrapolation of image A measurements when the shifted data points fall in interseason gaps. We used linear interpolation for epochs falling in the middle of an observing season and permitted 10 days of linear extrapolation for points in interseason gaps, increasing the uncertainties for successive extrapolated points as described in Morgan et al. (2012). The second issue is that all trial light curves must have the same number of epochs, requiring the truncation of light curves for trial time delays shorter than 60 days. After all interpolations, extrapolations, and truncations, the shifted light curves used for our analysis each contained 158 epochs.
Next, for each trial time delay, we use the Monte Carlo method of Kochanek (2004, also see Poindexter & Kochanek 2010) to generate large numbers of light curves from microlensing magnification patterns for random combinations of effective velocity v e between quasar, lens galaxy, and observer, mean microlens mass, and lens galaxy stellar mass/dark matter fraction. We then fit the simulated light curves to the observed r-band light curve shifted by the trial time delay. To construct the magnification patterns, we must utilize physical models of the macroscopic (strong) lensing, microlens mass function, and accretion disk surface brightness profile. We used the lensmodel software package (Keeton 2001) to generate a sequence of strong lensing models for the SBS 0909 system over a range of mass contributions from the dark matter and stellar components of the lens galaxy. Each model is a sum of a de Vaucouleurs component, representing the stellar content of the lens galaxy, and a Navarro, Frenk, & White (NFW; 1996) component, representing the galaxy's dark matter halo and concentric with the de Vaucouleurs component. We ran our first realization of the Monte Carlo simulation using a model sequence whose coefficients best reproduce the lens galaxy model of Lehár et al. (2000), but we eventually used a range of model sequences at both the Lehár and Sluse et al. (2012) positions. A model sequence contains ten models spanning 0.1 ≤ f M/L ≤ 1.0 in steps of 0.1, where f M/L represents the ratio of the mass of the stellar component to its mass in a uniform massto-light ratio model. For the stellar (microlens) mass function, we use a power law, dN/dM ∝ M −1.3 , with a ratio of maximum-to-minimum mass ratio of 50; this function reasonably approximates the Galactic disk mass function of Gould (2000). We model the quasar's accretion disk as a face-on, thin disk radiating as a blackbody with a power-law temperature profile T ∝ R −3/4 . Our model matches the outer regions of the thin disk model of Shakura & Sunyaev (1973), but we neglect the drop in temperature in the center due to the inner edge of the disk and the correction factor from general relativity to avoid introducing additional parameters. Provided the disk sizes we obtain are significantly larger than the radius of the inner disk edge, these simplifications introduce insignificant uncertainties relative to those associated with other parameters (Dai et al. 2010). With these parameters we create 40 independent magnification patterns for each quasar image for each of the ten different strong lens models, using the method described in the Appendix of Kochanek (2004). The patterns are 8192×8192 images representing 20 R E ×20 R E , where R E is the Einstein radius for the mean microlens mass M projected into the source plane, yielding a pixel scale of 1.1 × 10 14 ( M /M ⊙ ) 1/2 cm. The outer dimensions and pixel scale are chosen to be sufficiently large to representatively sample the magnification distribution and sufficiently small to adequately resolve the accretion disk in the g-band simulations (see Section 3.2).
For the time delay analysis we carry out 10 7 realizations of the r-band light curve from each of the 400 sets of magnification patterns and for each of the 141 trial time delays. We randomly select an initial position and effective velocity for the source trajectory from their prior distributions under the assumption that these variables are independent and uniformly distributed. For computational simplicity, we neglect the motion of the stars within the lens galaxy and describe the observer's motion as the projection of the CMB dipole velocity onto the lens plane, as done by Kochanek (2004). We compare the simulated light curves to the observed light curves and calculate the goodness-of-fit (χ 2 ) statistics for each, discarding trials with a χ 2 statistic per degree of freedom (χ 2 /ν) greater than 3.
As we discussed in Section 1, the time delay measurement in SBS 0909 is made more challenging by the uncertainties in the lens galaxy photometric model. Since Lehár et al. (2000) and Sluse et al. (2012) disagree about the astrometric position and effective radius (r eff ) of the lens galaxy, we took steps to ensure the robustness of our time delay measurement in the resulting degenerate model parameter space. First, we created a series of additional model sequences for a range of lens galaxy effective radii (r eff ) at the astrometric positions as determined by both Sluse et al. (2012) and Lehár et al. (2000). Each of these model sequences employs a two component (deVaucouleurs -NFW) lens galaxy, where we vary the stellar mass fraction as described above, imposing ellipticity constraints on the lens galaxy from Sluse et al. (2012) in the former case, and a small ellipticity of 1 − b/a = 0 ± 0.08 in the latter case. Then, for each new model sequence, we generated a new set of magnification patterns with which we repeated the fullscale Monte Carlo time delay analysis described above. For completeness, we also repeated this procedure for a Singular Isothermal Ellipsoid (SIE) model at both the Sluse et al. (2012) and the Lehár et al. (2000) positions.
Following completion of the Monte Carlo light curve simulations for all model sequences, we aggregated the results from all trials and performed a Bayesian analysis of the χ 2 statistics of the light curve fits. We formally marginalize over the microlensing variables of source size, microlens mass, and stellar mass fraction used to construct the magnification patterns, as well as the effective velocities from the Monte Carlo simulations, to calculate the posterior probability density for the time delay of SBS 0909. By aggregating the results from the model sequences at both the Lehár et al. (2000) and Sluse et al. (2012) positions, we effectively marginalize over uncertainties in the lens galaxy photometric fits as well. The trials with high χ 2 which were thrown out would not contribute significantly to the Bayesian integrals, so removing them does not affect our results.
We present the posterior probability density for the time delay ∆t AB of SBS 0909 resulting from our Bayesian analysis of the full set of model sequences in Figure 4. For ease of viewing we show only the portion of the distribution for 0 days < ∆t AB < 70 days; the values of the probability density for negative time delays were essentially zero. The time delay distribution is narrowly peaked, with a median ∆t AB = 50 days (so image B leads image A) and a 68% confidence interval of 46 days < ∆t AB < 52 days. Our result is consistent with, but more precise than, the previous time delay measurement ∆t AB = 49 ± 6 days by Goicoechea et al. (2008), based on the first two seasons of r-band monitoring data used in this study. We use this revised time delay for the remainder of our analysis, while acknowledging that we were unable to estimate the influence of all possible systematic errors in our measurement technique. So, as is the case with any sophisticated measurement, it is possible that we underestimated the influence of systematic errors, and these unknown systematics may contribute to the discrepancy we describe in the next paragraph and explore in Section 4.
In Figure 5, we illustrate the influence of the discrepancies in the lens galaxy photometric fits from the literature. We display contour plots of the time delay predicted by the SBS 0909 lens models as a function of the effective radius r eff and fractional mass of the de Vaucouleurs component relative to constant M/L model, f M/L , for lens galaxy astrometry from both Sluse et al. (2012) and Lehár et al. (2000). As expected, the predicted time delay is longer for the more compact lens models with larger stellar mass components, but we also see that the delays are systematically longer for a galaxy located at the Lehár et al. (2000) position than at the Sluse et al. (2012) position. This difference arises because the lens galaxy in the Lehár et al. (2000) fit is closer to image A than in the Sluse et al. (2012) fit, yielding a larger gravitational delay in image A. In any case, neither the fiducial Lehár et al. (2000) model, nor the fiducial Sluse et al. (2012) model can reproduce our new measurement of the time delay ∆t AB = 50 +2 −4 days, but galaxy models at the Lehár position and on the small end of the band of uncertainty in the Lehár effective radius measurement (r eff = 1. ′′ 58 ± 0. ′′ 9) yield delays that are easily consistent with our new time delay measurement. On the other hand, reproducing our new delay measurement with models at the Sluse et al. (2012) position requires a lens galaxy that is signicantly smaller than the already much more compact galaxy in Sluse measurement (r eff = 0. ′′ 54 ± 0. ′′ 02). Fortunately, our Monte Carlo microlensing simulation is sufficiently realistic as to be sensitive to the differences between intrinsic and microlensing variability, despite the uncertainties in the macroscopic lens model, since for each model sequence we sample a wide range of the stellar-to-total convergence ratio κ * /κ, with significant overlap between the macro model sequences.

Simultaneous Dual-Band Microlensing Analysis
Our dual-band (r and g) microlensing analysis follows the method of the simultaneous optical and X-ray analyses by Dai et al. (2010) and Morgan et al. (2012). We first shift the r and g-band light curves by the new time delay, ∆t AB = 50 days, in the same manner as we employed to construct the difference light curves in Section 2.3. Using the same magnification patterns from the time delay analysis, we carry out 10 7 simulations of the r-band light curve for each of the sets of magnification patterns and discard solutions for which χ 2 /ν > 2.5. We saved all the physical parameters from the surviving r-band light curve fits. We then attempted to fit the g-band light curve using the trajectories from the best r-band fits for a new grid of source sizes, and we compute the joint χ 2 /ν for the combined r-and g-band fits. A sample simulated r-and g-band difference light curve which is a best fit to the observed data is shown in Figure 3.
We calculate posterior probability density distributions for the accretion disk sizes in the g-and r-bands and the lens galaxy stellar mass fraction (f M/L ) by performing a Bayesian analysis on the combined set of r-and g-band solutions. We note that our simulations are carried out in Einstein units, where source sizes and velocity are scaled by ( M/M ⊙ ) 1/2 and denoted byr s andv e , respectively. To obtain the probability density for the true, unscaled physical source size P (r s ) from that for the scaled source size, P (r s ), we combine P (r s ) with the probability density for the scaled effective velocity, P (v e ), and a statistical model (i.e., a prior) for the true effective source velocity, P (v e ), in our analysis. We construct P (v e ) using the method described in Kochanek (2004). For that purpose we use the peculiar velocity estimates for the redshifts of SBS 0909 and the lens galaxy from the models presented in  and estimate the velocity dispersion of the lens galaxy from its Einstein radius, assuming the galaxy is a singular isothermal sphere with relaxed dynamics, which Treu & Koopmans (2004) and Bolton et al. (2008) show is a good approximation. As a final step we must correct the scale radius for the disk's inclination i by multiplying by (cos i) −1/2 , which is necessary because we have assumed a face-on disk in our simulations and microlensing amplitudes depend on the projected area of a source rather than the shape.

RESULTS & DISCUSSION
In the top panel of Figure 6 we show the posterior probability distribution for the physical size of the quasar's accretion disk in the observed-frame r-band resulting from our two-band Bayesian microlensing analysis. Hereafter, we state all sizes in terms of the thin disk scale radius, r s , defined as the radius at which the disk temperature matches the rest-frame wavelength of the filter used in our monitoring observations, kT = hc/λ rest (for r-band monitoring of SBS 0909, λ rest = 2620Å; for gband, λ rest = 2020Å). This can be converted to a halflight radius using the relation r 1/2 = 2.44 r s . We note that the distribution shown and all numerical quantities in the discussion which follows have been corrected for disk inclination i through multiplication by a factor of cos i −1/2 , assuming i = 60 • (corresponding to the expectation value of a random distribution of disk inclinations). The median of the probability distribution for r-band source size is log(r s,r /cm) = 15.3 ± 0.3, where the error bar represents the bounds of the 68% confidence interval. In Figure 6 we also show for comparison the 1σ range of values for the accretion disk size obtained by Mediavilla et al. (2011, hereafter M11) from their chromatic microlensing analysis of SBS 0909 under the assumption of a logarithmic size prior and a disk temperature profile power-law index p = 1/α = 4/3. In order to make the most accurate comparison, we converted M11's half-light radius to a thin-disk scale size, scaled the result from its rest-frame wavelength of 1460Å to λ rest (r) = 2620Å assuming R λ ∝ λ 4/3 for thin disks, scaled once more to a mean microlens mass M = 0.3 M ⊙ , and corrected for inclination assuming the same i = 60 • we applied to our disk sizes. As can be seen in Figure 6, our microlensing source size for SBS 0909 is smaller, but marginally consistent with M11's. We suspect that the difference may arise from the effect of the magnification pattern pixel sizes on the size distribution (noting that our r-band disk size is similar to the size of the pixels in M11's magnification patterns), evidence of microlensing of the quasar's broad emission lines, and the different modeling approaches used in the different studies. It is also conceivable that M11's result has been affected by their use of single-epoch spectra and/or their combination of optical and near-IR spectra obtained at epochs separated by several years: the uncorrected time delay and intrinsic variability may alter the continuum and emission line flux ratios from their true values.
Because it is possible that the r-band flux which we observe from SBS 0909 and model in our simulations could be contaminated by UV or optical photons scattered by the broad line region, or higher energy continuum photons reprocessed by the broad line region and re-emitted as emission lines, our r-band accretion disk size may be an overestimate (see Morgan et al. 2010; Guerras et al. -Joint probability density for the ratio of the accretion disk sizes in observed-frame r-band and g-band (rs,r/rs,g) for SBS 0909. The vertical line highlights the location of rs,r/rs,g = 1. The distribution is very wide, reflecting the poor constraints we are able to place on the observed-frame g-band accretion disk size. The median and 1σ values for the size ratio distribution are log rs,r/rs,g = 0.5 +0.9 −1.0 , which are larger but not statistically inconsistent with the r/g-band size ratio expected for a thin accretion disk. 2013). In fact, the prominent Mg II emission line (restframe 2798Å) in the spectrum of SBS 0909 falls within the passband of our r-band filter (see Lubin et al. 2000, M11), so contamination from line emission is of particular concern. To investigate the possibility, we have repeated our microlensing simulations under the assumption that a fraction of the observed r-band flux should actually be attributed to unmicrolensed emission from scattered light or the broad line region. We find, however, that unmicrolensed contamination is not a significant factor in our accretion disk size determination: even when we assume that as much as 30% of the observed r-band flux is contributed by contamination from emission on large physical scales, the median of the r-band size probability distribution is essentially unchanged at log(r s,r /cm) = 15.3 +0.3 −0.4 . We show the posterior probability density for the observed-frame g-band accretion disk size for SBS 0909 resulting from our two-band microlensing analysis in the lower panel of Figure 6. We have included in this panel the 1σ range of the accretion disk size result from M11 as well, corrected as described in the paragraph above, except the result has now been scaled to the rest-frame wavelength of g-band (2020Å) instead of r-band. As might be expected due to the shorter time baseline of the g-band monitoring data and somewhat poorer data quality, the constraints we obtain for the g-band size are not nearly as tight as those for r-band. In light of the wide peak of the g-band probability distribution and the hints of secondary peaks, we regard the g-band result as a preliminary estimate. Despite the large uncertainty, it is encouraging that the median of the g-band physical size distribution, log(r s,g /cm) = 14.8 ± 0.9, is indeed smaller than the median of our r-band size distribution, consistent with the shorter wavelength of g-band and the values for the disk temperature slope in the literature. Like our r-band disk size, our g-band result is notably smaller than the scaled result from M11, although the significance of the discrepancy is low due to the large uncertainties. We also calculate the probability density for the ratio of the r-and g-band disk sizes, which we show in Figure 7. The distribution is very broad, with a median value and 1σ confidence level of log r s,r /r s,g = 0.5 +0.9 −1.0 . Using a standard thin-disk temperature profile (T ∝ R −3/4 ) to predict the observedframe r/g size ratio would result in log r s,r /r s,g = 0.15, which is smaller but statistically consistent with our observed value. Unfortunately, our best value for the r/g size ratio is too uncertain to provide any conclusive indication of the temperature profile of the accretion disk in SBS 0909. We expect that a future analysis utilizing a longer g-band time baseline for photometric monitoring will significantly improve the precision of the g-band size measurement.
An interesting point made by Figure 6 is that the accretion disk sizes we obtain from our two-band microlensing analysis are not consistent with the central supermassive black hole mass derived from SBS 0909's Hβ emission lines by Assef et al. (2011), 10 9.29 M ⊙ . The Hβ line measurement appears to be the most reliable for SBS 0909 in Assef et al. (2011), because the line profiles for C IV and Hα emission lines were difficult to model. By extension, we might expect the black hole mass calculation from Hβ to be the most reliable as well. However, the innermost stable orbit of a maximally rotating Kerr black hole and the Schwarzschild radius predicted for a black hole of mass 10 9.29 M ⊙ both fall within the 1σ bounds of our r-and g-band accretion disk sizes. Another possible clue that the Hβ black hole mass may be problematic is that the theoretical scale radius of a thin accretion disk at 2620Å surrounding a black hole of mass 10 9.29 M ⊙ (the "theory size", log[R 2620 /cm] = 15.73) is larger than the r-band microlensing size from our simulations; yet, in the quasar microlensing literature, the theory size is consistently smaller than the results of microlensing simulations (see, e.g., Morgan et al. 2010;Blackburne et al. 2011). The accretion disk sizes predicted by our microlensing simulations are more consistent with the C IV black hole mass estimate, 10 8.51 M ⊙ , highlighting the difficulties and uncertainties associated with estimating black hole masses from quasar emission lines.
Our Monte Carlo microlensing analysis of SBS 0909 gand r-band monitoring data has enabled us to estimate the r-to-g-band accretion disk size ratio and improve the precision of the system's time delay, despite uncertainties in the macroscopic lens galaxy model. We suggest that deep, high-resolution imaging of the SBS 0909 system will be necessary to completely resolve the lingering questions about the lens galaxy model. For the moment, our delay measurement leads us to favor the lens galaxy astrometry of Lehár et al. (2000) over that of Sluse et al. (2012). Our result for the r/gband size ratio is rather coarse and requires confirmation; to that end, we have begun monitoring SBS 0909 in g-band again. Additionally, we recently expanded the USNA/USNO lensed quasar monitoring campaign to near-infrared (NIR) wavelengths. We look forward to the improvements in our ability to constrain the temperature profiles of quasar accretion disks that will be enabled by size measurements across a significantly larger wavelength baseline.