EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 671 (March 2023) A&A, 671 (2023) A60 Full HTML
Open Access
Issue
A&A
Volume 671, March 2023
Article Number A60
Number of page(s) 19
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202244680
Published online 06 March 2023

© The Authors 2023

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

In the last few decades, several studies have suggested the presence of a significant amount of dark matter (DM), both in elliptical (e.g., Loewenstein & Mushotzky 2003) and spiral (e.g., Rubin 1983) galaxies. These studies have measured total mass-to-light ratios considerably higher than the stellar ones and observed a remarkable difference between the visible, that is, stellar plus gas, and total mass estimated from galaxy dynamics (e.g., Rubin et al. 1970; Gerhard et al. 2001; Cappellari et al. 2006). Since then, cosmological observations of the cosmic microwave background (Smoot et al. 1992; Hinshaw et al. 2013; Planck Collaboration VI 2020), of baryon acoustic oscillations (Efstathiou et al. 2002; Eisenstein et al. 2005; DES Collaboration 2022), and of Type Ia supernovae (Riess et al. 1998; Perlmutter et al. 1999; Scolnic et al. 2018) have all supported the currently accepted Λ cold dark matter (CDM) model, according to which ≈30% of the Universe is composed of baryons and CDM and the remaining ≈70% is a poorly understood component, known as dark energy, that is responsible for the accelerated expansion of the Universe. However, the ΛCDM model, which has been very successful at large scales (≳1 Mpc), cannot accurately predict the properties of structures at smaller scales, such as the value of the inner slope of DM halos (e.g., Gnedin et al. 2004; Newman et al. 2013a,b; Martizzi et al. 2012). Hence, the interplay between DM and baryons on galactic scales is still being intensively investigated.

N-body simulations of DM particles show that, at equilibrium, they are distributed following an almost universal mass density profile, ρ(r), first described by the Navarro-Frenk-White profile (Navarro et al. 1997). This profile has a characteristic slope, which varies from ρr ≪ rs ∝ r−1 at small radii, and thus has a central cusp, to ρr ≫ rs ∝ r−3 at large radii. Recent simulations with higher resolution (e.g., Golse & Kneib 2002; Graham et al. 2006; Navarro et al. 2010; Gao et al. 2012; Collett et al. 2017; Dekel et al. 2017) and with realistic models that describe the baryonic components predict mass density profiles that deviate from the Navarro-Frenk-White one. This mainly depends on processes such as gas cooling, which allows baryons to condense toward the center of a galaxy (e.g., Blumenthal et al. 1986; Gnedin et al. 2004; Sellwood & McGaugh 2005; Gustafsson et al. 2006; Pedrosa et al. 2009; Abadi et al. 2010; Sommer-Larsen & Limousin 2010), active galactic nuclei feedback (e.g., Peirani et al. 2008; Martizzi et al. 2013; Li et al. 2017), dynamical heating in the central cuspy region, due to infalling satellites and mergers (e.g., El-Zant et al. 2001, 2004; Nipoti et al. 2004; Romano-Díaz et al. 2008; Tonini et al. 2006; Laporte & White 2015), and thermal and mechanical feedback from supernovae (e.g., Navarro et al. 1996; Governato et al. 2010; Pontzen & Governato 2012).

In this context, the observation of very massive (and DM-rich) galaxies and the measurement of their inner total, DM, and baryonic mass profiles represent a key step in understanding how the different mass components are distributed, and inferring how galaxies formed and evolved over cosmic time. These studies can be addressed by exploiting very massive galaxies that act as gravitational lenses.

Strong gravitational lensing represents a powerful tool for studying galaxy evolution and for exploring the properties of the Universe (e.g., Bartelmann 2010; Treu 2010). Given the fact that the deflection of the light emitted by a source only depends on the total gravitational potential of the lens, gravitational lensing is sensitive to both luminous and DM mass and allows one to obtain some of the most precise total mass measurements, with relative errors on the order of a few percent, in extragalactic astrophysics (e.g., Grillo 2010; Zitrin et al. 2012). Together with the results obtained by exploiting the most massive lenses, such as galaxy clusters, about the slope value of the mass density of the DM halos in their inner regions (Sand et al. 2004; Grillo et al. 2015; Annunziatella et al. 2017; Bergamini et al. 2019), the measurements of the cosmological density parameter values (Jullo et al. 2010; Caminha et al. 2016; Grillo et al. 2020), and the value of the Hubble constant (Grillo et al. 2018), galaxy-scale strong gravitational lenses have provided several key results (e.g., Suyu et al. 2013, 2017). In particular, it has been possible to characterize the physical properties of the lens, for instance to reconstruct the total and DM mass distributions, and thus the DM over total mass fraction (Gavazzi et al. 2007; Grillo et al. 2009; Suyu & Halkola 2010; Sonnenfeld et al. 2015; Schuldt et al. 2019) and the mass density slopes (e.g., Treu & Koopmans 2002; Koopmans et al. 2009; Barnabè et al. 2011; Shu et al. 2015), to infer the most likely lens stellar initial mass function (IMF; e.g., Cañameras et al. 2017b; Barnabè et al. 2013; Sonnenfeld et al. 2019), and to identify DM substructures (e.g., Vegetti et al. 2012; Hezaveh et al. 2016; Ritondale et al. 2019). It has also been shown how the values of some cosmological parameters can be measured in strong lensing systems using kinematic data of the lenses (e.g., Grillo et al. 2008; Cao et al. 2012) or in systems where two or more sources are multiply imaged by the same lens galaxy (Tu et al. 2009; Collett & Auger 2014; Tanaka et al. 2016; Smith & Collett 2021), once the mass sheet degeneracy is broken (Schneider 2014).

Moreover, gravitational lensing offers the opportunity to study the lensed background sources: they can be fully reconstructed by taking their surface brightness (SB) distribution into consideration during the strong lensing modeling (e.g., Schuldt et al. 2019; Rizzo et al. 2021; Wang et al. 2022). They can also be analyzed in detail thanks to high magnification factors, which allow one to investigate the local feedback mechanisms driving the evolution of the faintest and smallest high-z galaxies (e.g., Förster Schreiber et al. 2009; Cañameras et al. 2017a; Cava et al. 2018; Iani et al. 2021).

In this paper we measure the total and baryonic mass distributions of the deflector SDSS J010049.18+181827.7 (hereafter, SDSS J0100+1818) included in the Cambridge And Sloan Survey Of Wide ARcs in the skY (CASSOWARY) survey (Belokurov et al. 2009; Stark et al. 2013)1. The exceptionally large Einstein radius of ≈42 kpc, together with the results from our follow-up observations with the Nordic Optical Telescope (NOT), the Very Large Telescope (VLT), and the Hubble Space Telescope (HST), presented in Sect. 2, suggests that the SDSS J0100+1818 deflector is an uncommonly strong lens and a candidate fossil system at intermediate redshift (Johnson et al. 2018).

This paper is organized as follows. In Sect. 2 we present the photometric and spectroscopic observations of SDSS J0100+1818. In Sect. 3 we focus on the main lens galaxy, showing the results of its stellar mass, kinematics, luminosity profile, and environment. In Sect. 4 we perform a strong lensing analysis by exploiting different sets of point-like multiple images and illustrate the best-fit models and the deflector total mass profile. We enhance the analysis in Sect. 5, where we model the deflector by considering the extended SB of the lensed sources. At the end of the section, we study the reconstructed sources assuming different models and comparing their reconstructed sizes. In Sect. 6 we summarize and discuss the main results. Throughout this work, we assume H0 = 70 km s−1 Mpc−1, Ωm = 0.3, and ΩΛ = 0.7. In this model, 1 arcsec corresponds to a linear size of 6.59 kpc at the deflector redshift of z = 0.581.

2. Observations and data reduction

2.1. HST imaging

We observed SDSS J0100+1818 with the HST Wide Field Camera 3 (WFC3) in June and August 2018 (program GO-15253; PI: R. Cañameras), during one orbit in each of the two F438W and F160W filters. We used astrodrizzle from the DrizzlePac software package (Fruchter 2010) to improve the correction for geometric distortion on the individual exposures calibrated with the standard HST pipeline. We also optimized both the rejection of cosmic rays with dedicated bad pixel masks, and the subtraction of the local sky background. Small regions with decreased sensitivity in the IR channel of WFC3 (blobs) were corrected using flat fields from the most recent version of the IR blob monitoring (Sunnquist 2018). The individual exposures were then re-drizzled and combined with inverse-variance weighting. For the F160W band, we adopted a pixel scale 0.066″ pix−1 to adequately sample the point spread function (PSF), and a value final_pixfrac = 0.8 optimized for this pixel scale. For the F438W band, we used 0.033″ pix−1 and final_pixfrac = 0.6. The PSF full widths at half maximum (FWHMs) measured on isolated stars in the field are 0.086″ and 0.187″ in F438W and F160W, respectively.

2.2. NOT imaging

We used the 2.5 m NOT (program 56-032, PI: R. Cañameras), to obtain color information for galaxies in the lens environment and to characterize their spectral energy distributions (SEDs). In October 2017, we obtained StanCam imaging in the B and R bands with good or average weather conditions, and total exposure times of 2400 and 3600 s, respectively. We added a 900 s long V-band exposure obtained with ALFOSC in October 2013 through the fast track service (program 47-426, PI: C. Grillo). Image reduction was conducted with standard IRAF routines. We then used the Scamp and SWarp software (Bertin 2006, 2010) to correct for geometric distortion, resample individual frames and align them with respect to the WCS, using the USNO-B12 catalog as reference, and to obtain relative alignments with a rms accuracy below 0.1″. This resulted in a joint field-of-view of 2.5′ × 2.5′, pixel sizes of 0.19″, and PSF FWHMs of 1.05″, 0.70″, and 0.95″, in B, V, and R, respectively. Photometric calibration was performed with respect to SDSS, using color corrections from Jester et al. (2005), and refined by fitting blackbody functions to stars in the field. The AB magnitudes of the main lens galaxy in B, V, and R bands are 22.57 ± 0.06, 21.62 ± 0.10, and 20.22 ± 0.12, respectively.

2.3. VLT/X-shooter spectroscopy

SDSS J0100+1818 was observed with VLT/X-shooter (Vernet et al. 2011) (program 091.A-0852, PI: L. Christensen), in order to measure the lens and source spectroscopic redshifts, and to infer the lens stellar velocity dispersion. The main lens galaxy was observed in August and September 2013 with a 1.2″ wide slit, seeing ≃1″, clear sky conditions, and with a total on-source integration time of 3.2 h aimed at reaching a S/N of 5−10 per spectral bin in the VIS arm. The five individual OBs used a generic object, sky, object observing sequence with a 21″ offset. As part of the same program, a separate one-hour long OB targeted two multiple images of the brightest lensed source, by nodding between the two images with a 19″ offset.

After rejecting cosmic rays with Astro-SCRAPPY (McCully & Tewes 2019), the OBs were reduced with the X-shooter pipeline (Modigliani et al. 2010) using the EsoReflex environment version 2.8.3 (Freudling et al. 2013). The OB nodding between the two lensed images was reduced in stare mode in order to measure the sky background directly from the frames. We used the observations of telluric standard stars taken within two hours of the science observations at similar airmass, to correct for telluric absorption using the atmospheric synthetic spectra from Molecfit (Smette et al. 2015; Kausch et al. 2015). The combination of individual exposures, the subtraction of residual sky emission, and the optimization of the wavelength calibration were then performed with separate scripts (see Selsing et al. 2019). The 0.1−1″ offsets between the UVB and VIS arms resulting from the lack of atmospheric dispersion compensation at the time of the observations were also corrected on the 2D spectra. Finally, we conducted an optimal, profile-weighted extraction of the 1D spectrum of the lens galaxy of SDSS J0100+1818, and we extracted 1D spectra of multiple images using fixed apertures.

3. The SDSS J0100+1818 deflector

SDSS J0100+1818, (RA, dec) = (01:00:49.18, +18:18:27.79), was introduced in a late version of the CASSOWARY catalog, and was not included in the spectroscopic confirmation program from Stark et al. (2013). We therefore used the X-shooter spectra to measure secure redshifts for the main lens elliptical galaxy and lensed source. The best-fit lens redshift of z = 0.581 was inferred from the most prominent rest-frame optical absorption lines, as part of the stellar kinematic analysis detailed in Sect. 3.3. We inferred a joint redshift of z = 1.880 for the two source components forming image families A and B visible in Fig. 1. Given the low S/N values of the lines detected in the 2D spectra, this source redshift estimate is obtained from a joint analysis of the multiple images targeted by the X-shooter OBs. The width of the lines detected at about 10 740 Å in the binned 2D spectra of A1/B1, A3/B3, and A4/B4 is consistent with the [OII]λλ3727 doublet. Together with a faint detection of [OIII]λ5007 in A3/B3 (Appendix A) and with the lack of additional line detections over the spectral range covered by X-shooter, the detection of [OII]λλ3727 for the separate images ensures that the source redshift is robust. The 1D spectra of these multiple images are shown in the Appendix A.

thumbnail Fig. 1.

F160W image of the strong lensing system analyzed in this work. The multiple images corresponding to each of the three background sources considered in the analysis are labeled A1–A4, B1–B4, and C1–C2. Later in this study we show that sources A and B represent two emission peaks of a single extended source.

By combining observations in the HST F160W and F438W filters with the first strong lensing models of the system, we identified another candidate background source, whose two multiple images are labeled with C1 and C2 in Fig. 1. Even though we are lacking spectroscopic confirmation, we include source C in our strong lensing models, with its redshift value as a free parameter, because the observation of two images with similar colors, robustly predicted at the same positions by the strong lensing models makes its multiply imaged nature highly likely. However, without a spectroscopic redshift measurement for source C, we cannot determine whether AB or C is the most distant source, so we will not be able to perform a multiplane lensing analysis here, and thus the light emitted by each source will be deflected only by the total mass distribution of the main lens, and not by any other background source. Then, thanks to the available photometric data and lensing model predictions, we also hypothesize the presence of an additional background source, with four multiple images. These images are barely detected in the HST data and two of the multiple images appear angularly very close to candidate group members. Thus, the one-component mass models we develop in this study might not be able to properly reconstruct these images. For these reasons, we do not consider this extra background source in the strong lensing models presented in this work, where we focus on the multiple images A, B, and C. We postpone to a future study a more detailed analysis of this fourth multiple-image family, when new data will possibly become available (see the last point in Sect. 6).

In the following we present the physical properties of the main deflector, such as its stellar mass, luminosity profile and stellar kinematics, measured from the full photometric and spectroscopic data set. The last subsection is dedicated to the study of the lens environment, discussing the possible group nature of this system.

3.1. Stellar mass

The SED of the main lens galaxy is modeled with the Code Investigating GAlaxy Emission (CIGALE; Burgarella et al. 2005; Noll et al. 2009; Boquien et al. 2019). We use the multiband photometry from the Panoramic Survey Telescope & Rapid Response System (Pan-STARRS), NOT, and HST, taking the F160W flux corrected from contamination by neighboring sources. The magnitudes are measured with SExtractor (Bertin & Arnouts 1996), using ISOCORR down to 3σ isophotes in all bands for the main lens galaxy, except for the magnitude in HST F160W filter, measured from its luminosity profile, as described in the next subsection. The grid of models to conduct SED fitting relies on Bruzual & Charlot (2003) single stellar population templates, and assume delayed star formation histories with exponential bursts and e-folding times in the range 0.1–1 Gyr, ages between 0.5 and 8 Gyr, and the modified Charlot & Fall (2000) extinction law. We assume a Salpeter (Salpeter 1955) stellar IMF and fix the metallicity to solar values, as expected for massive early-type galaxies at z ∼ 0.6 (Sonnenfeld et al. 2015; Conroy et al. 2013; Gallazzi et al. 2014). The best-fit SED shown in Fig. 2 results in a total stellar mass value of (1.5 ± 0.3)×1012 M⊙.

thumbnail Fig. 2.

Spectral energy distribution of the main lens galaxy of SDSS J0100+1818. The observed flux densities from NOT (red triangles), Pan-STARRS (red circles), and HST (red square) are used to fit the SED with the CIGALE software. Error bars are smaller than the symbols. The best-fit stellar continuum template (blue curve) corresponds to the model fluxes plotted in orange and is used to infer the total stellar mass value of the deflector. The bottom panel shows the relative residuals of the fit.

3.2. Luminosity profile

We perform the photometric modeling of the system using the public software GALFIT (Peng et al. 2002, 2010). Considering the HST image in the F160W band, we extracted a cutout centered on the deflector with size of approximately 30″ × 30″ (see Fig. 3). Then, we decided to model the main lens and two other bright galaxies, located between the main lens and the multiple images/arcs. We masked out the external regions, where some bright objects are present, and the arcs, including the multiple image C2 (see the masked pixels, in white, in the right-hand panel of Fig. 3). In particular, C2 is too close to the main lens and radially distorted to be separately and adequately modeled. The very low values of the normalized residuals, without any clear pattern, confirm the goodness of our modeling choices.

thumbnail Fig. 3.

Results of the photometric modeling of the system with GALFIT. Left: considered cutout of the HST image in the F160W band. Center: best-fit model, consisting of a combination of two Sérsic profiles for the main lens in the middle, and a single Sérsic profile for each of the neighboring galaxies. Right: normalized residuals on a color scale from approximately −2σ to +2σ. The white pixels, with values of 0.0, show the regions masked during the luminosity profile modeling.

We assume Sérsic (Sérsic 1963) profiles for the galaxies, and a uniform distribution for the background. The parameter values are optimized by minimizing a standard χ2 function, defined as

χ 2 = ∑ i N pix ( f i data − f i model ) 2 σ i 2 , $$ \begin{aligned} \chi ^2 = \sum _i^{N_{\rm pix}} \frac{\left( f_{i}^\mathrm{data} -f_{i}^\mathrm{model} \right) ^2}{\sigma ^2_i} , \end{aligned} $$(1)

where the index i runs over all the non-masked pixels, up to Npix, f i data $ f_{i}^\mathrm{data} $ is the intensity value of the ith pixel in the input data image and f i model $ f_{i}^\mathrm{model} $ is the intensity value in the same pixel predicted by the model, after the convolution with the estimated PSF. The PSF we use was obtained by combining several non-saturated bright stars in the field-of-view, stacking them, subtracting the background, and then normalizing the resulting image.

To better reproduce the observations, and thus to reduce the residuals in the inner regions, the main lens galaxy has been modeled with a combination of two Sérsic profiles. The centers of the two components are not forced to coincide, but the best-fit centers differ by less than 1 pixel. In Fig. 3, the considered cutout is shown in the left panel, the best-fit model in the middle one, and the normalized residuals in the right one. The number of degrees of freedom is evaluated from the number of considered (non-masked) pixels Npix and the number of free parameters of the model Npar, as ndof =Npix − Npar. The minimum χ2 value is 13 899, which, with a number of degrees of freedom ndof = 46752, corresponds to a reduced χ2 value of ∼0.3.

We made use of the best-fit model of the main lens galaxy for three purposes: (1) We worked with lens-subtracted images in the strong gravitational lensing modeling, explained in the following sections, to avoid the lens light contaminating the multiple images and the arcs. This will be particularly important when considering source C, for which the multiple image C2 is angularly very close to the lens center.

(2) We built the cumulative luminosity profile of the deflector, from which we measure a total F160W magnitude value of mF160W = 17.06 ± 0.05, and an effective (i.e., half-light) radius value of (1.42 ± 0.02)″, corresponding to (9.32 ± 0.12) kpc at z = 0.581. The value of the total magnitude was considered in the SED fitting (see Fig. 2), from which we measured the total stellar mass value.

(3) We converted the luminosity profile into a stellar mass profile, shown in Fig. 4, by assuming a constant stellar mass-to-light-ratio. In detail, we considered the galaxy cumulative luminosity profile, normalized it with the total luminosity value, and then multiplied this adimensional profile by the value of the total stellar mass of the central galaxy, obtained from its SED fitting. The figure also shows the 1σ uncertainty and the value of the effective radius.

thumbnail Fig. 4.

Cumulative stellar mass (M*) profile (solid blue line) and ±1σ uncertainties (shaded area) of the main lens. It is measured from the total stellar mass measured through the SED fitting and assuming a constant stellar mass-to-light ratio. The vertical black line shows the half-light radius of the lens galaxy.

3.3. Stellar kinematics

We used the penalized Pixel Fitting software (pPXF; Cappellari & Emsellem 2004; Cappellari 2017) that follows a maximum penalized likelihood method to determine the stellar velocity dispersion σ* of the foreground lens galaxy. The X-shooter spectrum has a resolution of ≃6500 in the VIS arm, where the main absorption lines fall. We masked the strong telluric absorption windows that could affect the analysis, and re-binned by two spectral pixels in order to obtain a S/N > 10 per bin and ensure reliable measurements of the stellar kinematics. We used the SYNTHE high-resolution template library (Munari et al. 2005) that covers the entire 2500–10 500 Å wavelength range, with 1 Å pix−1 sampling and a constant resolution FWHM of 2 Å pix−1. Given the large size of the library, we extracted a subset of spectra representing the overall range of effective temperatures. The template spectra were then matched to the resolution of the binned X-shooter spectrum. The pPXF fit was conducted in the rest-frame wavelength range 3800–5300 Å that covers the main absorption lines, Ca H and K at 3935 and 3970 Å, the G band at 4305 Å, Hβ at 4863 Å, the MgI triplet at 5167, 5173, and 5184 Å, as well as the strong continuum break at 4000 Å. This results in a stellar velocity dispersion value of (451 ± 37) km s−1, confirming that the main lens galaxy of SDSS J0100+1818 is among the rarest, most massive elliptical galaxies (Loeb & Peebles 2003).

3.4. Environment

Given the large angular separation between the multiple images, we must consider the possibility that the main lens galaxy is located within a rich and overdense environment, as for other high-mass lens early-type galaxies at intermediate redshift (e.g., Newman et al. 2015; Johnson et al. 2018; Wang et al. 2022). We estimate photometric redshifts using the template-fitting BPZ package (Benítez 2000) to study the source distribution over the 2.5′ × 2.5′ field-of-view covered by our multiband images. In particular, we use the grizy bands from Pan-STARRS and BVR from NOT. At this stage, to avoid the possible introduction of systematics, we do not include the F160W band from HST, which only covers the central part of the considered field. The combination of optical and near-IR bands cover the 4000 Å break for early-type galaxies up to z ∼ 1. We obtain redshift estimates for a total of 142 galaxies over the field, with a typical uncertainty of ±0.2, and 53 of the most reliable redshifts are consistent at the 2σ level with z = 0.581. Figure 5 shows the difference between the redshift probability distribution functions of sources within 0.5′ of the main lens (about 200 kpc at z = 0.581), and over the rest of the field, both normalized to the same area. Our analysis is broadly consistent with the presence of an overdensity around the main deflector, and possible galaxy members mostly have R-band AB magnitudes of 22–24 mag. Definite confirmation of this structure nonetheless requires spectroscopic follow-up and, given the large spread of candidate members over the field, we simply model its effect on the light deflection with an external shear component.

thumbnail Fig. 5.

Difference between the redshift distribution functions of galaxies located within 0.5′ of the main ultra-massive lens galaxy at z = 0.581, and over the rest of the 2.5′ × 2.5′ field, after normalizing to the same area (orange curve). The black histogram is obtained from the peak of the redshift distribution of each source. A total of 14 galaxies within 0.5′ of the main deflector have redshifts consistent, at the 2σ confidence level, with z = 0.581, and 39 over the rest of our analyzed field.

4. Point-like source modeling

We used the Gravitational Lens Efficient Explorer (GLEE; Suyu & Halkola 2010; Suyu et al. 2012) software to perform our strong lensing modeling of the system. GLEE supports several types of mass and light profiles and uses Bayesian analyses, such as simulated annealing and Markov chain Monte Carlo (MCMC), to infer the probability distributions of the parameter values and thus characterize the deflector. It also employs the Emcee package, developed by Foreman-Mackey et al. (2013), for sampling the model parameter posterior.

We started our modeling by making use of point-like images, that is, considering each multiple image position as that of its brightest pixel, with an uncertainty of one pixel. The best-fit values of the parameters of each model were estimated through a multistep simulated annealing technique, first minimizing the χ2 on the source plane, and then on the deflector plane. Furthermore, the median values and the uncertainties (mainly at the 68% confidence level) were extracted from MCMC chains of 106 total steps, with acceptance rates between 20% and 30%, and the first 10% of the burn-in steps rejected. The final results are extracted from the final chains of a sequence, in which each intermediate chain is used to estimate the covariance matrix of the model parameters and to extract the starting point for the following one.

We model the total mass distribution of the lens by assuming different mass profiles. The goodness of each model is evaluated by comparing the number of degrees of freedom (d.o.f.) and the value of the minimum χ2, computed by comparing the observed, θ i obs $ {\boldsymbol{\theta}}_i^\mathrm{obs} $, and model-predicted, θ i pred $ {\boldsymbol{\theta}}_i^\mathrm{pred} $, multiple image positions with the following statistics:

χ 2 = ∑ i = 0 N | θ i obs − θ i pred | 2 σ i 2 , $$ \begin{aligned} \chi ^2 = \sum _{i=0}^N \frac{|\boldsymbol{\theta }_i^\mathrm{obs} -\boldsymbol{\theta }_i^\mathrm{pred} |^2}{\sigma ^2_i} , \end{aligned} $$(2)

where N is the total number of multiple images and σi is the positional uncertainty relative to the ith image. To compare different models, we considered three statistical estimators, often employed in similar strong lensing studies (see, e.g., Acebron et al. 2017; Mahler et al. 2018; Caminha et al. 2022): (1) the root-mean-square (rms) between the observed and model-predicted image positions, defined as

rms = 1 N | θ i obs − θ i pred | 2 ; $$ \begin{aligned} \mathrm{rms} = \sqrt{\frac{1}{N}|\boldsymbol{\theta }_i^\mathrm{obs} -\boldsymbol{\theta }_i^\mathrm{pred} |^2} ; \end{aligned} $$(3)

(2) the Bayesian information criterion (BIC; Schwarz 1978), given by

BIC = k ln n + χ 2 , $$ \begin{aligned} \mathrm{BIC} = k \ln {n} + \chi ^2 , \end{aligned} $$(4)

where k is the number of free parameters; and (3) the corrected Akaike information criterion (AICc; Akaike 1974; Cavanaugh 1997), defined as

AICc = 2 k + χ 2 + 2 k ( k + 1 ) n − k − 1 . $$ \begin{aligned} \mathrm{AICc} = 2 k + \chi ^2 + \frac{2k(k+1)}{n-k-1} . \end{aligned} $$(5)

The BIC and AICc estimators penalize models with increasing number of free parameters, to contrast overfitting. Thus, models with lower BIC and AICc values are preferred.

4.1. Modeling with A and B

At the beginning, we consider in the strong lensing modeling only the positions of the spectroscopically-confirmed eight multiple images of sources A and B.

First, we assume a pseudo-isothermal total mass distribution (PIEMD; Kneib et al. 1996) for the main lens. In GLEE, it is described by six parameters: the x and y coordinates of the center, the semimajor (a) to semiminor (b) axis ratio q = b/a, the position angle θ, the Einstein radius θE, and the core radius rcore. Here, we fix rcore = 0. We note that the value of θE is defined for a source at z = ∞ and does not correspond to that of the Einstein radius of the system, which should be nearly independent of the mass modeling details. The value of θE is a parameter that describes the lens strength and enters the dimensionless surface mass density κPIEMD as

κ PIEMD ( x , y ) = θ E 2 r core 2 + x 2 ( 1 + e ) 2 + y 2 ( 1 − e ) 2 , $$ \begin{aligned} \kappa _\mathrm{PIEMD} (x,{ y}) = \frac{\theta _{\rm E}}{2 \sqrt{r_{\rm core}^2+\frac{x^2}{(1+e)^2}+\frac{{ y}^2}{(1-e)^2}}} , \end{aligned} $$(6)

where the ellipticity e = 1 + q 1 − q $ e=\frac{1+q}{1-q} $. The value of the physical Einstein radius, as well as being estimated from this equation, would be clearly estimated from the total mass profiles presented in the following.

Then, we assume a singular power law elliptical mass distribution (SPEMD; Barkana 1998). In GLEE, it is described by seven parameters: the first six are in common with the PIEMD profile, and the slope g, which is related to the 3D logarithmic density slope γ′=dlog[ρ(r)]/dlog(r) (i.e., ρ∝ r − γ ′ $ \rho \propto r^{-\gamma^\prime} $) through γ′ = 2g + 1 (i.e., an isothermal profile corresponds to γ′ = 2 and g = 0.5). In the following we refer to the physical parameter γ′. Similarly to the PIEMD case, θE is a parameter of the mass distribution introduced in the dimensionless surface mass density κSPEMD as

κ SPEMD ( x , y ) = θ E ( x 2 + y 2 q 2 + 4 r core 2 ( 1 + q ) 2 ) − γ ′ − 1 2 . $$ \begin{aligned} \kappa _\mathrm{SPEMD} (x,{ y}) = \theta _{\rm E} \left( x^2+\frac{{ y}^2}{q^2} + \frac{4r^2_{\rm core}}{(1+q)^2}\right) ^ {-\frac{\gamma ^\prime -1}{2}}. \end{aligned} $$(7)

In both total mass distributions, we included an external shear component. In GLEE, this is described by two parameters: the shear strength, γext, and its position angle, ϕext (where 0 means that images are stretched horizontally, along the x-axis, and π/2 means images are stretched vertically, along the y-axis).

The best-fit values of the parameters for these first two models are reported in the two upper entries of Table 1. The values of the x and y coordinates are in arcseconds, relative to the center of light of the elliptical lens galaxy, identified as its brightest pixel. For each model, we list the best-fit parameter values, the number of observables (Nobs), the number of degrees of freedom (d.o.f.), and the value of the minimum ( χ min 2 $ \chi^2_{\min} $).

Table 1.

Modeling of the deflector total mass distribution using the observed multiple image positions.

4.2. Modeling with A, B, and C

Then, we add the two multiple images of source C, and use the same two mass density profiles defined in the previous section to fit the positions of all ten multiple images. The introduction of C requires an additional free parameter, its unknown redshift. In GLEE, this is parametrized with the value of Dds/Ds, which is the ratio between the distances between the deflector and the source and between the observer and the source. Below, we convert the values (and the posterior probability distribution) of Dds/Ds into those corresponding to the redshift of the source, zC, and we always refer to this parameter. The best-fit parameter values are reported in the third and fourth entries of Table 1.

Since the position of the multiple image C2 is angularly very close to the main lens galaxy, we decide to try also PIEMD and SPEMD mass density profiles with the value of rcore free to vary. We refer to these models as PIEMD+rc and SPEMD+rc. We expect that C2 can affect the reconstructed lens total mass distribution in the inner region, depending mainly on the values of rcore and γ′. In fact, from strong lensing data, one can accurately measure the total mass distribution of a lens at the radii where the multiple images are observed. The presence of source C allows us to expand the radial interval from approximately (30–50) to (15–63) kpc. The best-fit values of the parameters of the models PIEMD+rc and SPEMD+rc are reported in the two lower entries of Table 1.

By comparing the best-fit parameter values of these different models, we notice that the center of the total mass is always shifted with respect to the luminosity center, along the positive x-axis and the negative y-axis. Furthermore, the values of the angles θ and ϕext are similar and approximately orthogonal in almost all the models. When source C is included, the total mass of the deflector becomes more elliptical (i.e., the value of q decreases), unless a nonvanishing core radius is considered. Moreover, the introduction of source C and of a possible core radius changes significantly the best-fit values of the Einstein radii, suggesting a degeneracy between the values of the θE, rcore, γ′, and zC parameters. Finally, the best-fit SPEMD profiles without a core radius are shallower than the corresponding isothermal profiles (γ′≈1.5 vs. γ′ = 2); the best-fit values of the core radius and slope are very similar when the former is allowed to vary.

The low values of the minimum χ2 suggest that we can accurately reproduce the positions of the observed multiple images by assuming a one-component total mass distribution for the deflector. In fact, the average distance between the observed and model-predicted positions of the multiple images is only of 0.02″. In Fig. 6, we show the observed (filled circles) and the model-predicted (crosses) positions of the multiple images for the PIEMD+rc model. All the considered models would result in very similar figures.

thumbnail Fig. 6.

HST/F160W image of SDSS J0100+1818 and comparison between the observed multiple images of sources A (blue), B (red), and C (green), represented with circles, and those predicted by the best-fit PIEMD+rc model (fifth entry of Table 1), represented with crosses in the corresponding colors.

Furthermore, on the top of Table 4 we report the median and the 68% confidence level uncertainty values of the parameters of the different models. They are extracted from MCMC posteriors with 106 steps, with acceptance rates between 20% and 30%, after excluding the first 10% steps (considered as the burn-in phase). To correctly compare different models with different numbers of observables and degrees of freedom, for each model we rescaled the errors on the observed multiple images so that the value of χ min 2 $ \chi^2_{\min} $ is approximately equal to that of d.o.f..

In Fig. 7, we show the posterior probability density distributions of the parameters of the PIEMD+rc (in blue) and SPEMD (in red) models, with sources A, B, and C included in the modeling. In order to compare them, only the common parameters are plotted, and thus γ′ and rcore are not included. The x, y, q and θ distributions of the PIEMD+rc and SPEMD models are consistent. The Einstein radius distributions are centered on different values, but this is not surprising: its value is degenerate with those of the slope γ′ (see Fig. 11) and of the core radius rcore. Finally, we remark that the SPEMD model shows a unimodal probability density distribution for zC, as opposed to the bimodal one of the PIEMD+rc model.

thumbnail Fig. 7.

Probability density distributions of the parameters of the PIEMD+rc (in blue) and SPEMD (in red) models. The marginalized 1D histograms of each parameter are shown along the diagonal, while the other panels show the joint 2D probability distributions of the two parameters reported on the horizontal and vertical axes. The parameters plotted are those introduced in Sect. 4.1. The vertical dashed lines in the 1D histograms represent the 16th, 50th, and 84th percentiles, while the solid lines in the 2D distributions represent the 0.68 and 0.95 contour levels.

Finally, we show in Fig. 8 the comparison between two models with the same total mass parameterization (i.e., SPEMD), optimized first with the eight multiple images of sources A and B (in purple) and then with the ten multiple images of sources A, B, and C (in green). We can conclude that the results obtained with only two sources are consistent with those obtained with three sources and that the introduction of source C, at a different redshift, significantly reduces all but one confidence intervals.

thumbnail Fig. 8.

Probability density distribution of the parameters of the SPEMD model considering the multiple image positions of sources A and B only (in purple) and of sources A, B, and C (in green). The panels show the probability distributions as described in Fig. 7.

4.3. Total mass profile of the deflector

For each model, we estimate the total mass (MT) distribution of the deflector by extracting randomly 1000 models from the 106 steps of the last MCMCs described above. First, we convert the convergence maps provided by GLEE into the corresponding total mass maps. Then, we sum up the contribution of all the pixels within circular apertures centered on the brightest pixel of the main lens galaxy, with a step of 0.5 pixels, to obtain the cumulative total mass profiles that are presented in the following. Within each fixed aperture, we consider the distribution of the total mass values of all the 1000 random models, from which we measure the median value and the 16th and 84th percentiles, as the uncertainties at the 1σ confidence level.

We started comparing the cumulative total mass profiles of each different mass model reconstructed by using the multiple images of sources A and B only and A, B, and C. We find that the profiles show consistent median values and that the introduction of source C significantly reduces the statistical uncertainties in the inner and outer regions. For example, for the SPEMD model, when source C is considered, the relative uncertainties are reduced from about 32% to 17% and from 9% to 5% in the inner (R ≈ 15 kpc) and outer (R ≈ 63 kpc) regions, respectively. Thus, in the following discussion and figures, we refer only to the models optimized by taking the observed multiple image positions of sources A, B, and C into account.

In Fig. 9 we show the cumulative projected total mass profiles of the deflector for the PIEMD+rc and SPEMD models. Despite the different total mass density parameterizations, and the slightly different optimized values of the parameters, the profiles are very similar and agree on the same mass value, with the smallest statistical errors, at the mean distance (R ≈ 42 kpc) of the multiple images from the lens center, which represents approximately the physical Einstein radius of the system. The relative uncertainties at the distances of the innermost (C2, R ≈ 15 kpc) and outermost (C1, R ≈ 63 kpc) multiple images are, respectively, about 17% and 5% for the SPEMD model, while its minimum of approximately 1% is reached at R ≈ 42 kpc. The total projected mass value measured within this radius is (9.06 ± 0.04)×1012 M⊙ for the PIEMD+rc model and of (9.09 ± 0.05)×1012 M⊙ with the SPEMD model. In the same figure, we also include the stellar mass profile of the system we described above.

thumbnail Fig. 9.

Cumulative projected total mass profiles for the PIEMD+rc (blue) and SPEMD (red) models with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C as point-like objects. The cumulative projected stellar mass profile from Fig. 4 is shown in dashed black. The vertical lines close to the x-axis locate the distances from the lens galaxy center of the different multiple images, color-coded following Figs. 1 and 6. The black line shows the effective radius of the main lens galaxy.

Finally, in Fig. 10, we plot the cumulative projected stellar-over-total mass fraction profiles for the PIEMD+rc and SPEMD models. The two profiles are consistent, given the uncertainties, and differ mainly in the inner regions. At the lens galaxy effective radius, we estimate values of (70 ± 29)% and (49 ± 12)% for the PIEMD+rc and SPEMD models, respectively. According to both models, the fraction values are of (15 ± 3)% and (10 ± 2)% at R ≈ 42 kpc and R ≈ 63 kpc, respectively.

thumbnail Fig. 10.

Cumulative projected stellar-over-total mass fraction profiles for the PIEMD+rc (blue) and SPEMD (red) models with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C as point-like objects. The vertical lines are the same as in Fig. 9.

5. Extended source modeling

In this section we describe the reconstruction of the extended SB of the lensed sources that we performed with GLEE. At this stage, we merge sources A and B into a single, double-peaked source, at z = 1.88. Thus, in the following, we refer to “single source” models as the analogous of modeling with A and B in the point-like source modeling, and to “two source” models as of modeling with A, B, and C.

To model the SB of extended images, GLEE needs as input a mask with the pixels containing the observed multiple images and possible arcs of the lensed sources. These pixels are mapped onto the source plane through the ray-tracing equation, where the source SB is reconstructed on a regular grid of pixels. The ray-tracing equation depends on the total mass distribution of the lens and on the angular diameter distances between the lens and the sources. When we model the observed images of a single source (i.e., A and B), with its spectroscopic redshift value measured and considered fixed, GLEE finds the best-fit model by varying the values of the parameters of the mass distribution of the lens, and optimizing the posterior probability distribution, which is proportional to the product of the likelihood and the prior of the lens mass parameters (Suyu & Halkola 2010). Instead, when we model the images of the two sources, the mapping also depends on the relative distances between the observer, the deflector, and source C, and thus its redshift is introduced as an additional free parameter.

In the optimization process, the pixelated SB distributions of the sources are mapped back onto the image plane to obtain the model-predicted images, with flux intensity f ij pred $ f_{ij}^\mathrm{pred} $, and compared with those observed, f ij obs $ f_{ij}^\mathrm{obs} $. The χ SB 2 $ \chi^2_\mathrm{SB} $ is thus defined as

χ SB 2 = ∑ i , j ( f ij obs − f ij pred ) 2 σ ij 2 , $$ \begin{aligned} \chi ^2_\mathrm{SB} = \sum _{i,j} \frac{\left( f_{ij}^\mathrm{obs} -f_{ij}^\mathrm{pred} \right) ^2}{\sigma ^2_{ij}} , \end{aligned} $$(8)

where i and j identify the position of each pixel. During the optimization, we impose a curvature regularization on the values of the source SB pixels, weighted by the value of the parameter λ ̂ $ \hat{\lambda} $, which is included in the χ SB 2 $ \chi^2_\mathrm{SB} $ evaluation, as described in Suyu et al. (2006). This parameter penalizes physically less plausible solutions with very irregular SB distributions for the sources. The number of degrees of freedom is here computed as

d . o . f . = n mask − n lens − n source , $$ \begin{aligned} \mathrm{d.o.f.} = n_\mathrm{mask} - n_\mathrm{lens} - n_\mathrm{source} , \end{aligned} $$(9)

where nmask is the number of pixels on the image plane included in the mask as observables, nlens is the number of free parameters related to the total mass distribution of the lens, and nsource is the effective number (i.e., smaller than the adopted one) of pixels onto which the source is reconstructed on the source plane. At the beginning, we consider the SB of a source on a regular grid of 25 × 25 pixels. Then, we refine (and show here) the results adopting a 40 × 40 grid on the source plane. We note that GLEE, after mapping the pixels included in the mask onto the source plane, reconstructs the source SB on a regular grid of pixels, with dimensions chosen by the user. Given that the maximum physical distance of the mapped pixels along the two coordinate axes is not known a priori, the SB source reconstructions usually have different pixel scales on the x and y axes, as can be observed in the plots presented in the following.

5.1. Modeling the deflector and the arcs

We started from the best-fit models achieved with the point-like source modeling, described in Sect. 4. In Table 1, we observe that the χ min 2 $ \chi^2_{\min} $ value decreases significantly, from approximately 4 to 1, when the number of degrees of freedom is reduced from 6 to 5 (from the PIEMD to the SPEMD and PIEMD+rc models). On the contrary, when d.o.f. is further reduced to 4 (for the SPEMD+rc model), the value of χ min 2 $ \chi^2_{\min} $ shows only a minor decrement. Furthermore, the SPEMD+rc model has the largest values for both the BIC and AICc estimators. This suggests that this model is the least favored one and, thus, that its extra complexity is statistically not fully justified. For these reasons, in the extended source modeling with two sources (AB and C), we decide to consider only the SPEMD and PIEMD+rc lens mass models, both with d.o.f. = 5 in the point-like source approximation. The free parameters of each mass distribution are the same as described previously. In Table 2, we report the best-fit values from the χ SB 2 $ \chi^2_\mathrm{SB} $ minimization with a 40 × 40 pixelated source. As before, we can observe some trends in the values of the reconstructed parameters. The center of the total mass distribution is slightly shifted from the light peak, along the positive x and negative y directions. The value of the magnitude of the external shear is high in the PIEMD model with a single source, but it significantly decreases when considering a SPEMD model or when source C is introduced. The value of θ is approximately the same for all models, while the value of q decreases for the SPEMD models, where the deflector is found to be more elliptical. Finally, both the SPEMD models suggest a lens total mass profile shallower than an isothermal one, and the value of γ′ is almost identical, when a single or two background sources are considered.

Table 2.

Modeling of the deflector total mass distribution using the extended source modeling.

A comparison with the best-fit values of the point-like models shows that the center of the total mass distribution is always shifted in the same direction from the light distribution center. Moreover, the values of axis ratio, position angle, Einstein and core radii do not differ significantly when comparing each extended model with its corresponding point-like model.

The main differences for the lens reconstruction between the two ways of modeling the background sources can be observed in Table 4, where we list the median values and the statistical uncertainties extracted from MCMCs of 106 steps, tuned to have χ min 2 $ \chi^2_{\min} $ ≈ d.o.f. The larger number of observables (approximately 7200, instead of 20) makes the marginalized probability distributions of the parameters much narrower when the sources are modeled as extended, as visible from the uncertainty values reported in the table. We show in Fig. 11 the probability density distribution of the SPEMD model with two sources. There, we can observe the expected degeneracies between the values of θE, γ′ and zC, found also in the point-like modeling. The values of the shear parameters show some significant degeneracies with those of several other parameters, like q and θ.

thumbnail Fig. 11.

Probability density distribution for the SPEMD model with two sources (source AB and source C). The marginalized 1D histograms of each parameter are shown along the diagonal, while the other panels show the joint 2D probability distributions of the two parameters reported on the x and y axes. The parameters are those described in Sect. 4.1. The vertical dashed lines in the 1D histograms represent the 16th, 50th and 84th percentiles, which are also reported on top of each 1D histogram, while the solid lines and shaded areas in the 2D distributions represent the 0.68, 0.95, and 0.99 contour levels. In the floating panel, we zoom in the marginalized 1D histogram of the γ′ parameter (in black), and compare it with its distribution for the model with the AB source only (in red). The introduction of source C provides tighter constraints on the value of the γ′ parameter.

Similarly to Fig. 6, we show the comparison between the observations and the model predictions for the extended modeling. It is now possible to compare directly, pixel by pixel, the observed and predicted arcs. We show the results for the PIEMD (Fig. 12) and SPEMD (Fig. 13) models when considering the single source AB, and for the SPEMD model when including also source C (Fig. 14). In the first panels from the left, we show the pixels of the observed images included in the masks as observables. In the second panels, we plot the extended images predicted by the best-fit models and, in the third panels, the normalized residuals in the (−5σ, +5σ) interval. When a single source is reconstructed, we observe that the arcs are very well reproduced with the PIEMD model and the residuals do not exceed in absolute value 2.5σ. The results are even better with the SPEMD model, which is able to further reduce the residuals, in particular of the A2B2 and of A4B4 images. Besides the expected increase of the residuals (up to about 3.5σ around the images A3B3) when source C is added, the model can predict remarkably well the radial arc (C1). As already discussed, this multiple image is the most sensitive to the deflector inner total mass density profile and its good reconstruction supports the reliability of the strong lensing modeling results and of the lens SB modeling and subtraction. In the panels on the right, we show the reconstructed SB of the sources, which are discussed in Sect. 5.3.

thumbnail Fig. 12.

PIEMD model with one source. From the left to the right: observed SB in the F160W band of the multiple images considered in the extended source modeling; model-predicted SB; normalized residuals in the range from −5σ to +5σ; reconstructed SB of the source on a 40 × 40 pixel grid. Angular scales of 1″ and of 0.5″ are represented on the lens and source plane, respectively.

thumbnail Fig. 13.

As in Fig. 12 for the SPEMD model with one source.

thumbnail Fig. 14.

As in Fig. 12 for the SPEMD model with two sources. The top row shows sources A and B and their SB is reconstructed on a 40 × 40 pixel grid, while the bottom row shows source C and its SB is reconstructed on a 25 × 25 pixel grid.

5.2. Total mass profiles of the deflector

We measure the lens projected total mass distribution from the extended source modeling using the same procedure described for the point-like source modeling. In Fig. 15, we show the cumulative projected total mass profile of the deflector measured with the PIEMD+rc model and two extended sources (solid red line, the 1σ uncertainties are smaller than the linewidth), compared to its corresponding point-like source model (solid black and 1σ uncertainties as shaded area). The very small statistical uncertainties on the lens total mass profile from the extended source modeling, which are a factor of ≳10 smaller than those presented in Sect. 4.3, are related to the small uncertainties on the values of all the model parameters discussed above. In addition to the statistical uncertainties, we can estimate the systematic uncertainties by comparing the results of the different models, which fit the data almost equally well. Within R ≈ 42 kpc, we measure a projected total mass value of 8.76, 9.00, 9.08, and 9.08 × 1012 M⊙ for the PIEMD and SPEMD (single source), SPEMD and PIEMD+rc (two sources) models, respectively.

thumbnail Fig. 15.

Cumulative projected total mass profiles for the PIEMD+rc model with point-like (solid black) and extended source (red) modeling, with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C. For the extended source modeling, the uncertainties are smaller than the linewidth. The vertical lines close to the x-axis locate the distances from the lens galaxy center of the different multiple images, color-coded following Figs. 1 and 6. The black line shows the effective radius of the main lens galaxy.

On the contrary, the statistical uncertainties on the cumulative projected stellar-over-total mass fraction profiles for corresponding point-like and extended models are very similar, because the uncertainty on the stellar mass profile dominates over that on the total mass profiles. In Fig. 16, we show the profiles for the PIEMD+rc models. They present slightly different slopes, but they both result in a stellar-over-total mass fraction value of (15 ± 3)% within 42 kpc. From a comparison with the PIEMD+rc point-like model, we find that the extended one predicts a higher value of (108 ± 22)% consistent with one, at the effective radius, and a lower value of (9 ± 2)% in the outermost regions (R ≈ 63 kpc).

thumbnail Fig. 16.

Cumulative projected stellar-over-total mass fraction profiles for the PIEMD+rc model with point-like (black) and extended source (red) modeling, with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C. The vertical lines are the same as in Fig. 9.

5.3. Source SB reconstruction

The source SB reconstructions are shown in the fourth panels of Figs. 12–14. In each image, there are two different rods along the x and y directions, representing a scale of 0.5″, as discussed above. By observing them, we note that, despite recovering the same morphology, the sizes of the background sources vary by a factor of ∼2, when the different total mass profiles for the lens are assumed.

We investigate this difference by measuring the angular separation between the sources A and B, at redshift z = 1.880, and by estimating the half-light radii of the sources A, B, and C. In the first case, we compare the angular separations predicted by the different point-like and extended models. In particular, the positions of the point-like sources are optimized directly by GLEE, while, when considered as extended, we measure the angular separation between the brightest pixels associated with the sources. The results are shown in Table 3. First, we observe that each model predicts approximately the same angular distance, when the point-like or extended modeling is adopted. This test further proves the consistency of the two methods. Then, we confirm quantitatively what we observed from the reconstructed images: both the PIEMD and SPEMD models can reproduce well the observed multiple images, but they predict angular separations between the sources A and B that differ by a factor of approximately 2. This shows that the choice on the lens total mass distribution has a significant impact on the inferred properties of the background sources. This effect might be very relevant for lens clusters, where hundreds of background lensed sources are detected and the total mass distribution of the lens is more complex. However, in this strong lensing system, even if both the PIEMD and SPEMD models can reproduce successfully the observed multiple images and have approximately the same number of degrees of freedom, the SPEMD model results in a lower χ min 2 $ \chi^2_{\min} $ value. We remark that the PIEMD+rc and SPEMD models, with sources A, B, and C, have the lowest χ min 2 $ \chi^2_{\min} $ values and predict consistent angular distances. Considering the limited number of different models that can be tested in a strong lensing analysis, we highlight that this effect should be quantified and quoted as a possible source of systematic uncertainty on the reconstructed sources.

Table 3.

Angular separation and sizes of the reconstructed sources.

Table 4.

Median and 68% confidence level uncertainty values of the parameters of the different models, based on the point-like (top) and extended (bottom) modeling of the sources, as described in Tables 1 and 2.

Next, we measured the effective radii of the two components A and B of the first source, and of source C. We computed the luminosity growth curve of each reconstructed source, measuring the flux included within concentric circular apertures, centered on the brightest pixel, with a step of 0.5 pixels (which corresponds to approximately 0.03″ for all the sources). The apertures were corrected to take into account the different pixel sizes along the x- and y-axes. For the peaks of the composite AB source, we consider semicircular apertures, to avoid the contribution of source B being included in the measurement of the half-light radius of source A, and vice versa. Hence, we build the luminosity cumulative profiles and infer the half-light radii by considering only the halves of the plane less affected by the other peaks.

We observe that the cumulative luminosity profiles do not converge to a plateau, because of the noise in which the reconstructed sources are immersed. Because of that, we had to choose a maximum radius for the two sources, on the visual inspection of the images, of 4 and 3 pixels for the sources A and B, respectively. For source C instead, the cumulative luminosity profile converges to a plateau. The results are summarized on the left part in Table 3.

We check our results with another test. We analyze the multiple images A1 and B1, which are the less distorted and better resolved, and consider a rectangular cutout of 2.1″ × 1.3″ around them. We fit the SB of the background arc with an extended Sérsic component, and the SB of the two images A1 and B1 with two additional Sérsic profiles. From the best-fit model, we measure the values of the effective radii re A1 = 0.23″ and re B1 = 0.17″. Then, we estimate the local magnification factor around the multiple images A1 and B1 to infer their intrinsic sizes on the source plane. For each model, we compute the median magnification map from 1000 models randomly extracted from the final MCMCs (see Sect. 4.3). For each pixel of coordinates (i, j) of each model, we compute the local magnification factor as

μ i , j = 1 ( 1 − κ i , j ) 2 − γ 1 i , j 2 − γ 2 i , j 2 , $$ \begin{aligned} \mu _{i,j} = \dfrac{1}{(1-\kappa _{i,j})^2-\gamma _{1\, i,j}^2 - \gamma _{2\, i,j}^2} , \end{aligned} $$(10)

where κi, j, γ1 i, j and γ2 i, j are the values of the convergence and shear components, respectively. Then, from the 1000 values of magnification estimated on each pixel, we consider the median value and build the median magnification map μ ¯ i , j $ \overline{\mu}_{i,j} $. With the 16th and 84th percentile values of the same distribution, we quantify the ±1σ uncertainties. Finally, we measure the local magnification factors for the multiple images A1 and B1 by taking the median value within circles centered on them and with radii equal to re A1 and re B1. The measured local magnification factors around A1 and B1 are reported in Table 3. If we correct the observed effective radii on the lens plane according to

r e source = r e image μ ¯ , $$ \begin{aligned} r_{\rm e}^\mathrm{source} = \frac{r_{\rm e}^\mathrm{image} }{\sqrt{\overline{\mu }}} , \end{aligned} $$(11)

we find that these results (see the last two columns of Table 3) are consistent, given the uncertainties, with those obtained on the source plane with the luminosity growth curve of the SB reconstructions.

6. Discussion and summary

In this paper we have studied SDSS J0100+1818, a strong lensing system consisting of a massive early-type galaxy with a spectroscopic redshift of z = 0.581, which acts as a gravitational lens on two background sources, AB and C. One of the sources (AB), spectroscopically confirmed at z = 1.880, has four multiple images, visible around the deflector at a projected distance of about 7″, and presents two components. The other source (C), instead, does not have a spectroscopic redshift measurement, and its two multiple images represent the closest and the most distant images from the deflector center. Thus, the introduction of this source is key to the reconstruction of the total mass profile of the lens in its inner and outer regions, approximately from 15 to 63 kpc, and to the reduction of some degeneracies among the model parameters.

We have developed several strong lensing models of the deflector with the software GLEE, combining (cored or singular) PIEMD and SPEMD total mass profiles with or without an external shear component. At the beginning, we considered the sources as point-like objects and modeled the lensed source positions. We took advantage of the eight observed positions of source AB (four multiple images for each of the two components) and of the two observed positions of source C. The redshift of source C was also included as a free parameter in the models.

Then, we considered the sources as extended objects and reconstructed their SB distributions. With this improvement, we were able to exploit the image structure and the extended arcs in which they are distorted, over ∼7200 HST pixels. We finally used the reconstructed sources to measure their sizes and have discussed how much they can be affected by the adopted total mass profile for the deflector.

The main results of this work can be summarized as follows:

  • We combined the available multiband photometry from Pan-STARRS, NOT, and HST to model the SED of the main lens galaxy. The best-fit model results in a stellar mass value of (1.5 ± 0.3)×1012 M⊙. By using the public software GALFIT on the HST image in the F160W band of the system, we modeled the light distribution of the lens with a combination of two Sérsic profiles. Starting from them, we measured the cumulative luminosity profile and then converted it into a stellar mass profile by assuming a constant stellar mass-to-light-ratio.

  • We used the public software pPXF to estimate the value of the stellar velocity dispersion, σ, of the main lens galaxy from its X-shooter spectrum. We find σ = (451 ± 37) km s−1, which is consistent with the very large values of the galaxy stellar mass and of the mean distance between the observed multiple images, confirming that SDSS J0100+1818 is among the rarest, most massive elliptical galaxies (Loeb & Peebles 2003).

  • With the point-like source modeling, we have found a total mass value projected within the Einstein radius (of approximately 42 kpc) of (9.1 ± 0.1)×1012 M⊙, consistent for the PIEMD+rc and SPEMD models (Fig. 9). Source C is predicted at zC = 1.72 (1.98) for the PIEMD+rc (SPEMD) model (Table 1). The best-fit value of the logarithmic slope of the SPEMD model is shallower than for a singular isothermal profile. However, we remark the observation of the expected degeneracies between the values of the γ′, θE, and zC parameters (Fig. 7).

  • With the extended source modeling, we have confirmed a projected total mass value enclosed within the Einstein radius of 9.1 × 1012 M⊙ with both the PIEMD+rc and SPEMD models. Source C is predicted at zC = 1.69 (2.04) for the PIEMD+rc (SPEMD) model (Table 2). By considering the extended structure of the sources and of their multiple images, the number of observables increases to approximately 7200 (compared to the 20 observables of the point-like source modeling). As a consequence, the statistical uncertainties on the parameter values and on the derived quantities are strongly reduced (see Fig. 15).

  • In Dutton & Treu (2014) and Newman et al. (2015), the value of γtot was defined as the 3D mass-weighted mean value of the density slope within Re. They investigated how this value varies by considering a sample of 59 galaxy-scale lenses from the Sloan Lens ACS (SLACS; Auger et al. 2009, 2010) survey (⟨z⟩  =  0.20), 10 group-scale lenses (⟨z⟩  =  0.36), and 7 central galaxies of massive clusters (⟨z⟩ = 0.25; Newman et al. 2013b,a). We measured the value of γtot from the projected total mass profiles, M(R), presented in Sects. 4 and 5. Thus, we could only directly compare the results for the models with a constant slope, the SPEMD profiles with rcore = 0 in our study. We obtained a γtot value of 1 . 89 − 0.09 + 0.22 $ 1.89^{+0.22}_{-0.09} $ for the point-like source SPEMD model and 1 . 68 − 0.01 + 0.01 $ 1.68^{+0.01}_{-0.01} $ for the extended source SPEMD model, considering the A, B, and C background sources. In both cases, we observed that the introduction of source C slightly increases the value of γtot, by 0.03 and 0.05, respectively. These values are consistent with those found by Newman et al. (2015) for group-scale systems with Re ≈ 10 kpc, which lie between the values obtained for galaxy-scale (2.09 ± 0.03) and cluster-scale ( 1.18 ± 0 . 07 − 0.07 + 0.05 $ 1.18 \pm 0.07^{+0.05}_{-0.07} $) systems.

  • Newman et al. (2015) exploited the samples introduced above to compare the values of the projected stellar-over-total mass fraction within Re, defined as

    f ∗ , Salp = M ∗ ( R e ) M T ( R e ) , $$ \begin{aligned} f_{*, \, \mathrm{Salp} } = \frac{M_*(R_{\rm e})}{M_T(R_{\rm e})} , \end{aligned} $$(12)

    where the subscript clarifies that the stellar mass values were estimated assuming a Salpeter stellar IMF. They also showed that the value of the projected stellar-over-total mass fraction decreases with increasing Einstein radius and halo mass, and that galaxy-, group-, and cluster-scale systems populate different regions in these parameter spaces. They measured a mean projected stellar-over-total fraction value of 0.60 for galaxy-scale lenses, 0.17 for group-scale lenses, and 0.06 for cluster-scale lenses. In this work, for SDSS J0100+1818, we measured fractions of 0.70 ± 0.29 and 0.49 ± 0.12 for the PIEMD+rc and SPEMD models, respectively (see Fig. 10).

  • We have explored the possible presence of other group members and estimated photometric redshifts over the 2.5′ × 2.5′ field-of-view covered by our multiband images. It has not been possible to gather quantitative conclusions, due to a lack of spectroscopy. However, we qualitatively observed an overdensity of candidate members in the northeast direction, which corresponds to a shear position angle, ϕext, of about 20° (note the different orientation of the compass in Fig. 1 and of the GLEEx-axis). This is in relatively good agreement with the best-fit values of ϕext, which range from 16° to 71°. We note that the value of ϕext shows some expected degeneracies, mainly with the values of b/a, θ, and γext.

  • We reconstructed the SB distributions of the background sources and measured their half-light radii from the luminosity profiles. We have successfully recovered the two-peaked structure of the AB source, with a small physical separation between 1.4 and 2.9 kpc (at z = 1.880) from all the models (Table 3). We have found that, when considering different models with similar χ min 2 $ \chi^2_{\min} $/d.o.f. values, the sizes of the reconstructed sources can vary by a factor of about 2. These results can have important consequences on the strong lensing modeling of more complicated lens mass distributions (i.e., on cluster scales). However, we remark that in this study the models with a SPEMD or a PIEMD+rc lens total mass profile have the lowest χ min 2 $ \chi^2_{\min} $ values, and they predict reconstructed sources with consistent sizes (confirmed also with the extended source modeling).

  • We have measured values of the effective radius between 0.5 and 1 kpc at z = 1.880 for the A and B components, depending on the adopted model. Approximately 60% of z ∼ 2 galaxies show bright star-forming regions, dubbed “clumps” (e.g., Cowie et al. 1995; Elmegreen & Elmegreen 2005; Elmegreen et al. 2009, 2013; Guo et al. 2015, 2018; Zanella et al. 2019). Our measured sizes are in agreement with those of strongly lensed clumps in star-forming galaxies at z ∼ 1 − 3, with similar, moderately magnified sources (μ ∼ 3 − 10). Hydrodynamic simulations suggest that the measured clump sizes depend on the spatial resolution of the observations (e.g., Oklopčić et al. 2017; Behrendt et al. 2016, 2019; Tamburello et al. 2017; Faure et al. 2021), and thus that the measured size upper limit decreases with an increasing magnification factor (Meštrić et al. 2022; Claeyssens et al. 2023). In high amplification regimes, recent observations have been able to explore clumps with sizes of ∼ 100 pc (Dessauges-Zavadsky et al. 2019; Livermore et al. 2015; Cava et al. 2018), down to ≲ 20 pc in some extremely magnified cases (Rigby et al. 2017; Vanzella et al. 2020; Meštrić et al. 2022; Claeyssens et al. 2023; Messa et al. 2022). Hence, considering the measured magnification factors for the multiple images of the SDSS J0100+1818 system, it remains unclear whether A and B are monolithic, isolated clumps or blends of smaller clumps (or, in general, subcomponents) unresolved with HST.

This study could be improved in several ways with additional integral field spectroscopy. With these data, we could: (1) measure the redshift value of source C, which would be crucial for breaking the degeneracy between the values of θE and γ′. A measurement of the lens total mass enclosed within two (different) Einstein radii would indeed allow us to estimate the values of θE and γ′ with a precision of a few percent. The measured redshift of C would also allow us to perform a multiplane lensing analysis, in which the light emitted by source AB is also deflected by the total mass distribution of C (or vice versa, if zC > 1.88); (2) confirm or reject the hypothesis of the group nature of SDSS J0100+1818. If some neighbor galaxies were confirmed at redshifts similar to that of the main deflector, we would be able to include them individually in the strong lensing model (and not simply as an external shear component); (3) measure the velocity dispersion profile of the main lens galaxy and combine kinematics and strong lensing information. The previous points would also pave the way to the confirmation and inclusion of the additional background source we mentioned in Sect. 3. If confirmed, the introduction of another source, at a different redshift, strongly lensed into four multiple images, would greatly improve our strong lensing model and make SDSS J0100+1818 one of the few galaxy-scale systems known to date with three lensed sources at different redshifts (see, e.g., Collett & Smith 2020).

1

The ID of this system has changed in the CASSOWARY survey in the past years. It was called CSWA 15 when we first targeted it, but it is currently identified as CSWA 115. In the paper, to avoid confusion, we refer to it as SDSS J0100+1818.

2

http://tdc-www.harvard.edu/catalogs/ub1.html

Acknowledgments

This research is based on observations made with the NASA/ESA Hubble Space Telescope obtained from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. These observations are associated with program GO-15253. Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere under ESO programme 091.A-0852(A), and on observations with the Nordic Optical Telescope, operated by the Nordic Optical Telescope Scientific Association at the Observatorio del Roque de los Muchachos, La Palma, Spain, of the Instituto de Astrofisica de Canarias. The data presented were obtained with ALFOSC, which is provided by the Instituto de Astrofisica de Andalucia (IAA) under a joint agreement with the University of Copenhagen and NOTSA. This project is partially funded by PRIN-MIUR 2020SKSTHZ. RC and SHS thank the Max Planck Society for support through the Max Planck Research Group for SHS.

References

  1. Abadi, M. G., Navarro, J. F., Fardal, M., Babul, A., & Steinmetz, M. 2010, MNRAS, 407, 435 [NASA ADS] [CrossRef] [Google Scholar]
  2. Acebron, A., Jullo, E., Limousin, M., et al. 2017, MNRAS, 470, 1809 [Google Scholar]
  3. Akaike, H. 1974, IEEE Trans. Autom, Control, 19, 716 [NASA ADS] [CrossRef] [Google Scholar]
  4. Annunziatella, M., Bonamigo, M., Grillo, C., et al. 2017, ApJ, 851, 81 [NASA ADS] [CrossRef] [Google Scholar]
  5. Auger, M. W., Treu, T., Bolton, A. S., et al. 2009, ApJ, 705, 1099 [Google Scholar]
  6. Auger, M. W., Treu, T., Bolton, A. S., et al. 2010, ApJ, 724, 511 [NASA ADS] [CrossRef] [Google Scholar]
  7. Barkana, R. 1998, ApJ, 502, 531 [NASA ADS] [CrossRef] [Google Scholar]
  8. Barnabè, M., Czoske, O., Koopmans, L. V. E., Treu, T., & Bolton, A. S. 2011, MNRAS, 415, 2215 [Google Scholar]
  9. Barnabè, M., Spiniello, C., Koopmans, L. V. E., et al. 2013, MNRAS, 436, 253 [Google Scholar]
  10. Bartelmann, M. 2010, CQG, 27, 233001 [NASA ADS] [CrossRef] [Google Scholar]
  11. Behrendt, M., Burkert, A., & Schartmann, M. 2016, ApJ, 819, L2 [NASA ADS] [CrossRef] [Google Scholar]
  12. Behrendt, M., Schartmann, M., & Burkert, A. 2019, MNRAS, 488, 306 [NASA ADS] [CrossRef] [Google Scholar]
  13. Belokurov, V., Evans, N. W., Hewett, P. C., et al. 2009, MNRAS, 392, 104 [NASA ADS] [CrossRef] [Google Scholar]
  14. Benítez, N. 2000, ApJ, 536, 571 [Google Scholar]
  15. Bergamini, P., Rosati, P., Mercurio, A., et al. 2019, A&A, 631, A130 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. Bertin, E. 2006, in Astronomical Data Analysis Software and Systems XV, eds. C. Gabriel, C. Arviset, D. Ponz, & S. Enrique, ASP Conf. Ser., 351, 112 [Google Scholar]
  17. Bertin, E. 2010, Astrophysics Source Code Library [record ascl:1010.068] [Google Scholar]
  18. Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Blumenthal, G. R., Faber, S. M., Flores, R., & Primack, J. R. 1986, ApJ, 301, 27 [Google Scholar]
  20. Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A, 622, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
  22. Burgarella, D., Buat, V., & Iglesias-Páramo, J. 2005, MNRAS, 360, 1413 [NASA ADS] [CrossRef] [Google Scholar]
  23. Cañameras, R., Nesvadba, N. P. H., Kneissl, R., et al. 2017a, A&A, 600, L3 [Google Scholar]
  24. Cañameras, R., Nesvadba, N., Kneissl, R., et al. 2017b, A&A, 604, A117 [Google Scholar]
  25. Caminha, G. B., Grillo, C., Rosati, P., et al. 2016, A&A, 587, A80 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Caminha, G. B., Suyu, S. H., Grillo, C., & Rosati, P. 2022, A&A, 657, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Cao, S., Pan, Y., Biesiada, M., Godlowski, W., & Zhu, Z.-H. 2012, JCAP, 2012, 016 [CrossRef] [Google Scholar]
  28. Cappellari, M. 2017, MNRAS, 466, 798 [Google Scholar]
  29. Cappellari, M., & Emsellem, E. 2004, PASP, 116, 138 [Google Scholar]
  30. Cappellari, M., Bacon, R., Bureau, M., et al. 2006, MNRAS, 20, 127 [Google Scholar]
  31. Cava, A., Schaerer, D., Richard, J., et al. 2018, Nat. Astron., 2, 76 [Google Scholar]
  32. Cavanaugh, J. E. 1997, Stat. Probab. Lett., 33, 201 [CrossRef] [Google Scholar]
  33. Charlot, S., & Fall, S. M. 2000, ApJ, 539, 718 [Google Scholar]
  34. Claeyssens, A., Adamo, A., Richard, J., et al. 2023, MNRAS, 520, 2180 [NASA ADS] [CrossRef] [Google Scholar]
  35. Collett, T. E., & Auger, M. W. 2014, MNRAS, 443, 969 [NASA ADS] [CrossRef] [Google Scholar]
  36. Collett, T. E., & Smith, R. J. 2020, MNRAS, 497, 1654 [NASA ADS] [CrossRef] [Google Scholar]
  37. Collett, T. E., Buckley-Geer, E., Lin, H., et al. 2017, ApJ, 843, 148 [NASA ADS] [CrossRef] [Google Scholar]
  38. Conroy, C., Dutton, A. A., Graves, G. J., Mendel, J. T., & van Dokkum, P. G. 2013, ApJ, 776, L26 [NASA ADS] [CrossRef] [Google Scholar]
  39. Cowie, L. L., Hu, E. M., & Songaila, A. 1995, AJ, 110, 1576 [NASA ADS] [CrossRef] [Google Scholar]
  40. Dekel, A., Ishai, G., Dutton, A. A., & Maccio, A. V. 2017, MNRAS, 468, 1005 [NASA ADS] [CrossRef] [Google Scholar]
  41. DES Collaboration (Abbott, T. M. C., et al.) 2022, Phys. Rev. D, 105, 043512 [NASA ADS] [CrossRef] [Google Scholar]
  42. Dessauges-Zavadsky, M., Richard, J., Combes, F., et al. 2019, Nat. Astron., 3, 1115 [NASA ADS] [CrossRef] [Google Scholar]
  43. Dutton, A. A., & Treu, T. 2014, MNRAS, 438, 3594 [NASA ADS] [CrossRef] [Google Scholar]
  44. Efstathiou, G., Moody, S., Peacock, J. A., et al. 2002, MNRAS, 330, L29 [NASA ADS] [CrossRef] [Google Scholar]
  45. Eisenstein, D. J., Zehavi, I., Hogg, D. W., et al. 2005, ApJ, 633, 560 [Google Scholar]
  46. El-Zant, A., Shlosman, I., & Hoffman, Y. 2001, ApJ, 560, 636 [NASA ADS] [CrossRef] [Google Scholar]
  47. El-Zant, A. A., Hoffman, Y., Primack, J., Combes, F., & Shlosman, I. 2004, ApJ, 607, L75 [NASA ADS] [CrossRef] [Google Scholar]
  48. Elmegreen, B. G., & Elmegreen, D. M. 2005, ApJ, 627, 632 [Google Scholar]
  49. Elmegreen, B. G., Elmegreen, D. M., Fernandez, M. X., & Lemonias, J. J. 2009, ApJ, 692, 12 [CrossRef] [Google Scholar]
  50. Elmegreen, B. G., Elmegreen, D. M., Sánchez Almeida, J., et al. 2013, ApJ, 774, 86 [NASA ADS] [CrossRef] [Google Scholar]
  51. Faure, B., Bournaud, F., Fensch, J., et al. 2021, MNRAS, 502, 4641 [CrossRef] [Google Scholar]
  52. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  53. Förster Schreiber, N. M., Genzel, R., Bouché, N., et al. 2009, ApJ, 706, 1364 [Google Scholar]
  54. Freudling, W., Romaniello, M., Bramich, D. M., et al. 2013, A&A, 559, A96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Fruchter, A. S., Hack, W., Dencheva, N., Droettboom, M., & Greenfield, P. 2010, in 2010 Space Telescope Science Institute Calibration Workshop, 382 [Google Scholar]
  56. Gallazzi, A., Bell, E. F., Zibetti, S., Brinchmann, J., & Kelson, D. D. 2014, ApJ, 788, 72 [Google Scholar]
  57. Gao, L., Navarro, J. F., Frenk, C. S., et al. 2012, MNRAS, 425, 2169 [NASA ADS] [CrossRef] [Google Scholar]
  58. Gavazzi, R., Treu, T., Rhodes, J. D., et al. 2007, ApJ, 667, 176 [Google Scholar]
  59. Gerhard, O., Kronawitter, A., Saglia, R. P., & Bender, R. 2001, AJ, 121, 1936 [Google Scholar]
  60. Gnedin, O. Y., Kravtsov, A. V., Klypin, A. A., & Nagai, D. 2004, ApJ, 616, 16 [Google Scholar]
  61. Golse, G., & Kneib, J. P. 2002, A&A, 390, 821 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  62. Governato, F., Brook, C., Mayer, L., et al. 2010, Nature, 463, 203 [NASA ADS] [CrossRef] [Google Scholar]
  63. Graham, A. W., Merritt, D., Moore, B., Diemand, J., & Terzić, B. 2006, AJ, 132, 2701 [NASA ADS] [CrossRef] [Google Scholar]
  64. Grillo, C. 2010, ApJ, 722, 779 [NASA ADS] [CrossRef] [Google Scholar]
  65. Grillo, C., Lombardi, M., & Bertin, G. 2008, A&A, 477, 397 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Grillo, C., Gobat, R., Lombardi, M., & Rosati, P. 2009, A&A, 501, 461 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Grillo, C., Suyu, S. H., Rosati, P., et al. 2015, ApJ, 800, 38 [Google Scholar]
  68. Grillo, C., Rosati, P., Suyu, S. H., et al. 2018, ApJ, 860, 94 [Google Scholar]
  69. Grillo, C., Rosati, P., Suyu, S. H., et al. 2020, ApJ, 898, 87 [Google Scholar]
  70. Guo, Y., Ferguson, H. C., Bell, E. F., et al. 2015, ApJ, 800, 39 [NASA ADS] [CrossRef] [Google Scholar]
  71. Guo, Y., Rafelski, M., Bell, E. F., et al. 2018, ApJ, 853, 108 [NASA ADS] [CrossRef] [Google Scholar]
  72. Gustafsson, M., Fairbairn, M., & Sommer-Larsen, J. 2006, Phys. Rev. D, 74, 123522a [NASA ADS] [CrossRef] [Google Scholar]
  73. Hezaveh, Y. D., Dalal, N., Marrone, D. P., et al. 2016, ApJ, 823, 37 [Google Scholar]
  74. Hinshaw, G., Larson, D., Komatsu, E., et al. 2013, ApJS, 208, 19 [Google Scholar]
  75. Iani, E., Zanella, A., Vernet, J., et al. 2021, MNRAS, 507, 3830 [NASA ADS] [CrossRef] [Google Scholar]
  76. Jester, S., Schneider, D. P., Richards, G. T., et al. 2005, AJ, 130, 873 [Google Scholar]
  77. Johnson, L. E., Irwin, J. A., White, R. E., et al. 2018, ApJ, 856, 131 [NASA ADS] [CrossRef] [Google Scholar]
  78. Jullo, E., Natarajan, P., Kneib, J. P., et al. 2010, Science, 329, 924 [Google Scholar]
  79. Kausch, W., Noll, S., Smette, A., et al. 2015, A&A, 576, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Kneib, J. P., Ellis, R. S., Smail, I., Couch, W. J., & Sharples, R. M. 1996, ApJ, 471, 643 [Google Scholar]
  81. Koopmans, L. V. E., Bolton, A., Treu, T., et al. 2009, ApJ, 703, L51 [NASA ADS] [CrossRef] [Google Scholar]
  82. Laporte, C. F. P., & White, S. D. M. 2015, MNRAS, 451, 1177 [NASA ADS] [CrossRef] [Google Scholar]
  83. Li, Y., Ruszkowski, M., & Bryan, G. L. 2017, ApJ, 847, 106 [Google Scholar]
  84. Livermore, R. C., Jones, T. A., Richard, J., et al. 2015, MNRAS, 450, 1812 [Google Scholar]
  85. Loeb, A., & Peebles, P. J. E. 2003, ApJ, 589, 29 [NASA ADS] [CrossRef] [Google Scholar]
  86. Loewenstein, M., & Mushotzky, R. 2003, Nucl. Phys. B - Proc. Suppl., 124, 91 [NASA ADS] [CrossRef] [Google Scholar]
  87. Mahler, G., Richard, J., Clément, B., et al. 2018, MNRAS, 473, 663 [Google Scholar]
  88. Martizzi, D., Teyssier, R., Moore, B., & Wentz, T. 2012, MNRAS, 422, 3081 [NASA ADS] [CrossRef] [Google Scholar]
  89. Martizzi, D., Teyssier, R., & Moore, B. 2013, MNRAS, 432, 1947 [NASA ADS] [CrossRef] [Google Scholar]
  90. McCully, C., & Tewes, M. 2019, Astrophysics Source Code Library [record ascl:1907.032] [Google Scholar]
  91. Messa, M., Dessauges-Zavadsky, M., Richard, J., et al. 2022, MNRAS, 516, 2420 [Google Scholar]
  92. Meštrić, U., Vanzella, E., Zanella, A., et al. 2022, MNRAS, 516, 3532 [CrossRef] [Google Scholar]
  93. Modigliani, A., Goldoni, P., Royer, F., et al. 2010, in Observatory Operations: Strategies, Processes, and Systems III, eds. D. R. Silva, A. B. Peck, & B. T. Soifer, SPIE Conf. Ser., 7737, 773728 [NASA ADS] [CrossRef] [Google Scholar]
  94. Munari, U., Sordo, R., Castelli, F., & Zwitter, T. 2005, A&A, 442, 1127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  95. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1996, ApJ, 462, 563 [Google Scholar]
  96. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 [Google Scholar]
  97. Navarro, J. F., Ludlow, A., Springel, V., et al. 2010, MNRAS, 402, 21 [Google Scholar]
  98. Newman, A. B., Treu, T., Ellis, R. S., & Sand, D. J. 2013a, ApJ, 765, 25 [NASA ADS] [CrossRef] [Google Scholar]
  99. Newman, A. B., Treu, T., Ellis, R. S., et al. 2013b, ApJ, 765, 24 [NASA ADS] [CrossRef] [Google Scholar]
  100. Newman, A. B., Ellis, R. S., & Treu, T. 2015, ApJ, 814, 26 [NASA ADS] [CrossRef] [Google Scholar]
  101. Nipoti, C., Treu, T., Ciotti, L., & Stiavelli, M. 2004, MNRAS, 355, 1119 [NASA ADS] [CrossRef] [Google Scholar]
  102. Noll, S., Burgarella, D., Giovannoli, E., et al. 2009, A&A, 507, 1793 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  103. Oklopčić, A., Hopkins, P. F., Feldmann, R., et al. 2017, MNRAS, 465, 952 [CrossRef] [Google Scholar]
  104. Pedrosa, S., Tissera, P. B., & Scannapieco, C. 2009, MNRAS, 395, L57 [NASA ADS] [Google Scholar]
  105. Peirani, S., Kay, S., & Silk, J. 2008, A&A, 479, 123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  106. Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2002, AJ, 124, 266 [Google Scholar]
  107. Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2010, AJ, 139, 2097 [Google Scholar]
  108. Perlmutter, S., Aldering, G., Goldhaber, G., et al. 1999, ApJ, 517, 565 [Google Scholar]
  109. Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  110. Pontzen, A., & Governato, F. 2012, MNRAS, 421, 3464 [NASA ADS] [CrossRef] [Google Scholar]
  111. Riess, A. G., Filippenko, A. V., Challis, P., et al. 1998, AJ, 116, 1009 [Google Scholar]
  112. Rigby, J. R., Johnson, T. L., Sharon, K., et al. 2017, ApJ, 843, 79 [Google Scholar]
  113. Ritondale, E., Vegetti, S., Despali, G., et al. 2019, MNRAS, 485, 2179 [Google Scholar]
  114. Rizzo, F., Vegetti, S., Fraternali, F., Stacey, H. R., & Powell, D. 2021, MNRAS, 507, 3952 [NASA ADS] [CrossRef] [Google Scholar]
  115. Romano-Díaz, E., Shlosman, I., Hoffman, Y., & Heller, C. 2008, ApJ, 685, L105 [CrossRef] [Google Scholar]
  116. Rubin, V. C. 1983, Sci. Am., 248, 96 [NASA ADS] [CrossRef] [Google Scholar]
  117. Rubin, V. C., Ford, W., & Kent, J. 1970, ApJ, 159, 379 [NASA ADS] [CrossRef] [Google Scholar]
  118. Salpeter, E. E. 1955, ApJ, 121, 161 [Google Scholar]
  119. Sand, D. J., Treu, T., Smith, G. P., & Ellis, R. S. 2004, ApJ, 604, 88 [NASA ADS] [CrossRef] [Google Scholar]
  120. Schneider, P. 2014, A&A, 568, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  121. Schuldt, S., Chirivì, G., Suyu, S. H., et al. 2019, A&A, 631, A40 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  122. Schwarz, G. 1978, Ann. Stat., 6, 461 [Google Scholar]
  123. Scolnic, D. M., Jones, D. O., Rest, A., et al. 2018, ApJ, 859, 101 [NASA ADS] [CrossRef] [Google Scholar]
  124. Sellwood, J. A., & McGaugh, S. S. 2005, ApJ, 634, 70 [NASA ADS] [CrossRef] [Google Scholar]
  125. Selsing, J., Malesani, D., Goldoni, P., et al. 2019, A&A, 623, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  126. Sérsic, J. L. 1963, Boletin de la Asociacion Argentina de Astronomia La Plata Argentina, 6, 41 [Google Scholar]
  127. Shu, Y., Bolton, A. S., Brownstein, J. R., et al. 2015, ApJ, 803, 71 [NASA ADS] [CrossRef] [Google Scholar]
  128. Smette, A., Sana, H., Noll, S., et al. 2015, A&A, 576, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  129. Smith, R. J., & Collett, T. E. 2021, MNRAS, 505, 2136 [NASA ADS] [CrossRef] [Google Scholar]
  130. Smoot, G. F., Bennett, C. L., Kogut, A., et al. 1992, ApJ, 396, L1 [Google Scholar]
  131. Sommer-Larsen, J., & Limousin, M. 2010, MNRAS, 408, 1998 [NASA ADS] [CrossRef] [Google Scholar]
  132. Sonnenfeld, A., Treu, T., Marshall, P. J., et al. 2015, ApJ, 800, 94 [Google Scholar]
  133. Sonnenfeld, A., Jaelani, A. T., Chan, J., et al. 2019, A&A, 630, A71 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  134. Stark, D. P., Auger, M., Belokurov, V., et al. 2013, MNRAS, 436, 1040 [Google Scholar]
  135. Sunnquist, B. 2018, Instrument Science Report WFC3, 2018-6 [Google Scholar]
  136. Suyu, S. H., & Halkola, A. 2010, A&A, 524, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  137. Suyu, S. H., Marshall, P. J., Hobson, M. P., & Blandford, R. D. 2006, MNRAS, 371, 983 [Google Scholar]
  138. Suyu, S. H., Hensel, S. W., McKean, J. P., et al. 2012, ApJ, 750, 10 [Google Scholar]
  139. Suyu, S. H., Auger, M. W., Hilbert, S., et al. 2013, ApJ, 766, 70 [Google Scholar]
  140. Suyu, S. H., Bonvin, V., Courbin, F., et al. 2017, MNRAS, 468, 2590 [Google Scholar]
  141. Tamburello, V., Rahmati, A., Mayer, L., et al. 2017, MNRAS, 468, 4792 [NASA ADS] [CrossRef] [Google Scholar]
  142. Tanaka, M., Wong, K. C., More, A., et al. 2016, ApJ, 826, L19 [NASA ADS] [CrossRef] [Google Scholar]
  143. Tonini, C., Lapi, A., & Salucci, P. 2006, ApJ, 649, 591 [NASA ADS] [CrossRef] [Google Scholar]
  144. Treu, T. 2010, ARA&A, 48, 87 [NASA ADS] [CrossRef] [Google Scholar]
  145. Treu, T., & Koopmans, L. V. E. 2002, ApJ, 575, 87 [NASA ADS] [CrossRef] [Google Scholar]
  146. Tu, H., Gavazzi, R., Limousin, M., et al. 2009, A&A, 501, 475 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  147. Vanzella, E., Caminha, G. B., Calura, F., et al. 2020, MNRAS, 491, 1093 [Google Scholar]
  148. Vegetti, S., Lagattuta, D. J., McKean, J. P., et al. 2012, Nature, 481, 341 [NASA ADS] [CrossRef] [Google Scholar]
  149. Vernet, J., Dekker, H., D’Odorico, S., et al. 2011, A&A, 536, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  150. Wang, H., Cañameras, R., Caminha, G. B., et al. 2022, A&A, 668, A162 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  151. Zanella, A., Le Floc’h, E., Harrison, C. M., et al. 2019, MNRAS, 489, 2792 [Google Scholar]
  152. Zitrin, A., Rosati, P., Nonino, M., et al. 2012, ApJ, 749, 97 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Spectra of the multiple images

We show here the 1D spectra of the multiple images of the AB source that are covered by the X-shooter data (as described at the beginning of Sect. 3).

thumbnail Fig. A.1.

X-shooter 1D spectra of the multiple images A1 and B1. Zoom in of the wavelength range where we observe the main emission lines used to confirm the redshift measurement. The black lines are the observed spectra in units of 10−17 erg s−1 cm−2 Å−1, and the gray regions indicate the data variance.

thumbnail Fig. A.2.

X-shooter 1D spectra of the multiple images A3 and B3. See the caption of Fig. A.1 for further details.

thumbnail Fig. A.3.

X-shooter 1D spectra of the multiple images A4 and B4. See the caption of Fig. A.1 for further details.

Appendix B: Spectrum of the main deflector

We show here the 1D spectrum of the main lens galaxy.

thumbnail Fig. B.1.

1D X-shooter spectrum of the main lens galaxy. It is smoothed with a boxcar filter, and the most prominent spectral features are marked.

All Tables

Table 1.

Modeling of the deflector total mass distribution using the observed multiple image positions.

In the text
Table 2.

Modeling of the deflector total mass distribution using the extended source modeling.

In the text
Table 3.

Angular separation and sizes of the reconstructed sources.

In the text
Table 4.

Median and 68% confidence level uncertainty values of the parameters of the different models, based on the point-like (top) and extended (bottom) modeling of the sources, as described in Tables 1 and 2.

In the text

All Figures

thumbnail Fig. 1.

F160W image of the strong lensing system analyzed in this work. The multiple images corresponding to each of the three background sources considered in the analysis are labeled A1–A4, B1–B4, and C1–C2. Later in this study we show that sources A and B represent two emission peaks of a single extended source.
In the text
thumbnail Fig. 2.

Spectral energy distribution of the main lens galaxy of SDSS J0100+1818. The observed flux densities from NOT (red triangles), Pan-STARRS (red circles), and HST (red square) are used to fit the SED with the CIGALE software. Error bars are smaller than the symbols. The best-fit stellar continuum template (blue curve) corresponds to the model fluxes plotted in orange and is used to infer the total stellar mass value of the deflector. The bottom panel shows the relative residuals of the fit.
In the text
thumbnail Fig. 3.

Results of the photometric modeling of the system with GALFIT. Left: considered cutout of the HST image in the F160W band. Center: best-fit model, consisting of a combination of two Sérsic profiles for the main lens in the middle, and a single Sérsic profile for each of the neighboring galaxies. Right: normalized residuals on a color scale from approximately −2σ to +2σ. The white pixels, with values of 0.0, show the regions masked during the luminosity profile modeling.
In the text
thumbnail Fig. 4.

Cumulative stellar mass (M*) profile (solid blue line) and ±1σ uncertainties (shaded area) of the main lens. It is measured from the total stellar mass measured through the SED fitting and assuming a constant stellar mass-to-light ratio. The vertical black line shows the half-light radius of the lens galaxy.
In the text
thumbnail Fig. 5.

Difference between the redshift distribution functions of galaxies located within 0.5′ of the main ultra-massive lens galaxy at z = 0.581, and over the rest of the 2.5′ × 2.5′ field, after normalizing to the same area (orange curve). The black histogram is obtained from the peak of the redshift distribution of each source. A total of 14 galaxies within 0.5′ of the main deflector have redshifts consistent, at the 2σ confidence level, with z = 0.581, and 39 over the rest of our analyzed field.
In the text
thumbnail Fig. 6.

HST/F160W image of SDSS J0100+1818 and comparison between the observed multiple images of sources A (blue), B (red), and C (green), represented with circles, and those predicted by the best-fit PIEMD+rc model (fifth entry of Table 1), represented with crosses in the corresponding colors.
In the text
thumbnail Fig. 7.

Probability density distributions of the parameters of the PIEMD+rc (in blue) and SPEMD (in red) models. The marginalized 1D histograms of each parameter are shown along the diagonal, while the other panels show the joint 2D probability distributions of the two parameters reported on the horizontal and vertical axes. The parameters plotted are those introduced in Sect. 4.1. The vertical dashed lines in the 1D histograms represent the 16th, 50th, and 84th percentiles, while the solid lines in the 2D distributions represent the 0.68 and 0.95 contour levels.
In the text
thumbnail Fig. 8.

Probability density distribution of the parameters of the SPEMD model considering the multiple image positions of sources A and B only (in purple) and of sources A, B, and C (in green). The panels show the probability distributions as described in Fig. 7.
In the text
thumbnail Fig. 9.

Cumulative projected total mass profiles for the PIEMD+rc (blue) and SPEMD (red) models with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C as point-like objects. The cumulative projected stellar mass profile from Fig. 4 is shown in dashed black. The vertical lines close to the x-axis locate the distances from the lens galaxy center of the different multiple images, color-coded following Figs. 1 and 6. The black line shows the effective radius of the main lens galaxy.
In the text
thumbnail Fig. 10.

Cumulative projected stellar-over-total mass fraction profiles for the PIEMD+rc (blue) and SPEMD (red) models with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C as point-like objects. The vertical lines are the same as in Fig. 9.
In the text
thumbnail Fig. 11.

Probability density distribution for the SPEMD model with two sources (source AB and source C). The marginalized 1D histograms of each parameter are shown along the diagonal, while the other panels show the joint 2D probability distributions of the two parameters reported on the x and y axes. The parameters are those described in Sect. 4.1. The vertical dashed lines in the 1D histograms represent the 16th, 50th and 84th percentiles, which are also reported on top of each 1D histogram, while the solid lines and shaded areas in the 2D distributions represent the 0.68, 0.95, and 0.99 contour levels. In the floating panel, we zoom in the marginalized 1D histogram of the γ′ parameter (in black), and compare it with its distribution for the model with the AB source only (in red). The introduction of source C provides tighter constraints on the value of the γ′ parameter.
In the text
thumbnail Fig. 12.

PIEMD model with one source. From the left to the right: observed SB in the F160W band of the multiple images considered in the extended source modeling; model-predicted SB; normalized residuals in the range from −5σ to +5σ; reconstructed SB of the source on a 40 × 40 pixel grid. Angular scales of 1″ and of 0.5″ are represented on the lens and source plane, respectively.
In the text
thumbnail Fig. 13.

As in Fig. 12 for the SPEMD model with one source.
In the text
thumbnail Fig. 14.

As in Fig. 12 for the SPEMD model with two sources. The top row shows sources A and B and their SB is reconstructed on a 40 × 40 pixel grid, while the bottom row shows source C and its SB is reconstructed on a 25 × 25 pixel grid.
In the text
thumbnail Fig. 15.

Cumulative projected total mass profiles for the PIEMD+rc model with point-like (solid black) and extended source (red) modeling, with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C. For the extended source modeling, the uncertainties are smaller than the linewidth. The vertical lines close to the x-axis locate the distances from the lens galaxy center of the different multiple images, color-coded following Figs. 1 and 6. The black line shows the effective radius of the main lens galaxy.
In the text
thumbnail Fig. 16.

Cumulative projected stellar-over-total mass fraction profiles for the PIEMD+rc model with point-like (black) and extended source (red) modeling, with ±1σ uncertainties (shaded areas), obtained by modeling the multiple images of A, B, and C. The vertical lines are the same as in Fig. 9.
In the text
thumbnail Fig. A.1.

X-shooter 1D spectra of the multiple images A1 and B1. Zoom in of the wavelength range where we observe the main emission lines used to confirm the redshift measurement. The black lines are the observed spectra in units of 10−17 erg s−1 cm−2 Å−1, and the gray regions indicate the data variance.
In the text
thumbnail Fig. A.2.

X-shooter 1D spectra of the multiple images A3 and B3. See the caption of Fig. A.1 for further details.
In the text
thumbnail Fig. A.3.

X-shooter 1D spectra of the multiple images A4 and B4. See the caption of Fig. A.1 for further details.
In the text
thumbnail Fig. B.1.

1D X-shooter spectrum of the main lens galaxy. It is smoothed with a boxcar filter, and the most prominent spectral features are marked.
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.