Free Access
Volume 522, November 2010
Article Number A95
Number of page(s) 14
Section Extragalactic astronomy
Published online 05 November 2010

© ESO, 2010

1. Introduction

Refsdal (1964) was the first to state that gravitationally lensed quasars can be very useful for determining the parameters of our Universe: combined with a model of the mass distribution in the lensing galaxy, the time delay between different lensed images can lead to the determination of the Hubble constant, H0. This motivated many of the early lensed quasars studies and time delay measurements campaigns. Unfortunately, it quickly became clear that the final estimate of H0 is very sensitive to systematic errors in the lens modeling. One way to reduce these systematic errors is to derive accurate relative astrometry of gravitationally lensed images and lens galaxy light profiles based on high resolution frames. A good example of the consequences for H0 of accurate astrometry is the quadruply lensed quasar PG 1115+080. Courbin et al. (1997) found stronger constraints on the lensed images and lensing galaxy positions, which reduced the degeneracies and the range of acceptable models. This allowed them to measure H0 with an accuracy two times greater than in previous studies. Keeton & Kochanek (1997) also demonstrated that reducing the error in the lens galaxy position of PG 1115+080 from 50 mas to 10 mas could improve the constraints on the lens model and thus on H0. Lehár et al. (2000) highlighted that a poor knowledge on the position of the lens galaxy in two different systems, i.e. B0218+357 and PKS 1830-211, prevented them from being able to accurately measure the Hubble constant. Systematic studies of the effect of the astrometric accuracies on H0 are difficicult because they depend on the lens system configuration.

Strong lensing is also a promising tool for estimating the amount (possibly as a function of redshift) of dark matter clumps (which, following other authors, we hereafter call “substructures”) in distant galaxies and compare this estimated amount to predictions of numerical simulations (see e.g. Zackrisson & Riehm 2009; Koopmans et al. 2009a). The first piece of evidence that strongly lensed quasars are sensitive to substructures in galaxies comes from the so-called “anomalous flux ratios”: for many systems, the flux ratios between the lensed images deviate from those predicted by simple lens models (Kochanek 1991; Mao & Schneider 1998; Dalal & Kochanek 2002; Keeton et al. 2003). It has long been understood that substructures act only on the image flux ratios because of the dependence of the latter on the second derivative of the gravitational potential. However, recent work has explored two new routes to detect substructures in lensed quasars. One method involves using time delay measurements, as shown by Keeton & Moustakas (2009) who highlighted the small changes in time delays caused by substructures. Even if these delays are likely to be modified only by a few tenths of a percent, future large monitoring campaigns should allow the detection of the signature of substructures (Moustakas et al. 2009). Another method proposes detecting substructures in the lensing galaxy by means of their effects on the position of lensed images. The amplitude and probability we should expect for this phenomenon remains disputed. On the one hand, observable astrometric perturbations should be produced by the most massive substructures. Because of the scarcity of high-mass dark matter clumps, Metcalf & Madau (2001), however, derived a low probability and, on average, low astrometric perturbations. On the other hand, Chen et al. (2007) showed that lower mass substructures also play a role. Including a large range of sub-halo masses, they find that substructures could induce astrometric perturbations as large as 10 mas (see Zackrisson & Riehm 2009,for a more complete review). Observationally, astrometric perturbations caused by substructures were detected in a few systems, the most remarkable ones being MG2016+112 (Koopmans et al. 2002; More et al. 2009) and B0128+437 (Biggs et al. 2004). In both cases, the anomalies have been unveiled thanks to high resolution radio images.

Although we cannot yet reach the spatial resolution of the Very Long Baseline Array (VLBA) in the optical range, it has been shown by some of us (Chantry & Magain 2007) that a sophisticated deconvolution technique (ISMCS1) applied to Hubble Space Telescope (HST) images could lead to relative positions of lensed quasar images with milliarcsecond (mas) accuracy, reducing the error bars by a factor  > 2 relative to other techniques. In the present paper, we apply this technique to a sample of seven gravitational lenses without measured time delays. All these systems are photometrically monitored by the COSMOGRAIL2 collaboration and their time delays should be measured in the near future. The goals of this paper are twofold. First we wish to provide shape parameters for the lensing galaxy and accurate relative astrometry for these systems. Then we wish to determine simple lens models which reproduce the observed configuration. From these models, prospective time delays are also calculated, complementing time delays predicted with non-parametric modeling (Saha et al. 2006). Second, we systematically investigate, for quadruply imaged quasars, the ability of simple smooth models to reproduce the image configuration within a few milliarcsec. From this systematic and uniform approach, we wish to test whether the data show evidence of astrometric perturbations caused by substructures.

The lens sample studied in this paper is composed of 7 different systems without time delay measurements and currently monitored by the COSMOGRAIL collaboration: 4 doubly imaged quasars for which no detailed modeling and/or relative astrometry has ever been published and 3 quadruply imaged quasars, amongst which 2 have not yet been studied in detail either (no modeling, no time delay, and/or no lens redshift). The ISMCS deconvolution of the gravitational lenses for which time delays have already been measured will be presented in another paper (Chantry et al., in prep.).

The studied sample is detailed in Sect. 2 and the data in Sect. 3, while the image processing technique is explained in Sect. 4 along with the results. The modeling strategy is explained in Sect. 5. A discussion of the models is presented in Sect. 6. We then conclude in Sect. 7.

2. An overview of our sample

We present the seven gravitationally lensed quasars of our sample, the right ascension and declination being expressed in the J2000 coordinates system:

  • HE 0047-1756 (a). This object (RA = 00h50m27ṣ82 and Dec = ) was discovered by Wisotzki et al. (2000) in the framework of the Hamburg/ESO Survey (HES) for bright quasars, covering the Southern sky. It was later identified by Wisotzki et al. (2004) as a doubly imaged quasar at a redshift of zs = 1.68. The lens is an elliptical galaxy with a spectroscopic redshift of zl = 0.407 ± 0.001 (Eigenbrod et al. 2006b; Ofek et al. 2006).

  • RX J1131-1231 (b). This quadruply imaged quasar (RA = 11h31m55ṣ39 and Dec = ) was discovered in 2003 by Sluse et al. They found a redshift of zl = 0.295 ± 0.002 for the lens while the source lies at zs = 0.657 ± 0.001. Preliminary time delays have been proposed by Morgan et al. (2006) and revised estimates will be published in Kozłovski et al. (in prep.). The system was characterized in details in terms of astrometry and photometry by Sluse et al. (2006). Claeskens et al. (2006) modeled it and also reconstructed the source, which appears to be a type 1 Seyfert spiral galaxy. A similar analysis was performed by Brewer & Lewis (2008) using a Bayesian approach. Substructures in the main lens were searched for by Sugai et al. (2007), while Sluse et al. (2007) and Dai et al. (2010) used microlensing to study the quasar source.

  • SDSS J1138+0314 (c). This quadruply imaged object (RA = 11h38m03ṣ70 and Dec = ) was discovered in 2008 in the Sloan Digital Sky Survey (SDSS) by Inada et al. The redshifts of the quasar and the lens were measured by Eigenbrod et al. (2006b) and are respectively equal to zs = 2.438 and zl = 0.445 ± 0.001. No detailed modeling nor time delay have ever been published for this system.

  • SDSS J1155+6346 (d). This doubly imaged quasar (RA = 11h55m17ṣ34 and Dec = ) was discovered by Pindor et al. (2004) in the SDSS data set. They measured the redshifts of the quasar and the lens to be zs = 2.888 and zl = 0.176, respectively. They also found that one of the two images of the quasar is very close to the lensing galaxy (at around 10% in effective radius off the center of the lens) and is the brightest. That configuration cannot be reproduced by a simple model of mass distribution.

  • SDSS J1226-0006 (e). This system (RA = 12h26m08ṣ02 and Dec = ) is a doubly imaged quasar discovered in the framework of the SDSS by Inada et al. (2008). The quasar is located at a redshift of zs = 1.125. According to Eigenbrod et al. (2006b), the lens is likely to be an early-type galaxy with a spectroscopic redshift of zl = 0.516 ± 0.001. This system has no measured time delay, no published relative astrometry and no detailed modeling study.

  • WFI J2026-4536 (f). Morgan et al. (2004) discovered this quadruply imaged quasar (RA = 20h26m10ṣ43 and Dec = ) during an optical survey using the WFI camera mounted on the MPG/ESO 2.2 m telescope operated by the European Southern Observatory (ESO). The redshift of the source is zs = 2.23 while the redshift of the lens is unknown, although the latter is clearly visually detected on high resolution images. No time delay has ever been measured but according to Morgan et al. (2004), the longest one might be of the order of at most a week or two. According to them, the lensed images are likely affected by microlensing.

  • HS 2209+1914 (g). This system (RA = 22h11m30ṣ30 and Dec = ) is a doubly imaged quasar, with zs = 1.07, discovered during the Hamburg-Cfa Bright Quasar Survey (HS) by Hagen et al. (1999). They clearly detected the lensing galaxy. Nothing else is available for this system: no time delay, no lens redshift, and no modeling.

3. Observational material

The images we analyze were acquired with the camera 2 of NICMOS, i.e. the Near-Infrared Camera and Multi-Object Spectrometer (hereafter NIC2) mounted on the HST. They were all obtained in the framework of the CASTLES project (Cfa-Arizona Space Telescope LEns Survey3, PI: C.S. Kochanek), and are available in the HST archives. The filter used was the F160W, which is very close to the H-band filter. It was selected for several reasons. On the one hand, the PSF is well-sampled and only slightly variable across the field. On the other hand, the lensed-quasar flux ratios are less affected by microlensing effects, given the size of the quasar at these wavelengths and the dust extinction by the lensing galaxy. Details about the image acquisition are summarized in the first columns of Table 1: the name of the object, the date of observation, the number of frames, and the total exposure time. All the frames were obtained after the installation of the NICMOS Cooling System, or NCS, in 2002. Every image was acquired with dithering in the MULTIACCUM mode, each one of them being a combination of about twenty subframes. As these objects were all observed between October and December 2003, the pixel size of the detector on the sky does not change from one target to the other, also because the plate scale of NICMOS has become very stable since the installation of the NCS. The values that we use were measured during part b of the third Servicing Mission Observatory Verification, SMOV3b, and are x = 0.075948″ and y = 0.075355″ (STScI NICMOS Group 2007).

Table 1

General information about the acquisition of HST/NIC2 images and the application of ISMCS.

4. ISMCS on HST/NIC2 images

4.1. Technique

To extract accurate spatial and shape parameters from our data, we need a method capable of separating the contributions of the lensed point sources from those of the more diffuse components (galaxies, halos, arcs, rings, ...). This is exactly what the MCS deconvolution algorithm (Magain et al. 1998) provides. One of the advantages of this deconvolution method relative to other techniques is that it does not violate the sampling theorem. In practice that means that we do not try to fully deconvolve an image to obtain an infinite resolution. We instead choose a resolution for the final deconvolved image, in our case a Gaussian with a Full-Width-at-Half-Maximum (FWHM) of 2 pixels, and we deconvolve our images with a partial point spread function or PSF (which provides the total PSF when reconvolved with our 2 pixel FWHM Gaussian). To achieve this task, we need to know very well the shape of the PSF. As the NIC2 field is only 19″̣2  ×  19″̣2, we do not have the possibility of using field stars to determine the PSF. Since the lensed quasar images are contaminated by the lensing galaxy or partial Einstein rings beneath them, we cannot use these images directly to improve our PSF. We instead use ISMCS (see Chantry & Magain 2007,for further details), a special iterative strategy coupled with the MCS algorithm. The HST PSF is quite complex (including spike-like features and an intense first Airy ring), hence we start the deconvolution process using a PSF created by the Tiny Tim software (Krist & Hook 2004) as a first guess of the true PSF. We improve the Tiny Tim function by adjusting it simultaneously for all the point sources of a frame using a technique described in Magain et al. (2007), which allows us to add a numerical component to the input PSF so that it is better adapted to the actual frame. We then obtain a set of modified PSFs that we use to simultaneously deconvolve all the frames at our disposal, with a sampling step two times smaller than the original one. In doing so, we obtain a first approximation of the diffuse background and after reconvolving it to the inital resolution, we subtract it from the original frames and obtain new ones, partially cleaned from the extended structures. In these modified frames, which contain point sources that are less contaminated by smooth structures, we once again improve our PSFs. This iterative process has to be repeated until the reduced chi-square, (Eq. (1)), reaches a value close to unity in an area determined by the maximum extension of the residual structures after the first deconvolution, and until the residuals are sufficiently flat (no sharp structure). In practice, we stop when an additional iteration no longer improves the (typically, ). For an image with N pixels, the latter is defined to be (1)where ℳ(x) is the model reconvolved by the partial PSF, 𝒟(x) is the observed signal, and σ(x) is the standard deviation associated with that signal. In practice and as an improvement relative to the process applied to the Cloverleaf in the original paper of Chantry & Magain (2007), we note that convergence is reached more rapidly when performing a first simultaneous deconvolution of all the frames with Tiny Tim PSFs instead of first trying to improve these PSFs on the unmodified images, i.e. containing the untouched background structures. With this first deconvolution, we obtain a map of the diffuse structures. This map has then to be cleaned from some artificial ring structures present only to compensate for inaccurate PSFs. It can then be subtracted from the original frames, which can then be used to improve for the first time the PSFs. That first deconvolution is accounted for as one iteration.

4.2. Results

The original frame (combination of all observations), the deconvolved image4 obtained from the last iteration of ISMCS, and the mean residual maps5 after the first and last iterations are displayed for each system in Fig. 1. Both residual maps are expressed in units of σ and their color scale ranges from  − 5 in black to +5 in white. The black rectangle delimits the zone used to estimate the reduced χ2, the orientation of these rectangles being identical to that of the original frames. We emphasize that the PSF used for the first iteration is the one created by the Tiny Tim software. When we examine the residual maps, the improvement brought by ISMCS is undeniable. Moreover, in most cases (5 amongst 7), the remnant structures underneath the point sources (on the residual map from the last iteration) disagree with each other, which is the sign of a variable PSF throughout the detector, even on small spatial scales. The number of iterations necessary to reach convergence is shown in the last columns of Table 1 along with the values of the after the first and the last iterations.

The astrometry, corrected from the X / Y scale difference and the distortions of NIC2, and the photometry (Vega system) are shown in Table 2. The  ± 1σ error bars were calculated in deconvolving each frame individually at the last iteration and in determining the dispersion around the mean. They are very small because they are inherent to the deconvolution technique: no external systematic error is included in these error bars. To estimate the total error, we compare the spatial extension of each object on the detector to the one of the Cloverleaf (H1413+117). The latter was used as a test of ISMCS in Chantry & Magain (2007): in comparing the astrometry of the point sources obtained in two different filters and at two different orientations on the sky, they could estimate the total error to be 1 mas, accounting e.g. for a possible remnant distortion in the images. The estimated total errors based on the Cloverleaf are displayed in the fifth column of Table 2. Since they are based on the maximum extension of the object regardless of the direction, they should of course be considered upper limits. That is why they are called “MTE”, which stands for “maximum total error”.

Since the total error derived in Chantry & Magain (2007) for H1413+117 was based on a comparison of the relative positions of the lensed images obtained at different NIR wavelengths and image orientations with the same instrument, we attempted to derive independent estimates based on the comparison of the relative astrometry derived with HST and with high resolution radio data. In a future paper treating the lenses with already measured time delays (Chantry et al., in prep.), we will present HST astrometry for the radio quad JVAS B1422+231 (Patnaik et al. 1992). To estimate the error affecting our results, we chose one lensed image as an astrometric reference and calculated the distance between it and every other lensed image. We then measured the difference between the distances obtained with our positions and those calculated with the radio astrometry of Patnaik et al. (1999). The scatter in these differences of distance around the mean is about 2.6 mas. Assuming the uncertainty is the same in any direction, we derived an error in the relative astrometry of 1.8 mas in both RA and Dec. This value is larger then the MTE of 1.05 mas derived for B1422+231 from our standard method. This is expected as the radio emission in B1422+231 is slightly extended and is likely not to originate in the accretion disk (as observed in the optical range) but rather from a nearby region at the basis of the radio jet. This effect, known as core-shift, is observed between two different radio-bands (Porcas 2009; Kovalev et al. 2008) and may induce astrometric perturbations as large as a few mas on the relative astrometry of lensed quasar images (Mittal et al. 2006). Thus, it appears that the use of radio data as an independant calibrator of the systematic errors is difficult at such an accuracy. It requires the comparison of the relative astrometry, for several objects, between radio and optical wavelengths, a task that is beyond the scope of the present paper.

Table 2

Relative position, maximum total error (“MTE”), magnitude, and flux ratio of the lensed images and lensing galaxy (see Fig. 1 for the labels).

In addition to the point-source deconvolution, we used an analytical model to characterize the lensing galaxy light distribution. To ensure in this case that the maximum amount of light of the galaxy is included in the profile, the deconvolution was performed with no numerical component. Since most of the lensing galaxies are ellipticals, we used a de Vaucouleurs light profile (de Vaucouleurs 1948). This procedure allowed us to extract the galaxy shape parameters summarized in Table 3: the position angle or “PA” (orientation in degrees east of north) of the galaxy, its ellipticity (e = 1 − b/a), the effective semi-major and semi-minor axis (resp. aeff and beff). The effective radius Reff is calculated to be the geometrical mean between the two effective semi-axis (Kochanek 2002). These three quantities are expressed in arcseconds. The  ±1σ error bars were also calculated in deconvolving each frame individually and determining the dispersion around the mean. We note that the luminosity of the galaxies displayed in Table 2 is measured in an aperture equal to Reff.

4.3. Discussion

A few remarks can be made about the results of the deconvolution:

  • HE 0047-1756 (a). A faint Einstein ring, which is a stretched image of the quasar host galaxy, is revealed by the deconvolution.

  • RX J1131-1231 (b). Sluse et al. (2006) reports astrometric measurements on the same frames with the MCS deconvolution algorithm but with no iterative strategy. Their results agree within the error bars with those presented here. An offset of up to 3 mas between both results is observed. This difference is probably caused by the significant brightness of the Einstein ring. The different amount of recovered background under the PSF can lead to a small shift in position. Photometry is also affected by the presence of the ring. We derive an absolute photometry about 0.4 mag brighter than Sluse et al. (2006) but obtain compatible flux ratios. The remnant systematic structures in the final residual map are also caused by this very bright ring, which affects the PSFs and degrades their quality. Indeed, because a part of the background is identical for the three brightest lensed images (and thus with more weight in the determination of the PSF), it is impossible to completely disentangle the flux contribution of the ring from the one really included in the point sources. This is a limitation of the ISMCS method: for it to work properly, the background has to be different for each point source.

  • SDSS J1138+0314 (c). A faint Einstein ring is revealed by the deconvolution process.

  • SDSS J1155+6346 (d). Our astrometry is not in agreement with Pindor et al. (2004) especially concerning the lens. The difference for source B amounts to 0″̣1 in RA and 0″̣04 in Dec, while for the lens the offset is much larger and amounts to 1″̣55 in RA and 0″̣28 in Dec. However, our astrometry is in agreement with the one listed in the CASTLES database. The remnant systematic structures in the final residual map are due to the presence of the very bright and extended lensing galaxy (which is not clear on the presented frames, the cuts being chosen so that the two lensed images appear clearly). As in the case of RX J1131-1231, it is not possible to completely disentangle the background flux contribution from the one of the point sources, a part of the background being identical for the two lensed images. This degrades the quality of both the PSFs and the deconvolution.

  • WFI J2026-4536 (f). The astrometry we obtain, except for the right ascension of the lens (), is in agreement, within the error bars, with the results of Morgan et al. (2004), who used the same frames but a different image processing technique.

  • HS 2209+1914 (g). A bulge is clearly observable but an additional extended structure is also visible. It could either be some spiral arms, in which case the lens would be a late-type galaxy, or even a distorted Einstein ring. A spectrum of the lens and higher resolution imaging would help us to differentiate between these hypotheses. As we do not know what this structure is related to, we fit a de Vaucouleurs model to the bulge only, in using a special feature of the MCS algorithm: a mask encircling the lens galaxy to avoid the model fitting this extended structure. Moreover, the residual map contains many intense structures. However, as these structures do not have the same shape beneath both point sources, we cannot recover them with another iteration of ISMCS. Saturation is unlikely to be responsible for that phenomenon, as it is corrected by the NICMOS reduction pipeline. It may be a consequence of differential extinction by the lensing galaxy, resulting in a different color for both lensed images and thus a different shape of the PSF.

We note that all our results agree with those found in the CASTLES database, within their error bars (ours being smaller).

thumbnail Fig. 1

HST/NIC2 original and deconvolved frames (resp. top left and top right), mean residual maps from the first and from the last iteration of ISMCS (resp. bottom left and bottom right).

Open with DEXTER

Table 3

Measured shape parameters for the lensing galaxy.

5. Parametric modeling

Our goals are twofold. First, we aim to provide simple models and prospective time delays for the lensed quasars monitored by COSMOGRAIL. Second, we aim to test whether simple smooth lens models are able to reproduce the mas relative astrometry of quadruply imaged quasars in 3 systems without measured time delays. Our strategy consists of using the LENSMODEL software package v1.99o (Keeton 2001) to model the mass distribution of our seven systems. For a chosen model of the mass distribution, the code minimizes a χ2 defined as the square of the difference between observable quantities and their “model counterparts”, weighted by the observational errors on these quantities. Two different lens models are considered. First, an isothermal profile, which is the standard mass distribution used to model gravitational lenses (Kassiola & Kovner 1993), and second, a de Vaucouleurs profile, for which we assume that the light perfectly traces the mass in the inner regions of lensing galaxies. These two models should provide a good approximation of the extreme slopes of the mass distribution at the location of the lensed images and the expected time delays (Kochanek 2002; Kochanek & Schechter 2004). In addition, the study of the galaxy-galaxy lensing sample from the Sloan Lens ACS6Survey (SLACS, Bolton et al. 2006) revealed that the massive elliptical lensing galaxies are nearly kinematically indistinguishable from isothermal ellipsoids (see e.g. Koopmans et al. 2009b). This supports the use of an isothermal gravitational potential as a fiducial model to test the ability of smooth lens models to reproduce quadruply imaged quasars with mas accuracy. Since lensing galaxies are never isolated, we model the effect of the environment using an external shear term characterized by an amplitude γ and a position angle θγ (pointing towards the mass at the origin of the shear). All the models are computed for a flat universe with the following cosmological parameters: H0 = 70 km s-1 Mpc-1, Ωm = 0.3, and ΩΛ = 0.7.

To model our systems, we use every constraint at our disposal: the relative astrometry of the lensed images, with the MTE, i.e., the uncertainties displayed in the fifth column of Table 2, the position of the main lens, with the error inherent to the deconvolution when it is larger than the MTE, and, in the case of doubly imaged quasars, the flux ratio of the two point sources. In principle, fluxes can be contaminated by different effects such as microlensing by stars in the lens galaxy, dust extinction, and also by the time delay itself. As the flux ratios are measured in the near-infrared, all these effects should be small (Yonehara et al. 2008). We thus assume a 1σ error of 10% for the flux ratios. In summary, we have 10 constraints for the quads, while we have 8 for the doubles. For the de Vaucouleurs model, we assume that the total mass profile follows the light profile. We thus add three constraints to the model: the PA of the galaxy, its ellipticity e, and its effective radius Reff (see Table 3). Because of the limited number of observational constraints, isothermal mass profiles are assumed to be spherically symmetric (SIS, i.e. Singular Isothermal Sphere) when modeling doubles. This is not a strong assumption as the quadrupole term of the potential modifies only slightly the time delays of doubly imaged quasars (Kochanek 2002; Wucknitz 2002). For quads, we allow the ellipticity of the isothermal mass distribution (SIE, i.e. Singular Isothermal Ellipsoid) to deviate from the ellipticity of the light profile. This enables to account for dark matter halos that may be rounder than the light distribution (Ferreras et al. 2008). The position angle of the total mass distribution can be constrained to be that of the light profile because these two distributions might only be slightly misaligned (Keeton et al. 1997; Ferreras et al. 2008). Finally, we also assume that the center of the total mass distribution and that of the light profile are identical within the error bars. This is supported by the work of Yoo et al. (2006) who found, for 4 lensed quasars with an Einstein ring, that the offset between the light and the total mass distribution is limited to a few mas. Calculating the number of degree(s) of freedom (d.o.f.), which is the difference between the number of model parameters and observable quantities, we find 0 d.o.f. when modeling doubly imaged quasars and 2 (resp. 3) d.o.f. when modeling quads with SIE (resp. de Vaucouleurs) + external shear.

The search for the best model and estimate of uncertainties is performed in two steps. First, we generate an initial sample of 2000 different models with parameters distributed over the whole parameter space and optimize them. This method is efficient for finding the best models and identifying local minima. To estimate the model uncertainties, we then sample the posterior probability distribution of the parameter space using an adaptive Metropolis Hastings Monte Carlo Markov Chain (MCMC) algorithm. This technique is implemented in LENSMODEL and described in Fadely et al. (2010). In practice, an ensemble of 15 different chains are run, each chain consisting of a sequence of trial steps drawn from a multivariate Gaussian distribution of width estimated during the first step of the process. The sampling of the parameter space is optimized by using the covariance matrix. In 5% of the steps, the covariance matrix is diagonal, allowing us to use a large step along one of the axis and to escape local minima in the χ2 surface. We use the same criterion as Fadely et al. (2010) to assess that any MCMC run has converged. Finally, for each point of the MCMC, we calculate the relative likelihood of a parameter p based on the χ2 statistics (i.e. L(D | p) = exp( − χ2 / 2)), and calculate a 68% confidence interval for each parameter.

The parameters of the best fit models are displayed in Table 4. The columns display the following items: the name of the object, the type of mass distribution used (“DV” stands for de Vaucouleurs profile), the mass scale parameter (the angular Einstein radius REin in arcseconds), the mass distribution ellipticity e and its orientation θe in degrees positive east of north, the effective radius Reff in arcseconds in the case of a de Vaucouleurs model, the intensity of the shear γ and its orientation θγ in degrees (east of north), the corresponding flux ratios of the lensed images, the number of degree(s) of freedom (d.o.f.), the χ2 of the fit, and the predicted time delays in days when the lens redshift is known. For the quads and in the same column as the χ2, we also give the , which is the contribution of the lensed images position to the χ2, and , which is the contribution of the lens galaxy position to the χ2. We note that ΔtAB >  0 means that the flux of A varies before that of B. The median value of each parameter and its 68% confidence level is shown in Table 5.

Table 4

Results of the parametric modeling.

Table 5

Median value of the model parameters and 68% confidence interval.

6. Discussion

6.1. Doubles

For the doubly imaged quasars, both SIS+shear and DV+shear models can reproduce the image configuration as well as the flux ratio, even with our constraints on the shape of the galaxy in the case of a de Vaucouleurs profile. Two systems require a large amount of shear (i.e. γ > 0.1 for both mass models) to reproduce the lens configuration: SDSS J1226-0006 and SDSS J1155+6346. For SDSS J1226-006, the HST/NIC2 images show a galaxy G2 at and from image A (3″̣4 from the main deflector), about 15° off the direction of θγ. This galaxy, whose type is unknown, is likely not the only source of shear. The luminosity ratio of G2 to the lens is Llens/LG2 = 4.8. Assuming we can use the Faber-Jackson relation (L ∝ σ4, Faber & Jackson 1976), this ratio implies that σlens/σG2 = 1.5, σlens and σG2 being respectively the velocity dispersion of the lens and of G2. The isothermal model allows us to translate REin of the lens into σlens. We find σlens = 212 km s-1 and thus σG2 = 141 km s-1. Using formula A.20 of Momcheva et al. (2006) and assuming G2 is at the same redshift as the lens, this induces a shear of γ = 0.039, more than 2 times smaller than the one predicted by the SIS model. Other galaxies in the field are probably responsible for the remaining shear.

A more dramatic case is SDSS J1155+6346, for which models predict a shear as large as 0.4 to reproduce the observed configuration. This is one of the largest shears needed to reproduce a lensed quasar system. On some larger field images of this object (obtained with ACS onboard HST, PI: C.S. Kochanek), we do not see any bright galaxy in its vicinity. We thus suspect that a massive galaxy cluster lies outside the ACS field, though nothing is clearly visible on the SDSS data7. Deeper images would be necessary to exclude or confirm the existence of this cluster.

In the case of HE 0047-175, a diffuse component lies at and from image A ( from the lens), in the direction of the shear (see Fig. 1) independent of the employed model. Although very faint (about 2 mag fainter than the lens), this galaxy likely represents the major contribution to the shear in this system. A SIS model with σ = 88  kms-1 would indeed produce the observed amount of shear, if located at the position of this faint companion (assuming zcomp = zlens = 0.407).

6.2. Quads and astrometric anomalies

The quadruply imaged quasars allow us to test the ability of simple smooth models to reproduce a relative astrometry with mas accuracy. Only the relative astrometry of SDSS J1138+0314 is easily reproduced with our models. In contrast, for WFI J2026-4536 and RX J1131-1231, we find that a very large χ2 is associated with our models. In the first case, the main contribution to the χ2 comes from the difference between the PA of the model and the PA of the light distribution. In the second case, the large χ2 is mainly caused by the impossibility of the model to recover the positions of both the lensed images and the lensing galaxy. In any case, an underestimate of the error bars on the quasar lensed images is unlikely.

For RX J1131-1231, a reduced χ2 ~ 1 can only be obtained if we increase the error bars in the positions of the lensed images by a factor 10. Alternatively, we also derive a good fit if we allow more freedom in the position of the lensing galaxy (i.e. error = 0.02). Following this procedure, we find that the offset between the light and mass distribution centroid amounts to 88  pc. This value is marginally consistent with the upper limit of 70  pc derived for B1938+666 by Yoo et al. (2006). It is however inconsistent with the maximum offset value of  ~ 20  pc found for the 3 other systems they analyzed, suggesting that the offset between light and mass distributions is not the cause of the astrometric perturbation we observe. Claeskens et al. (2006) also found that models (simple or more complex) fitting simultaneously the Einstein ring of RX J1131-1231 and the quasar lensed images lead to a poor χ2. Brewer & Lewis (2008) were also unable to reproduce the lensed quasar relative positions to mas accuracy using the lens model based on the Einstein ring.

For WFI J2026-4536, an acceptable χ2 cannot be obtained in enlarging the error bars in the positions of the lensed images or of the lens galaxy but only in relaxing the constraint on the position angle of this latter. This could be due to a galaxy located at and (Morgan et al. 2004) from image B, which we ignored in the modeling.

For the three quads, the flux ratios predicted by our best-fit models differ from those measured on the HST/NIC2 frames. This discrepancy may have several origins such as invalid assumption about the lens model, i.e. a need for multipole components (Evans & Witt 2003), microlensing, and/or massive substructures in the lensing galaxies (see Keeton et al. 2003, 2005, for an exhaustive discussion).

The amount of shear needed in our models to reproduce the configuration of the three quads is quite high (γ ≳ 0.1). Our model estimates can be compared to the minimum amount of shear required to reproduce the image position, following the methodology described in Witt & Mao (1997). Using Eq. (20a) of Witt & Mao (1997), we find a minimum shear γmin of respectively 0.035, 0.004, and 0.062 for RX J1131-1231, SDSS J1138+0314, and WFI J2026-4536. The values found for the last two systems are much smaller than the one predicted by our models. This is not surprising as we find θγ and θe to be orthogonal, that implies that γ is in fact poorly known and strongly degenerate with the internal shear. For RX J1131-1231, γmin is only 3 times smaller than the value derived for our SIE+γ model. The difference between γmin and γobs is due to θe being nearly aligned with θγ (about 15° offset), a situation which also leads to a significant underestimate of γ using γmin. By using Eq. (22) from Witt & Mao (1997), we can also derive a range of allowed values for θe based on the image configurations. For WFI J2026-4536, we find that the observed θe is far from the allowed values, confirming the results of our models, which indicate a likely offset betwen the mass and the light matter distribution. For the two other systems, the observed θe falls at the limit of the allowed range, suggestive of a significant ellipticity of the lens, as observed.

The previous results suggest that it is quite common for simple lens models to fail in reproducing mas astrometry of quadruply imaged quasars. To investigate the question, we searched the literature for lensed quasars having images with mas astrometric error bars (i.e. up to 0″̣002 on the lensed image positions) and published simple models. We found eleven systems gathering these conditions8: B0128+437 (Biggs et al. 2004), MG0414+0534 (Ros et al. 2000), HE0435-1223 (Morgan et al. 2005; Kochanek et al. 2006), SDSS0924+0219 (Keeton et al. 2006; Eigenbrod et al. 2006a), H1413+117 (Chantry & Magain 2007; MacLeod et al. 2009), B1422+231 (Bradač et al. 2002), B1608+656 (Koopmans et al. 2003), B1933+5039 (Cohn et al. 2001), MG2016+112 (Chen et al. 2007; More et al. 2009), WFI2033-4723 (Vuissoz et al. 2008), and B2045+265 (McKean et al. 2007). Two of these systems (H1413+117, B1933+503) are easily reproduced by simple models because of the large uncertainty affecting the position of the lensing galaxy (). Out of the six systems for which VLA, VLBA, or VLBI data are available10 (B0128+437, MG0414+0534, B1422+231, B1608-656, MG2016+112, B2045+265), only B1422+231 shows convincing evidence that smooth models allow one to reproduce the relative astrometry, although substructures are needed to reproduce the flux ratios (Bradač et al. 2002). For MG0414+0534, B1608+656, MG2016+112, and B2045+265, complex models including a bright susbtructure (MG0414+0534) or a companion galaxy (B1608+656, MG2016+112, B2045+265) are needed to obtain acceptable fits. The case of B0128+437 is a bit peculiar as the lensing galaxy was not detected at the time of the modeling paper11 and so, only the lensed images were used for the models. As this system contains more than four lensed images, the models can be constrained. Based only on the position of the lensed images, an astrometric anomaly is detected as long as an astrometric accuracy smaller than 1 mas (but still compatible with the data error bars) is considered. Among the last three systems (HE0435-1223, SDSS0924+0219, WFI2033-4723), constrained by relative astrometric positions derived from HST images, a good fit is obtained only for HE0435-1223. A SIE+γ model leads to for SDSS0924-0219 (Keeton et al. 2006) and to for WFI2033-4723. In both cases, relaxing the constraint on the lens galaxy centroid allows one to derive a perfect fit to the astrometry. For WFI2033-4723, more complex models including a nearby group allow Vuissoz et al. to reproduce the astrometry but they also need the light and mass distribution to be misaligned, which is not totally satisfactory (Vuissoz et al. 2008). Kochanek & Dalal (2004) compile seven quads (some of them also compiled here) for which they fitted SIE+γ models. Unfortunately, only sparse information about the astrometric error used is available and we cannot infer any trend from this study.

To conclude, out of the nine usable systems amongst the eleven quads, at least four show astrometric perturbations with respect to predictions of simple lens models (B0128+437, MG0414+0534, SDSS0924+0219, WFI2033-4723). For three of the remaining systems (B1608+656, MG2016+112, B2045+265), conclusions are difficult to draw because the need to include a companion galaxy comes naturally from deep near-IR imaging. The last two systems (HE0435-1223, B1422+231) are reproduced well by simple smooth models.

Although the considered sample of quads represents heterogeneous data sets and analyses, it indicates that the relative astrometry of quads often deviates from simple model expectations when trying to reproduce it to mas precision. The considered sample suggests that the situation is less critical for “central quads” (i.e. with the source lying close to the center of the central astroid caustic) than for fold systems (i.e. sources lying close to a fold caustic). This might be a normal geometrical effect (image positions vary more slowly when moving the source in the central region of the astroid caustic) but it might also be caused by substructures producing severe deformations of the caustics (Bradač et al. 2002). It remains to be seen how significant this effect is with respect to the relative astrometric uncertainty in the image positions or the amount of shear.

Substructures are not the only explanation of the frequent inability of simple lens models to fit the configuration of quads. Other possible explanations are astrometric perturbations due to the lens environment, asymmetries in the mass distribution, disky/boxy projected mass profiles, and offsets between the galaxy light centroid and mass centroid. The last two possibilities seem however to be ruled out by Kochanek & Dalal (2004) and Yoo et al. (2006). The evidence that bright substructures/nearby satellite galaxies explain astrometric perturbations of some systems suggest that substructures may be one of the major contributors to the astrometric perturbations of quads.

All this motivates a systematic study of the ability of simple models to reproduce the configuration of quads, with good control on the error estimates and uniform modeling. This work is beyond the scope of this paper and delayed to a forthcoming paper, when the iterative deconvolution method will have been applied to a larger number of quadruply imaged quasars.

7. Conclusions

In applying ISMCS, i.e. the MCS deconvolution algorithm combined with an iterative strategy, to HST/NIC2 images of seven lensed quasars, we have obtained accurate relative positional constraints on both the lensed images and lensing galaxy, and the shape parameters of this galaxy. We have achieved an accuracy of around 1–2 mas on the lensed image positions. We have also detected for the first time a partial Einstein ring in two cases, the double HE 0047-1756, and the quadruple SDSS J1138+0314, and we highlight the already known ring in RX J1131-1231. Deeper images are needed to perform clear source reconstruction of HE 0047-1756 and SDSS J1138+0314. In the case of HS 2209+1914, the deconvolved frame reveals a structure around the bulge of the lens galaxy that cannot be clearly identified. This structure may be either an Einstein ring or some poorly resolved spiral arms of a large late-type galaxy. This question probably deserves further study: a spectrum of the surrounding structures could give extra information about their true nature.

We also obtain simple mass models for every system. In the case of doubles, both the isothermal and de Vaucouleurs profiles can reproduce the observed configuration. For SDSS J1155+6346, a good fit can only be reached with an extremely and anomalously high external shear, 0.392 in the case of a SIS and 0.453 in the case of a DV, which is indicative of there being a galaxy group or cluster located outside the field of view of the ACS. For SDSS J1226-006, the large shear (γ = 0.1) is probably partially due to a nearby galaxy located 3″̣4 from the main deflector. In the case of quads, a good χ2 can only be obtained for one object: SDSS J1138+0314. The two other quads of our sample, RX J1131-1231 and WFI J2026-4536, need more complicated models to account for their observed configuration. For RX J1131-1231, the offset between the light and mass distribution cannot account for the astrometric perturbation we observe.

Our study of the literature allowed us to conclude that most of the quads cannot be modeled with simple profiles when the astrometic accuracy reaches around 1 mas: some require the presence of companion galaxies, some others need substructures. This finding motivates the acquisition of mas astrometry for all the quads and their mass modeling.

In the framework of the COSMOGRAIL collaboration, the next step for these seven systems is the acquisition of well-sampled light curves to extract time delays. Then, if the redshift of the lens is known, our astrometric constraints will help us to reduce the systematic errors in the Hubble constant.


Iterative Strategy combined with the MCS deconvolution alogrithm.


COSmological MOnitoring of GRAvItational Lenses;


The labels of the lensed images are the same as in previous studies if any.


The residual map is the image of the difference between the model and the original frame in units of sigma.


ACS stands for Advanced Camera for Surveys.


Although B0712+472 has accurate astrometry, the model published by Jackson et al. (1998) is provided without information about the χ2, thus not allowing us to estimate the quality of the fit.


B1933+503 is actually a ten-image lens. Only the position of four of them is known with a precision of 2 mas or less.


B1933+503 also has VLA data but has already been ruled out, its configuration being easily reproduced with simple models because of the large uncertainty in the lens position.


The latter was discovered in AO images by Lagattuta et al. (2010).


This work is supported by the Belgian Federal Science Policy (BELSPO) in the framework of the PRODEX Experiment Arrangement C-90312. V.C. thanks the Belgian National Fund for Scientific Research (FNRS) as well as ESO Vitacura (Chile) for the use of their facilities. A fellowship from the Alexander von Humboldt Foundation to D.S. is gratefully acknowledged. We also thank Ross Fadely and Chuck Keeton for advices on the use of the MCMC implementation of LENSMODEL. And finally we thank the anonymous referee for many precious comments which allowed to improve this paper.


All Tables

Table 1

General information about the acquisition of HST/NIC2 images and the application of ISMCS.

Table 2

Relative position, maximum total error (“MTE”), magnitude, and flux ratio of the lensed images and lensing galaxy (see Fig. 1 for the labels).

Table 3

Measured shape parameters for the lensing galaxy.

Table 4

Results of the parametric modeling.

Table 5

Median value of the model parameters and 68% confidence interval.

All Figures

thumbnail Fig. 1

HST/NIC2 original and deconvolved frames (resp. top left and top right), mean residual maps from the first and from the last iteration of ISMCS (resp. bottom left and bottom right).

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.