First robotic monitoring of a lensed quasar: intrinsic variability of SBS 0909+532
To go into the details about the variability of the double quasar SBS 0909+532, we designed a monitoring programme with the 2 m Liverpool Robotic Telescope in the Sloan filter, spanning 1.5 years from 2005 January to 2006 June. The –band light curves of the A and B components, several cross–correlation techniques and a large number of simulations (synthetic light curves) lead to a robust delay = 49 6 days (1 interval) that agrees with our previous results (the B component is leading). Once the time delay and the magnitude offset are known, the magnitude– and time–shifted light curve of image A is subtracted from the light curve of image B. This difference light curve of SBS 0909+532 is consistent with zero, so any possible extrinsic signal must be very weak, i.e., the observed variability in A and B is basically due to observational noise and intrinsic signal. We then make the combined light curve and analyse its statistical properties (structure functions). The structure function of the intrinsic luminosity is fitted to predictions of simple models of two physical scenarios: accretion disc instabilities and nuclear starbursts. Although no simple model is able to accurately reproduce the observed trend, symmetric triangular flares in an accretion disc seems to be the best option to account for it.
keywords:Gravitational lensing, Galaxies: quasars: general, Galaxies: quasars: individual (SBS 0909+532)
, Corresponding author. , , , ,
SBS 0909+532 consists of two components (two–image gravitationally lensed quasar) separated by about 1.1 (Kochanek et al., 1997; Lehár et al., 2000; Lubin et al., 2000; Oscoz et al., 1997a). In optical frames taken at normal seeing conditions with relatively short exposure times, the lensing elliptical galaxy is undetectable (e.g. Ullán et al., 2006). Thus, a simple photometric model (with only two close point–like sources) is able to describe the whole crowded region associated with the quasar components. The first resolved light curves of SBS 0909+532 were presented by Ullán et al. (2006), who derived accurate fluxes in the Johnson–Cousins–Bessel filter from a multisite observing campaign in 2003. These first -band records were used to obtain the time delay between both components. From two different techniques and 1000 repetitions of the experiment (synthetic light curves based on the observed records), Ullán et al. (2006) reported a delay ranging from 41 to 56 days ( 90–95% confidence), where the minus sign means that the intrinsic signal is observed first in the faintest component (B) and later in the brightest component (A). However, the first variability study had some weak points. The –band light curves did not permit to rule out positive delays fairly, and a negative interval [ 90, 0] days was considered in the estimation of uncertainties (component B leading component A). Moreover, there was a relatively poor overlap between the A and B records, when the A light curve was shifted by the best solutions of the time delay and the magnitude offset.
The light curves of the two components of a double quasar at a given wavelength, provide extremely valuable astrophysical information. As the ray paths are different for different components, the corresponding traveltimes will not agree with each other: it appears a time delay between the observed components, which is directly related to the present expansion rate of the Universe (Hubble constant) and the mean surface density of the lensing galaxy (e.g. Kochanek et al., 2004; Refsdal, 1964). To tackle the determination of and , one needs to measure the time delay and the basic parameters of the gravitational mirage, i.e., redshifts and image positions, where the relevant image positions are the positions with respect to the centre of the lensing galaxy (the matter/energy content of the Universe plays also a role). Apart from the time delay, a magnitude offset can be also derived from the comparison of both brightness records. Thus, the possible extrinsic variability (due to intervening objects, e.g., microlenses or dusty clouds in the lensing galaxy) may be detected through the difference light curve (): the difference between the light curves of the components, when the time delay and the magnitude offset are taken into account properly (Schmidt & Wambsganss, 1998). While the s of QSO 0957+561 do not show evidence for extrinsic variability (Gil-Merino et al., 2001; Schmidt & Wambsganss, 1998; Wambsganss & Schmidt, 1998), the s of other lens systems seem to indicate the existence of extrinsic gradients and fluctuations (e.g. Paraficz et al., 2006).
When the difference signal is consistent with zero, the variability observed in both components can be attributed to observational noise and intrinsic phenomena, i.e., physical processes in the source quasar. In this hypothetical case there is a unique opportunity to study the intrinsic signal of a distant quasar (e.g. Kawaguchi et al., 1998). The light curves of the two components can be combined to produce one better–sampled record, and this combined light curve () can be used to carry out statistical analyses. Kawaguchi et al. (1998) analysed the logarithmic slope of the structure function (e.g. Simonetti et al., 1985) of a of QSO 0957+561, which (slope) is in reasonable agreement with the predictions by an accretion disc–instability model. They remarked one important caveat about the methodology to discriminate between different scenarios, which is based on the logarithmic slope of the structure function of the intrinsic luminosity for nuclear starbursts and accretion disc instabilities. The analysis only incorporated two simple models, i.e., a standard starburst model (Aretxaga et al., 1997) leading to relatively high slopes and a cellular–automaton disc–instability model (Mineshige et al., 1994) producing the smallest slopes, so the conclusions could be biased if the simple models do not describe the actual behaviour of the two physical scenarios or the intrinsic signal is caused by several independent mechanisms. However, the comparison with simple analytical or numerical models represents an important first step to reveal the origin of the intrinsic signal.
In order to get high–quality information about the variability of SBS 0909+532 at a red wavelength, we designed a monitoring programme with the world’s largest fully robotic telescope. The 2 m Liverpool Robotic Telescope (Steele et al., 2004) at the Roque de los Muchachos Observatory (Canary Islands, Spain) was used from 2005 January to 2006 June to obtain nightly frames in the Sloan filter. In this paper we present the robotic programme and the corresponding light curves (Sect. 2). Section 3 is devoted to the time delay estimation from the records in the band. Unfortunately, the lens galaxy of SBS 0909+532 has a large effective radius, with a correspondingly low surface brightness (Lehár et al., 2000). These photometric properties complicate the goal of determining an accurate galaxy astrometry, and thus, accurate values of and . Even using the more recent optical frames taken with the Hubble Space Telescope (HST), the position of the centre of the galaxy has a large uncertainty. Consequently, we do not discuss the and values, and focus on a careful analysis of the nature of the observed variability. In Sect. 4, we make the and look for possible extrinsic signal. The and its statistical properties (structure function) are analysed and interpreted in Sect. 5. Finally, in Sect. 6 we summarise our conclusions and put the results in perspective.
2 Observations and light curves
The monitoring programme with the 2 m Liverpool Robotic Telescope (LRT) began in 2005 January. This robotic project was carried out with the RATCam optical CCD camera. The field of view and the pixel scale (binning 22) were and 0.278 arcsec, respectively. Although we concentrate here on the observations in the red arm of the optical spectrum, i.e., the frames in the Sloan passband, the gravitationally lensed quasar was also monitored in the blue arm (via LRT in the Sloan passband). A chromatic () multisite study will be presented further on. Our red subprogramme was optimized to get frames all nights when SBS 0909+532 is visible and there are no technical/atmospheric problems. Each useful night we usually obtained one –band frame of the lens system (exposure time of 120 s). However, many nights in the 2005 October–December period we obtained two 120 s frames, so we are able to check formal photometric errors by comparing them with the artificial (untrue) intranight variabilities. In fact, the intranight scatters are the best non–biased estimators of the typical photometric errors. In Table 1 we include details about the whole –band monitoring from 2005 January to 2006 June (month, frames and observing nights, sampling rate, and total exposure time). Note that this robotic project is quite cheap in observation time, since the total science time is only 18 ks = 5 hours. However, this kind of 1.5–year monitoring (see Table 1) is organizationally very complex or impossible with a conventionally scheduled and operated telescope (Steele et al., 2000).
|Month||Framesnights||Sampling rate (nights/week)||Total exposure time (s)|
|2005 October||11 + 26||1.7||1560|
|2005 November||15 + 211||4.0||3240|
|2005 December||14 + 212||4.0||3360|
|2006 January||11 + 41||0.5||600|
|2006 February||18 + 21||2.2||1200|
|2006 March||116 + 21||4.2||2160|
|2006 April||111 + 21||3.0||1560|
Basic instrumental reductions are applied to all RATCam frames before the data are passed to users. This incorporates bias subtraction, trimming of the overscan regions and flat fielding. After the basic pre–processing, cosmic ray rejection was applied to the frames. There are also bad pixel masks, which are kindly made available by the Angstrom project (Kerins et al., 2006), which is another gravitational lensing programme underway on the LRT. Additional details on the CCD, the Sloan filters and the data pipeline (LRT automatic pre–processing) can be found at http://telescope.livjm.ac.uk/.
A simple photometric model works well for SBS 0909+532 (see Introduction and Ullán et al., 2006). Due to the faintness of the lensing galaxy, the whole crowded region can be described through only two close point–like sources (two components of the lensed quasar). Thus, we apply a PSF fitting method to a representative subset of LRT optical frames verifying some elemental conditions (e.g., the telescope pointing was accurate enough so that the lens system is included in the field of view, there is no a strongly degraded signal, etc). The subset contains 92 frames, so 60% of the robotic exposures in Table 1 are used to make the light curves. In Fig. 1 we show subframes corresponding to three successful exposures: 2005 March 23 (top panels), 2005 December 4 (middle panels) and 2006 March 18 (bottom panels). The double quasar (left panels) appears as an extended structure similar to two tangent stars (see the right panels including the ”c” field star). Following the notation of Ullán et al. (2006), we obtain measurements of , , and , i.e., we compute relative magnitudes using the ”b” star as the reference object (the lens system is inside the triangle defined by the ”a–c” stars, e.g. Kochanek et al., 1997). Both ”a” and ”b” are non–variable field stars, and we expect a constant behaviour of . To check the trend of , we focus on the 69 relative fluxes corresponding to the 2005/2006 season: continuous monitoring from 2005 October to 2006 June. Using the formal photometric errors, the trend is not consistent with a constant because the reduced chi–square value is large ( = 6.27). The formal uncertainties do not seem to include all the sources of error, and we need a non–biased estimator of uncertainties. This is not a problem in our project, since we have several nights with two exposures at different times. Taking into account the deviations in for nights with two data (intranight deviations: , 8 hours, ), it is easy to obtain a standard intranight deviation = 6.8 mmag. Assuming this standard deviation as a typical error (i.e., the intranight variability is not real, but the tool to estimate the true photometric uncertainty) and re–doing the fit to a constant, we derive a very reasonable value of 1.12 (see the filled squares and the discontinuous line after day 3600 in Fig. 2). Therefore, the record in the 2005/2006 season is clearly consistent with a constant behaviour and there is a fair way to infer non–biased uncertainties of the quasar fluxes and (typical errors from standard intranight deviations).
From intranight variabilities in the 2005/2006 season, we derive self–consistent uncertainties: 7 mmag 10.5 mmag 14 mmag 18 mmag. We remark that the ”c” star and the A component have a similar brightness, but A is placed in a crowed region and reasonably is larger than . While the stellar relative fluxes are inferred with true uncertainties less than or equal to 10 mmag, the quasar relative fluxes have true errors in the 10–20 mmag range. Once non–biased errors in and are determined from intranight variations, we group the quasar fluxes each night including more than one exposure. In the 2005/2006 season, the quasar light curves are characterized by a mean sampling rate of two points every week. The LRT was offline in 2005 September (for two weeks) for engineering work which included realuminisation of the primary mirror. Thus, our whole original light curves (incorporating data before 2005 September) are affected by this maintenance work, and we must correct the realuminisation offsets, i.e., the instrumental offsets between 2004/2005 and 2005/2006 seasons. There is no problem with , since we can easily estimate the difference between the average fluxes before and after day 3600. However, the quasar components are variable objects, and we need independent brightness records to correct the offsets in both and . A simultaneous monitoring in the band at Mt. Maidanak, Uzbekistan, is used to derive the realuminisation biases in the LRT records. Maidanak (1.5 m AZT–22 Telescope) quasar fluxes are compared with LRT brightnesses at separations less than or equal to 3 days, which leads to instrumental offsets of about 70 and 10 mmag in and , respectively. As the and colours of the ”b” star are close to the colours of the B component (similar spectra), it is expected a very small correction in the differential curve . For this record , we just obtain an offset (slightly less than 10 mmag) that agrees with theoretical predictions. The whole final brightness records in the band are depicted in Fig. 2. In this figure, we plot the A light curve (filled circles), shifted by + 0.50 mag, and the B light curve (open circles). To compare quasar fluxes against stellar brightnesses, relative fluxes of the ”a” star (filled squares), shifted by + 3.13 mag, are also shown in Fig. 2. The LRT brightnesses in the -band (78 data of each component) have errors above 10 mmag, which could complicate some analyses. Thus, to improve the situation we also group fluxes within 3–day intervals. The new data set only contains 41 fluxes of each component, but both mean errors are reduced to the 10–mmag level: 10 mmag and 13 mmag. In Fig. 3 we show the new grouped light curves, i.e., the A curve (filled circles) and the B curve (open circles). Now, to easily compare both trends, the A record is properly shifted in time and magnitude ( 49 days and + 0.65 mag; see sections 3–4). The individual calibrated quasar fluxes are presented in Table 2. We use the SDSS flux of the ”b” star in the band ( = 14.87 mag) to get apparent magnitudes of both components111The SDSS Web site is http://www.sdss.org/..
3 Confirmation of the time delay
If we concentrate our attention in Fig. 3, the A and B light curves show a decline larger than 50 mmag as well as a 50–mmag event around day 3800 (A is shifted in time and magnitude). Thus, although there is an important gap in the LRT monitoring (due to occultation of the lens system in 2005 July–September), these two features are promising tools to confirm our previous time delay estimation (see Introduction). To calculate the time delay between both components of SBS 0909+532, we use different cross–correlation techniques. First, we focus on the grouped light curves and the discrete cross–correlation function (). The was introduced by Edelson & Krolik (1988) and has been extensively used in delay studies (e.g. Gil-Merino et al., 2002; Ofek & Maoz, 2003; Oscoz et al., 1997b). Here, the of SBS 0909+532 is evaluated every 1 day in the region from 200 to + 100 days, so we analyse a wide range of lags including both positive and negative values. The is binned in 2 day intervals centered at the lags. To work with a reasonable time resolution, we only take into account values less than or equal to 10 days, i.e., bins with width 20 days. For 4 days, the is very noisy, whereas for = 8–10 days, the main peak of the is significantly reduced with respect to the expected value of 1. The more interesting results are derived from = 5–7 days, and we choose the intermediate bin ( = 6 days) as the most suitable one. For = 6 days, there is a maximum at 50 days ( 0.9). Apart from the main peak around the maximum (delay–peak), there are other secondary peaks at negative and positive lags. These secondary structures have an amplitude of about 0.5, i.e., they are clearly smaller than the main feature. No correlation ( 0) or anticorrelation ( 0) is also found at the edges of the time lag–interval. For other values, we also find maxima at or around 50 days.
As it was discussed by Lehár et al. (1992), when the main fluctuations in the light curves have an intrinsic origin, the irregular delay–peak of the AB cross–correlation function should be closely traced by the symmetrical central peak (around a lag equal to zero) of the AA (or BB) autocorrelation function. Moreover, other features of the cross–correlation function around lags , ,… will be closely reproduced in the autocorrelation function around lags , ,…, respectively (with relation to A, we assume that B is delayed in ). Therefore, if the shifted discrete autocorrelation function () is matched to the discrete cross–correlation function (), one derives the time delay in a self–consistent way (note that the lag corresponding to the maximum of the is only a rough estimation of the delay). This self–consistent methodology is called the technique, and it has been successfully applied to some golden data sets (see next paragraph) for QSO 0957+561 (e.g. Goicoechea et al., 1998; Serra-Ricart et al., 1999).
Before applying the technique, we need to make a golden data set. A golden data set contains a variable and free from long gaps light curve A in a period [, ], and a variable and free from long gaps light curve B in a period [, ], where . As we know a rough initial estimation of the SBS 0909+532 delay (through the maximum of the , see here above), we can use this value ( = 50 days) to correct the unsuitable 50–day edges. It is also clear that the records in Fig. 3 are variable. Therefore, the existence of long gaps is the only difficulty to make a golden data set for SBS 0909+532. Can we relax the condition on the absence of gaps?. Simulations by Serra-Ricart et al. (1999) indicated that the method works even in the presence of relatively long gaps. Moreover, the main problem with gaps is the possible bias between the differences and , which are the key pieces in the cross–correlation. In the periods of interest, the two light curves are sampled in different ways, so the direct averages of the fluxes and could lead to a bias between the differences and . Thus we can obtain a quasi–golden data set provided that and are carefully determined. Instead of the direct averages of the grouped fluxes and in the periods of interest, the mean values and are inferred from a method that corrects for sampling bias.
Once we have a quasi–golden data set, a comparison between the corresponding (filled circles) and (open circles) is plotted in Fig. 4. The is the average of the AA and BB autocorrelation functions, and it is shifted by 50 days, i.e., the rough estimation of the delay. There is no important distortions in the features of the as compared with the features in the , which strengths the realibility of a delay close to 50 days. Using = 6 days, possible values of the time delay () versus the associated values normalised by its minimum value , are also plotted in Fig. 5. The function is defined in Eq. (1) and Eq. (7) of Goicoechea et al. (1998) and Serra-Ricart et al. (1999), respectively. In Fig. 5, a relatively narrow peak centered on 49 days (best value of the delay) is derived. To obtain uncertainties, we follow an approach similar to that described by Ullán et al. (2006). We make 1000 repetitions of the experiment, apply the minimization ( = 6 days) to each quasi–golden synthetic data set, and thus obtain 1000 best values of the delay. Through the distribution of delays, our 1 measurement is = 49 7 days (69.8% confidence interval). This result confirms the previous delay determination from Calar Alto and Maidanak frames in 2003. Moreover, the new measurement is very robust, since it is not a pre–conditioned estimation (a wide range from 200 to + 100 days is tested) and there is a significant overlap between the A and B records, when the A light curve is shifted by the best solution of the time delay (see Fig. 3).
We also apply the modified cross–correlation function () technique (Beskin & Oknyanskij, 1995; Oknyanskij, 1997). The combines properties of both standard cross–correlation functions: the interpolated cross-correlation function by Gaskell & Sparke (1986) and the by Edelson & Krolik (1988). This time we analyse the individual (non–grouped) fluxes in Fig. 2, and after doing some tests, we use 10–day bins in the component A. When the is applied to our data in the lag interval [ 200, + 200] days, the maximum correlation coefficient ( = 0.878) corresponds to a lag = 48 days. The correlation coefficient at different lags (days) is shown in Fig. 6. We also check if the removal of some fluxes from the B record can influence the estimation of (, ). The stability test is done in a blind way, i.e., we do not choose by eye the points to be dropped, but a systematic procedure is used to select those points. In a first step, the first 3 points are removed from the B light curve, then the second, third and fourth points are removed from, etc. From this 3–point systematic cleaning, almost all the values are close to 0.9 and 85% of the iterations lead to = 49 days, in good agreement with the analysis of the whole data set. In order to estimate delay errors, we carry out simulations taking into account all the properties of the observed records: kind of variability, sampling and photometric errors (Koptelova et al., 2006). About 1000 synthetic data sets leading to 0.8 are considered here. Our 1 measurement is = 49 5 days (70.6% confidence interval), which is practically identical to the 1 result from the method. We consider that the two techniques ( and ) have similar quality. Thus, there is no a fair way to choose either 5 or 7 days as error, and we adopt 49 6 days as the final 1 estimation, where 1 day is the uncertainty in the error.
4 Difference light curve
From now on we take the LRT grouped fluxes in Fig. 3 as basic tools for discussing variability properties in the filter. First, considering the time delay in section 3, we infer the time delay–shifted light curve of A (). Here, is the record shifted by = 49 days in the horizontal direction. Then a mean offset is computed in a direct way. This magnitude offset is close to + 0.65 mag. Second, the magnitude– and time–shifted light curve of image A (, where is the curve shifted by the magnitude offset in the vertical direction) is subtracted from the light curve of image B (). Both and are plotted in Fig. 3. We thus obtain the difference light curve () . To compute the mean offset as well as the , the dates in the time shifted curves and are taken as reference epochs. The and fluxes are then compared to the averaged values of in bins with semiwidth centred on the reference dates.
The of SBS 0909+532 in the years 2005–2006 appears in Fig. 7. We use two different values: 6 days (filled circles) and 15 days (open circles), and both trends are consistent with each other. The differences (in –magnitudes) do not exceed the 0.05 mag thresholds (discontinuous lines). Moreover, there is no evidence in favor of the existence of events or gradients. The difference signal is in apparent agreement with zero, i.e., Fig. 7 shows a noisy relationship = . From a quantitative point of view, we can also estimate the goodness of representing the data with = 0 (e.g. Schmidt & Wambsganss, 1998). The reduced chi–square values are 0.92 ( = 6 days) and 0.87 ( = 15 days), corroborating the good agreement between the light curves of the two images. Therefore, the observed variability in A and B is basically due to observational noise and intrinsic signal, although we cannot rule out the existence of a very weak extrinsic signal whose amplitude must be well below the noise level in the . Our constraints on the possible microlensing (extrinsic) variability can be used to obtain information on the granularity of the matter in the lensing galaxy and the size of the source (Gil–Merino et al., in preparation). Although Mediavilla et al. (2005) suggested the existence of an achromatic microlensing magnification of image B, a homogeneous microlensing pattern will not produce microlensing variations but a microlensing offset as part of the magnitude offset. Only an inhomogeneous microlensing pattern can be detected through the .
5 Combined light curve: origin of the intrinsic signal
The presented in Sect. 4 is in clear agreement with the absence or an extremely low (undetectable) level of extrinsic signal. Therefore, in this section, we make the combined light curve () and interprete it as due to the observational noise and intrinsic phenomena. The combined photometry consists of both light curves and , and this is plotted in Fig. 3 (filled and open circles). The data points in Fig. 3 show the global record (in –magnitudes) in a 500–day interval, with the data covering a period 450 days (there is an unavoidable gap of about 50 days that is related to the annual occultation of the quasar).
The is generated by intrinsic signal () and observational noise (), so . This simple relationship between the combined fluxes, the underlaying intrinsic signal and the noise permits to estimate the first–order structure function of (e.g. Gil-Merino et al., 2001; Simonetti et al., 1985). A structure function analysis is a method of quantifying typical flux variabilities at different lags. The structure function at lag is given by
where the sum only includes the (,) pairs verifying that (the number of such pairs is ). We take a normalization factor equal to 1/2, and thus, the asymptotic behaviour on long timescales is just the signal variance instead of 2 (e.g. Collier & Peterson, 2001). Fig. 8 shows the structure function of SBS 0909+532 obtained from the and Eq. (1), using independent 10–day bins (filled circles). To test the robustness of the derived , we also obtain extreme structure functions, i.e., using the extreme values of the delay range, = 43 days and = 55 days, to make the corresponding s. These two extreme trends (continuous lines in Fig. 8) are included in the error bars of the , so the typical structure function (from = 49 days) is a reliable tool. The discontinuous line in Fig. 8 corresponds to a 10 mmag threshold, and at lags 60 days, the typical fluctuations are below this 10 mmag level. This result explains why the flux ratio at the same time of observation basically coincides with the flux ratio corrected by the time delay, at least in the red arm of the optical spectrum (Ullán et al., 2006).
To discuss the origin of the intrinsic signal, we focus on the structure function of the intrinsic luminosity (e.g. Cid Fernandes et al., 2000). This can be directly compared to the predictions of different physical scenarios, e.g., nuclear starbursts and accretion disc instabilities (Kawaguchi et al., 1998). We remark that the observed flux at is emitted at a shorter wavelength , and we must consider the emission of near UV light ( 2600 Å), because the observations were made in the Sloan band ( 6200 Å). On the other hand, = = 2.5 log(/), where and denote –band magnitudes and monochromatic fluxes, respectively. We can also use the cosmological law , with , and being respectively monochromatic luminosity of the quasar, extinction–magnification factor and luminosity distance. From some rearrangement, it is inferred the relationship = 10, = . We initially take a set of units so that = 1 and = = 10, so the structure function at rest–frame lag = can be estimated through the averaged sum
where and the sum includes pairs verifying . As the values from Eq. (2) are 10, we decide to use more appropriate units: = 1/1000. In these new units, the structure function of the intrinsic luminosity is drawn in Fig. 9 (filled and open circles). In this figure, the dotted straight line represents the asymptotic behaviour , i.e., at long lags. From now on, we only consider the data points not exceeding the flat asymptotic behaviour (filled circles). There is an arrow to make mark on the last reliable point. This last valid point corresponds to 150 days /3. Moreover, at 150 days, the number of pairs in each bin is 100–200. Longer lags cannot be seriously considered in the analysis, since the associated information is probably biased and could not describe the typical behaviour of the signal.
The reliable trend in Fig. 9 (filled points) should be related to the variability scenario and may unveil the origin of the intrinsic fluctuations. First, we use three simple analytical models of accretion disc instabilities. In these three Poissonian models, the luminosity is due to the superposition of a variable component and a constant background (the nonflaring part of the accretion disc). The variable component is made by the superpositions of flares at random times, which are characterized by a typical timescale . Square (SQ), exponentially decaying (ED) and symmetric triangular (ST) flares are taken into account (see Appendix B of Cid Fernandes et al., 2000). When the observed is fitted to the analytical intrinsic structure functions, the SQ and ED models do not work well and lead to 4–5. However, the ST model leads to a more reasonable (but not sufficiently good) value of 2 (solid line in Fig. 9). For the best ST model, = 70 days and the lifetime is = 3/2 100 days. There are several physical timescales that might be associated with accretion disc instabilities. The main dynamical timescales are the free–fall time and the orbital time . The thermal timescale is given by , i.e., the time for vertical diffusion of heat. Apart from the black hole mass and the emission radius , the thermal timescale also depends on the disc viscosity parameter (this last parameter is usually named , but we rename it to avoid confusion). The viscous time ( is the disc thickness) is also relevant. This is linked to the radial inflow of a gas element, and represents the longest timescale (e.g. Czerny, 2004; Krolik, 1999, and references therein). Collier & Peterson (2001) estimated the values of , and for different values of at two rest–frame wavelengths: 1400 and 5000 Å (see Fig. 5 in that work). To do the estimations, they assumed a standard geometrically thin and optically thick disc Shakura & Sunyaev (1973). At 2600 Å (see above), the lifetime of ST flares is consistent with the orbital time for 10 M black holes and the thermal timescale for 10 M black holes.
Many previous studies focused on the logarithmic slope () of the growing part of the square root of structure function, which (slope) is directly related to the physical mechanism responsible for the variability (e.g. Kawaguchi et al., 1998). For example, for SQ flares of duration , at lags ( is flat at ). This leads to a constant slope = 0.5. For ED and ST flares, reaches its maximum value at the shortest lags. If the ED flares are characterized by a semi–lifetime , then 0.5 at . On the other hand, for ST flares with rise time , is approximately proportional to at . Thus, whereas should be less than or equal to 0.5 for the SQ and ED models, the ST model produces a steep slope 1 at relatively short lags. What about the observed structure function?. The observed slope is 1 using both the subset of data at 35 days (the first six points in Fig. 9) and the global data set until = 63 days. It is clear that the observed slope of about 1 favours the ST model of accretion disc instabilities. However, even ST flares with equal rise and decay times = 70 days (best solution) have not the ability to reproduce the observed slope at lags exceeding 50 days (see Fig. 9). In fact this is the reason to obtain a rough agreement ( = 2) instead of an accurate fit leading to 1. The cellular–automaton model (Mineshige et al., 1994) is another framework to describe accretion disc instabilities. In spite of its popularity, we cannot consider that physical mechanism as an alternative to explain the observations, since it generates a small slope: = 0.4–0.5 (Kawaguchi et al., 1998).
With respect to the nuclear starbursts, we use the standard model by Aretxaga et al. (1997). This simple model is characterized by a timescale , which is the time when the supernova remnant reaches the maximum of its radiative phase. The observed cannot be reproduced by the relationship in Eq. (17) of Aretxaga et al. (1997). Although the best timescale has an acceptable value of = 60 days (considering that the source is a very bright quasar), the observations–model comparison leads to a poor fit with = 4 (dashed line in Fig. 9). Standard nuclear starbursts are able to induce a steep slope at lags shorter than . However, the flat regime is only reached at lags significantly longer than (Aretxaga et al., 1997; Kawaguchi et al., 1998). This mechanism is thus unable to reproduce a steep slope from short lags to the flattening lags, and consequently it fails to trace the observed behaviour.
6 Conclusions and discussion
The 2 m Liverpool Robotic Telescope (Steele, 2001; Steele et al., 2004) at the Roque de los Muchachos Observatory (Canary Islands, Spain) is ideally suited to monitorize gravitationally lensed quasars (GLQs) and derive light curves of their components. In this paper we present the first GLQ monitoring campaign using the Liverpool Telescope. We have observed the double quasar SBS 0909+532 in 2005–2006, taking nightly frames in the Sloan filter. Our main results and conclusions are:
The –band frames and PSF photometry lead to variable light curves of the two quasar components A and B. Both brightness records show a decline larger than 50 mmag followed by a 50–mmag event. These features are used to confirm the recently reported time delay between the components (Ullán et al., 2006). From different cross–correlation techniques and a large number of repetitions of the experiment (synthetic light curves based on the observed records), we infer a delay = 49 6 days (1 interval), which agrees with our previous results. The new delay determination is robust, since a wide range of possible delays is tested and the two light curves overlap over a long time interval, when the A record is shifted by the best solution of the time delay.
To obtain the difference light curve of SBS 0909+532, the magnitude– and time–shifted light curve of image A is subtracted from the light curve of image B. There is no evidence in favor of the existence of extrinsic variability (e.g., fluctuations or gradients caused by microlensing in the lensing galaxy halo), since the difference curve is consistent with zero. Very recently, Paraficz et al. (2006) derived difference curves of five GLQs and found that two out of the five GLQs have significant extrinsic signal. Paraficz et al. (2006) also presented one difference curve that clearly agrees with zero (HE 2149–2745) as well as two doubtful cases. Here we report on another GLQ with flat difference curve, which joins the family of systems without important extrinsic variations (e.g. Burud et al., 2002; Gil-Merino et al., 2001; Schmidt & Wambsganss, 1998).
Taking into account the absence or the extremely low (undetectable) level of extrinsic signal, we make the combined light curve and interprete it as due to the observational noise and intrinsic variations. The combined photometry consists of both A and B records, where the A light curve is shifted by the best solutions of the time delay and the magnitude offset. In order to study the growth of intrinsic variability with rest–frame lag, we explicitly obtain the structure function of the intrinsic luminosity. This is then fitted to predictions of simple models of accretion disc instabilities and nuclear starbursts. With respect to the disc–instability scenario, we concentrate on the phenomenological models by Cid Fernandes et al. (2000), since the cellular–automaton model (Mineshige et al., 1994) is related to a relatively small logarithmic slope of the structure function (Kawaguchi et al., 1998). This small slope ( = 0.4–0.5) is not consistent with the structure function derived from observations ( 1). We also consider the standard nuclear starbursts by Aretxaga et al. (1997). Symmetric triangular flares in an accretion disc lead to the best (but not accurate enough) fit, whereas standard SN explosions cannot produce a large value of from short to long lags, and thus, cannot account for the observed variability.
We obtain a rough agreement between the observations of SBS 0909+532 and the production of 100–day symmetric flares in the accretion disc of the distant source quasar. Considering typical supermassive black holes, the lifetime of the flares seems to agree with the dynamical and thermal times at the emission radius corresponding to a standard disc (Collier & Peterson, 2001). This suggest the existence of disc instabilities having the local (emission ring) dynamical–thermal timescale. Local magnetorotational instabilities (e.g. Balbus & Hawley, 1991) would be possibly able to cause variability over this timescale. However, it is difficult to identify precisely the physical origin of the fluctuations, and we cannot rule out other processes. For example, the hottest temperature at the innermost radius of a standard gas disc is clearly smaller than X–ray energies. Therefore, as SBS 0909+532 is a bright X–ray source (Chartas, 2000; Page et al., 2004), inverse Comptonization in a hot inner corona may explain its X–ray emission. Instabilities in the corona, which is presumably unstable (e.g. Shakura & Sunyaev, 1976), could lead to variable X–ray emission as well as variable X–ray irradiation of the gas disc. After X–ray reprocessing in the disc, variability at optical/UV wavelengths is expected (e.g. Czerny, 2004). Manmoto et al. (1996) also proposed an unstable advection–dominated disc that might generate X–ray flares and subsequent optical/UV events. In both physical schemes, the shape and timescale of the optical/UV variations are determined by the shape and timescale of the X–ray fluctuations. Are these reasonable pictures?. X–ray symmetric flares were previously found by Negoro et al. (1994), but on a stellar scale. Cygnus X–1 is an accreting black hole system ( 10 M) that produces symmetric flares lasting about 1 second (see Fig. 1 of Negoro et al., 1994). The Cygnus X–1 variability was reproduced through theoretical calculations by Manmoto et al. (1996), so it may be associated with advection–dominated disc instabilities. If the same mechanism is responsible for the X–ray variability of both stellar and supermassive black hole systems, then the stellar timescale must be rescaled by 10–10 (assuming typical supermassive black holes with 10–10 M) to obtain 10–1000 days symmetric flares on a quasar scale. For 10 M, these hypothetical X–ray events would be good tracers of our near UV fluctuations. Thus, surprisingly, we find that advection–dominated disc instabilities could be able to account for variations on both stellar and quasar scales. An unstable hot corona is also a plausible source of X–ray flares, but a discussion on the involved timescale and (a)symmetry is out of the scope of this paper.
To date only one more GLQ with flat difference curve was analysed in some detail to determine the origin of its intrinsic variability. From light curves of QSO 0957+561 in two different optical bands ( filters), Kundić et al. (1997) inferred variability functions and their corresponding logarithmic slope 0.4. Here, are quasar magnitudes, i.e., the observational noise was not subtracted from the photometric measurements. Using only data in the band, Kawaguchi et al. (1998) also determined a shallow slope of : 0.3–0.4. Kawaguchi et al. (1998) presented a growing trend at lags 500 days ( 200 days), but they did not report on the asymptotic behaviour ( for the Kawaguchi et al.’s normalization) or the reliability of the results at long lags, so we cannot properly discuss the involved variability timescale. On the other hand, the observed rise in the variability is consistent with a cellular–automaton disc–instability model (Kawaguchi et al., 1998). However, in spite of the small errors in the light curves ( 10 mmag), the noise in the structure function ( 15 mmag) does have a significant effect on the measured variations at the shortest lags and the measured sope. Thus, the slope of should be steeper than the slope of . If we roughly construct a noise–less function in the band, the new slope is 0.7. Moreover, this new (and probably true) slope agrees with the slope in the band (Gil-Merino et al., 2001). Therefore, there are clear evidences that both QSO 0957+561 and SBS 0909+532 are = 1.4 bright quasars with 0.7–1.
We may also try to establish a relation between these GLQ results and variability studies of non–lensed quasars. While some works focused on ensemble structure functions of large samples (hundreds, thousands or even more) of quasars (e.g. de Vries et al., 2005; Hawkins, 2002; Vanden Berk et al., 2004), other works cocentrated on well–sampled individual quasars (e.g. Cid Fernandes et al., 2000; Collier & Peterson, 2001). The first ones did not include details on individual objects. Moreover, additional problems usually complicate a direct comparison with current GLQ data. For example, de Vries et al. (2005) stated that their study is insensitive (due to the measurement noise) to lags shorter than 1 year. Thus, we cannot put the behaviour of GLQs in perspective. Hawkins (2002) mostly used relatively faint objects showing coherent fluctuations of about 1 mag over a few years (see Fig. 4 in that work). However, both GLQs are relatively bright sources displaying coherent flux variations of about 0.1 mag over a few or several months. On timescales of a few years, the GLQs have a maximum scatter (difference between maximum and minimum flux) below 0.25 mag. Moreover, instead of a noise–less structure function having a good time resolution at lags 1–2 years, Hawkins (2002) studied a long term ensemble with time resolution of about one year. The most relevant ensemble variability was reported by Vanden Berk et al. (2004), who showed an vs. relationship having reasonable time coverage and resolution. This ensemble variability for quasars with very different luminosities (absolute magnitudes) and redshifts led to a slope 0.3 at 500 days. It is evident that the GLQ slope 0.5 disagrees with the ensemble slope ( 0.5), which could indicate the existence of different populations of intrinsically variable quasars (perhaps depending on luminosity and redshift) or even the existence of microlensing effects on a large amount of quasars (e.g. Hawkins, 2002). With respect to the well–sampled individual sources, Cid Fernandes et al. (2000) analysed the light curves of 42 0.4 quasars that were monitored by the Wise Observatory group (Giveon et al., 1999). Assuming the presence of square flares ( = 0.5), they derived lifetimes ranging from 124 days to a few years, so the shortest ones are close to the lifetime of the SBS 0909+532 flares. Cid Fernandes et al. (2000) also remarked the insensitivity of their fits to the flare shape (with the exception of exponentially decaying flares, since these asymmetric flares led to poor fits). These results (relatively long lifetimes and poor sensitivity to the flare shape) could be due, at least in part, to the procedure for fitting observations, which involved two free parameters, i.e., the flare lifetime and the variance. We, however, compute directly the variance from the light curve and then fit the lifetime. Our variance represents the typical variability in a period of about 200 days (rest–frame time), and taking into account the absence of important gradients over 1.5 years (from a comparison with previous Calar Alto light curves Ullán et al., 2006), this seems a good tracer of the variance over longer periods. Collier & Peterson (2001) studied 13 local active galactic nuclei with very good time coverage and resolution. The slopes of their optical/UV structure functions varied between 0.3 and 0.8 (similar to the GLQ slope). The flare lifetimes (using symmetric triangular flares and certain lag intervals) are 5–94 days, with the longest ones very close to the lifetime of the SBS 0909+532 events. In order to check the results based on a monitoring period of 200–300 days, Collier & Peterson (2001) presented the variability structure for 7 years of monitoring of NGC 5548. Using the short monitoring period, they found 40–60 days. However, from the additional analysis of the longer monitoring period, they inferred richer results. The asymptotic behaviour is reached at lags of 200 days, and the 50–day timescale seems to represent the lifetime of the shortest flares.
To accurately reproduce the observed variability of SBS 0909+532, new (more complex) models seem to be required. Moreover, a longer monitoring period would lead to a very precise description of the structure at 60 days and the appearance of reliable features at lags exceeding 60–70 days, which would significantly improve the observational constraints. We also note that the current optical data of SBS 0909+532 are insufficient to carry out more ambitious analyses, e.g., a test for non–linearity. Uttley et al. (2005) showed that a linear variability–flux relation in the light curve of an object is an indicator of non–linearity. This hypothetical linear relation would suggest that the variability process is multiplicative, or in other words, it is the result of several random subprocesses which multiply together. A process produced by multiplication of many independent processes has a lognormal distribution, so the hypothetical light curve would have a lognormal probability density function (PDF) (for details, e.g. Uttley et al., 2005). However, the observed fluctuations in the flux of SBS 0909+532 ( 10%) and the number of available data points (clearly insufficient to make highly populated bins in flux) do not permit to reliably obtain the variability–flux relation or the PDF of SBS 0909+532. Finally, we remark the unique advantages of using GLQs as a tool to study the origin of the intrinsic signal of quasars, because there is no way to disentangle intrinsic from extrinsic signal in a non–lensed quasar. New long term monitoring programmes of GLQs with flat difference curve will permit to fairly discuss the structure of the intrinsic variability on short and long timescales, and to compare between this behaviour and the structure of the fluctuations in both GLQs with non–flat difference curve and non–lensed quasars.
We thank an anonymous referee for several comments which significantly improved the discussion concerning the intrinsic variability. We are also indebted to C. Moss for guidance in the preparation of the robotic monitoring programmes. To check the monitoring campaign in Spain, we use observations at Mt. Maidanak (Uzbekistan). We acknowledge the observers and heads of the Maidanak collaboration for their maintenance and availability of the gravitational lenses database. This research has been supported by the Spanish Department of Education and Science grant AYA2004-08243-C03-02, University of Cantabria funds, grant for young scientists of the President of the Russian Federation (number MK-2637.2006.2), Deutscher Akademischer Austausch Dienst (DAAD) grant number A/05/56557, grant of Russian Foundation for Basic Research (RFBR) 06-02-16857 and the Science and Technology Center of Ukraine (STCU) grant U127k. We also acknowledge support by the European Community’s Sixth Framework Marie Curie Research Training Network Programme, Contract No.MRTN-CT-2004-505183 ”ANGLES”. Funding for the Sloan Digital Sky Survey (SDSS) has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the NASA, the NSF, the U.S. Department of Energy, the Japanese Monbukagakusho, and the Max Planck Society. The SDSS is managed by the Astrophysical Research Consortium (ARC) for the Participating Institutions. The Participating Institutions are The University of Chicago, Fermilab, the Institute for Advanced Study, the Japan Participation Group, The Johns Hopkins University, Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, University of Pittsburgh, Princeton University, the United States Naval Observatory, and the University of Washington.
- Aretxaga et al. (1997) Aretxaga, I., Cid Fernandes, R., Terlevich, R., 1997. MNRAS 286, 271.
- Balbus & Hawley (1991) Balbus, S. A., Hawley, J. F., 1991. ApJ 376, 214.
- Beskin & Oknyanskij (1995) Beskin, G. M., Oknyanskij, V. L., 1995. A&A 304, 341.
- Burud et al. (2002) Burud, I., et al., 2002. A&A 383, 71.
- Chartas (2000) Chartas, G., 2000. ApJ 531, 81.
- Cid Fernandes et al. (2000) Cid Fernandes, R., Sodré, Jr. L., Vieira da Silva, Jr. L., 2000. ApJ 544, 123.
- Collier & Peterson (2001) Collier, S., Peterson, B. M., 2001. ApJ 555, 775.
- Czerny (2004) Czerny, B., 2004. The role of an accretion disk in AGN variability, astro-ph/0409254.
- de Vries et al. (2005) de Vries, W. H., Becker, R. H., White, R. L., Loomis, C., 2005. AJ 129, 615.
- Edelson & Krolik (1988) Edelson, R. A., Krolik, J. H., 1988. ApJ 333, 646.
- Gaskell & Sparke (1986) Gaskell, C. M., Sparke, L. S., 1986. ApJ 305, 175.
- Gil-Merino et al. (2001) Gil-Merino, R., et al., 2001. MNRAS 322, 397.
- Gil-Merino et al. (2002) Gil-Merino, R., Wisotzki, L., Wambsganss, J., 2002. A&A 381, 428.
- Giveon et al. (1999) Giveon, U., Maoz, D., Kaspi, S., Netzer, H., Smith, P. S., 1999. MNRAS 306, 637.
- Goicoechea et al. (1998) Goicoechea, L. J., Mediavilla, E., Oscoz, A., Serra-Ricart, M., Buitrago, J., 1998. Ap&SS 261, 341.
- Hawkins (2002) Hawkins, M. R. S., 2002. MNRAS 329, 76.
- Kawaguchi et al. (1998) Kawaguchi, T., Mineshige, S., Umemura, M., Turner, E. L., 1998. ApJ 504, 671.
- Kerins et al. (2006) Kerins, E., et al., 2006. MNRAS 365, 1099.
- Kochanek et al. (1997) Kochanek, C. S., et al., 1997. ApJ 479, 678.
- Kochanek et al. (2004) Kochanek, C. S., Schneider, P., Wambsganss, J., 2004. In: Meylan, G., Jetzer, P., North, P. (Eds.), Proc. 33rd Saas-Fee Advanced Course, Part 2 of Gravitational Lensing: Strong, Weak & Micro, Springer-Verlag, Berlin.
- Koptelova et al. (2006) Koptelova, E. A., Oknyanskij, V. L., Shimanovskaya, E. V., 2006. A&A 452, 37.
- Krolik (1999) Krolik, J. H., 1999. In: Poutanen, J., Svensson, R. (Eds.), ASP Conf. Ser. Vol. 161, High Energy Processes in Accreting Black Holes, Astron. Soc. Pac., San Francisco, p. 315.
- Kundić et al. (1997) Kundić, T., et al., 1997. ApJ 482, 75.
- Lehár et al. (2000) Lehár, J., et al., 2000. ApJ 536, 584.
- Lehár et al. (1992) Lehár, J., Hewitt, J. N., Roberts, D. H., Burke, B. F., 1992. ApJ 384, 453.
- Lubin et al. (2000) Lubin, L. M., Fassnacht, C. D., Readhead, A. C. S., Blandford, R. D., Kundić, T., 2000. ApJ 119, 451.
- Manmoto et al. (1996) Manmoto, T., Takeuchi, M., Mineshige, S., Matsumoto, R., Negoro, H., 1996. ApJ 464, L135.
- Mediavilla et al. (2005) Mediavilla, E., et al., 2005. ApJ 619, 749.
- Mineshige et al. (1994) Mineshige, S., Ouchi, B. N., Nishimori, H., 1994. PASJ 46, 97.
- Negoro et al. (1994) Negoro, H., Miyamoto, S., Kitamoto, S., 1994. ApJ 423, L127.
- Ofek & Maoz (2003) Ofek, E. O., Maoz, D., 2003. ApJ 594, 101.
- Oknyanskij (1997) Oknyanskij, V. L., 1997. In: Hunt, G., Payne, H. E. (Eds.), ASP Conf. Ser. Vol. 125, Astronomical Data Analysis Software and Systems VI, Astron. Soc. Pac., San Francisco, p. 162.
- Oscoz et al. (1997b) Oscoz, A., et al., 1997b. ApJ 479, L89.
- Oscoz et al. (1997a) Oscoz, A., Serra-Ricart, M., Mediavilla, E., Buitrago, J., Goicoechea, L. J., 1997a. ApJ 491, L7.
- Page et al. (2004) Page, K. L., Reeves, J. N., O’Brien, P. T., Turner, M. J. L., Worrall, D. M., 2004. MNRAS 353, 133.
- Paraficz et al. (2006) Paraficz, D., Hjorth, J., Burud, I., Jakobsson, P., Elíasdóttir, Á., 2006. A&A 455, L1.
- Refsdal (1964) Refsdal, S., 1964. MNRAS 128, 307.
- Schmidt & Wambsganss (1998) Schmidt, R., Wambsganss, J., 1998. A&A 335, 379.
- Serra-Ricart et al. (1999) Serra-Ricart, M., et al., 1999. ApJ 526, 40.
- Shakura & Sunyaev (1973) Shakura, N. I., Sunyaev, R. A., 1973. A&A 24, 337.
- Shakura & Sunyaev (1976) Shakura, N. I., Sunyaev, R. A., 1976. A&A 175, 613.
- Simonetti et al. (1985) Simonetti, J. H., Cordes, J. M., Heeschen, D. S., 1985. ApJ 296, 46.
- Steele (2001) Steele, I. A., 2001. New AR 45, 45.
- Steele et al. (2000) Steele, I. A., Newsam, A. M., Mottram, C. J., McNerney, P., 2000. In: Kibrick, R., Wallander, A. (Eds.), Proc. SPIE Vol. 4010, Advanced Global Communications Technologies for Astronomy, p. 133.
- Steele et al. (2004) Steele, I. A., et al., 2004. In: Oschmann, J. M. (Ed.), Proc. SPIE Vol. 5489, Ground–based Telescopes, p. 679.
- Ullán et al. (2006) Ullán, A., et al., 2006. A&A 452, 25.
- Uttley et al. (2005) Uttley, P., McHardy, I. M., Vaughan, S., 2005. MNRAS 359, 345.
- Vanden Berk et al. (2004) Vanden Berk, D. E., et al., 2004. ApJ 601, 692.
- Wambsganss & Schmidt (1998) Wambsganss, J., Schmidt, R., 1998. New AR 42, 101.