A COMPREHENSIVE VIEW OF A STRONGLY LENSED PLANCK-ASSOCIATED SUBMILLIMETER GALAXY

We obtained a 3440 s K_S-band (hereafter K) image on 2011 April 13 (UT) and a 2100 s J-band image on 2011 June 30 (UT) with the Keck II laser guide star adaptive optics system (LGSAO; Wizinowich et al. 2006). An R = 15.8 mag star 48'' SW of HATLAS12−00 served as the tip-tilt reference star. The estimated Strehl ratios at the source position are ∼23% and 5% in the K and J bands, respectively. We used the NIRC2 camera at 004 pixel⁻¹ scale for both filters (40'' field) and dithered with 2''–3'' steps. The atmospheric seeing at 0.5 μm was ∼04 and 05 during the K- and J-band imaging, respectively.³⁹

We used our IDL (Interactive Data Language) programs to reduce the images. After dark subtraction and flat-fielding, sky background and object masks are updated iteratively. For each frame, after subtracting a scaled median sky, the residual background is removed with B-spline models. In the last iteration, we discard the three frames of the poorest image quality and correct the NIRC2 geometric distortion using the solution of P. B. Cameron⁴⁰ before combining the aligned frames. The resolutions of the final K- and J-band images are 016 and 027 in FWHM, respectively. We measure the FWHMs from the most compact source in the field located 10'' SE of HATLAS12−00 (labeled "PSF" in Figure 1(a)); we also use this object as the point-spread function (PSF) in the lens modeling (Section 3.2). The images are flux calibrated against UKIRT Infrared Sky Survey (UKIDSS; Lawrence et al. 2007) and reach depths of K = 25.6 and J = 25.0 AB for a 5σ detection with 01 and 02 radius apertures,⁴¹ respectively.

Zoom In Zoom Out Reset image size

Limited by the small field of NIRC2, a deep wide-field image is required for astrometry calibration. Optical imaging was obtained with the high-throughput auxiliary-port camera (ACAM) mounted at a folded-Cassegrain focus of the 4.2 m William Herschel Telescope (Benn et al. 2008) on 2011 April 26 (UT). We obtained four images of 200 s on a ∼2' field centered on HATLAS12−00, without any filter. The seeing was ∼09. The images were reduced and combined following standard techniques in IRAF.⁴² No accurate photometric calibration is possible because we did not use any broadband filter. But by comparing sources extracted from the ACAM image and the Sloan Digital Sky Survey (SDSS) i-band catalog in the same field, we find that our image reaches an equivalent i-band 5σ depth of 24.6 AB, or 2.3 mag deeper than the SDSS.

We solve the astrometry of the ACAM image using the on-sky positions of SDSS (Aihara et al. 2011) DR8 sources inside the field. We use the astrometry routines in Marc Buie's IDL library⁴³ to correct for offsets, rotation, and distortions, with four terms (constant, X, Y, and $R = \sqrt{X^2+Y^2}$). Sources that appear blended in the SDSS catalog are excluded. With 35 SDSS sources, we measure 1σ dispersions of δR.A. = 013 and δDecl. = 014 between the astrometry-calibrated ACAM image and the SDSS. Finally, we use the same routines to match the NIRC2 images to the ACAM image with 13 well-detected sources inside the 40'' NIRC2 field of view. The corrected NIRC2 images show 1σ dispersions of δR.A. = 004 and δDecl. = 005.

We obtained SMA interferometric imaging of HATLAS12−00 at 880 μm (339.58 GHz) in the compact array configuration with an on-source integration time (t_int) of 1 hr and at 890 μm (336.9 GHz) in the subcompact array configuration with t_int = 2 hr. The compact and subcompact observations took place on 2011 May 2 and 2012 January 14, respectively. During both nights, atmospheric opacity was low (τ_225 GHz ∼ 0.1) and phase stability was good. Both observations used an intermediate frequency coverage of 4–8 GHz and provide a total of 8 GHz bandwidth (considering both sidebands). The quasars 1229+020 and 1058+015 were used for time-variable gain (amplitude and phase) calibration. The blazar 3C 279 served as the primary bandpass calibrator. For the compact data, we used Titan as the absolute flux calibrator. For the subcompact data, we intended to use Callisto as the flux calibrator, but Jupiter might have fallen into one of the side lobes of the SMA primary beam while we observed Callisto. So we decided to use 3C 279 in lieu of Callisto as the flux calibrator. It is possible to use 3C 279 because we have reliable measurements of its flux both before and after the observation of HATLAS12−00.

We used the invert and clean tasks in the Multichannel Image Reconstruction, Image Analysis, and Display (MIRIAD) software (Sault et al. 1995) to invert the uv visibilities and deconvolve the dirty map, respectively. We used natural weighting to obtain the best sensitivity. For the compact data, the cleaned image has a synthesized beam with a FWHM resolution of 207 × 187 at a position angle (P.A.) of −23.6 deg east of north; for the subcompact data, the beam is 557 × 368 at P.A. = 65.7 deg. The primary beam of the SMA is ∼37''. The rms noise levels are 3.0 mJy beam⁻¹ and 3.6 mJy beam⁻¹ for the compact image and the subcompact image, respectively.

HATLAS12−00 is resolved into two components by the SMA (Figure 1). Taking into account the 10% flux calibration uncertainty, the total flux is 70 ± 10 mJy and 93 ± 12 mJy for the compact image and the subcompact image, respectively. The latter agrees well with the Large APEX BOlometer CAmera (LABOCA) bolometer array flux measurement at 870 μm (Section 2.5). The compact array data did not fully capture the source flux because of the sparser array configuration, i.e., ∼25% of the total flux is distributed on spatial scales larger than those accessible to the SMA in its compact array configuration. So we use the total flux from subcompact data for SED modeling (Section 4.2). We chose to use the compact image for lens modeling (Section 3.3) because of its higher spatial resolution. We find that using the subcompact image or the subcompact+compact combined image does not change the lens modeling result, but they give larger errors for the derived parameters.

We exploited the recent upgrade to the National Radio Astronomy Observatory's⁴⁴ Very Large Array (VLA; Perley et al. 2011), which includes the provision of Ka-band receivers (26.5–40 GHz), to observe the redshifted CO(1→0) emission from HATLAS12−00 at 27.06532 GHz (ν_rest = 115.27120 GHz; Morton & Noreau 1994).

Observations were carried out dynamically during excellent weather conditions on 2012 January 6 and 8. During this Open Shared Risk Observing period the available bandwidth from the new Wideband Interferometric Digital ARchitecture correlator consisted of two independently tunable output pairs of eight sub-bands each, with 64 × 2 MHz full-polarization channels per sub-band, giving a total bandwidth of 2048 MHz. At the redshift of HATLAS12−00, however, the CO(1→0) line could be reached by only the BD output pair, giving ∼11, 350 km s⁻¹ coverage and ∼22 km s⁻¹ resolution. We offset our tuning by 64 MHz to avoid noisier edge channels. The eight sub-bands of output pair AC were tuned to 32.5 GHz.

The bright compact calibration source, J1150−0023, was observed every few minutes to determine accurate complex gain solutions and bandpass corrections. 3C 286 (S = 2.1666 Jy at 27.06 GHz) was also observed to set the absolute flux scale, and the pointing accuracy was checked locally every hour. In total, around 2 hr of data were obtained for HATLAS12−00, with ∼1 hr of calibration.

The data were reduced using $\cal AIPS$ (31DEC12) following the procedures described by Ivison et al. (2011), though with a number of important changes: data were loaded using bdf2aips and fring was used to optimize the delays, based on 1 minute of data for 3C 286. The base bands were knitted together using the noifs task, yielding uv data sets with 512 × 2 MHz channels, which we then added together using the task dbcon. Finally, the channels were imaged over a 512 × 512 × 03 field, with natural weighting (robust = 5), to form a 512³ cube centered on HATLAS12−00. Integrating over those 55 channels found to contain line emission (so a fwzi of ∼1200 km s⁻¹) yielded an rms noise level of 27 μJy beam⁻¹.

The cleaned and velocity-integrated CO map is shown in Figure 1(c). The beam is 25×22 at P.A. = 85°. Similar to the SMA, the map resolves two components separated by ∼5''. The CO lines extracted from the two components show the same redshift and line profile, further confirming that they are lensed images of a single source. The best-fit Gaussian to the area-integrated spectrum gives a line width of ΔV_FWHM = 585 ± 55 km s⁻¹ and a line flux of S_COΔV = 1.52 ± 0.20 Jy km s⁻¹. In comparison, the CO(1→0) measurements reported by Harris et al. (2012) using Zpectrometer on the Green Bank Telescope are ΔV_FWHM = 680  ±  80 km s⁻¹ and S_COΔV = 1.18  ±  0.26 Jy km s⁻¹(corrected for the 20% difference in the absolute flux density of 3C 286). The reason for the discrepancy is unclear, but the two line flux measurements agree within the 1σ errors. So, hereafter, we use the weighted mean of the two measurements, S_COΔV = 1.40 ± 0.22 Jy km s⁻¹, to derive the molecular gas mass.

Photometry of HATLAS12−00 was obtained from the SDSS (u, g, r, i, and z), the UKIDSS (Y, J, H, and K), the Wide-Field Infrared Survey Explorer (WISE, 3.6 and 4.6 μm; Wright et al. 2010), the Herschel/PACS (100 and 160 μm; Program ID: OT1_RIVISON_1; Ibar et al. 2010), the Herschel/SPIRE (250, 350, and 500 μm; Pascale et al. 2011; Rigby et al. 2011), the LABOCA (870 μm; Siringo et al. 2009), the SMA (880 μm), the Max-Planck Millimetre Bolometer (MAMBO, 1.2 mm; Kreysa et al. 1999), the CARMA (2792 and 3722 μm; Bock et al. 2006), and the VLA FIRST survey (21 cm; Becker et al. 1995).

We obtained imaging at 870 μm with the LABOCA bolometer array at the Atacama Pathfinder EXperiment (APEX) telescope in 2011 November (D. L. Clements et al., in preparation). LABOCA observed a 114 diameter field with a resolution of FWHM = 186. The observations have a total integration time of ∼30 hr reaching a 1σ sensitivity of ∼2 mJy.

We obtained 1.2 mm imaging with MAMBO at the IRAM 30 m telescope (FWHM ∼ 107) in 2011 January and February (H. Dannerbauer et al., in preparation). Observing time in the on–off mode is 24 minutes, achieving a 1σ sensitivity of ∼1 mJy.

We obtained continuum observations at 81.2 and 108.2 GHz (3722 and 2792 μm; covering rest-frame CO(3→2) and CO(4→3) lines) on 2011 March 18 and September 1 as part of our CO follow-up campaign of bright, lensed H-ATLAS SMGs with CARMA in D array (D. A. Riechers et al., in preparation). Observations were carried out for 0.9 and 1.4 hr on source, respectively, using the 3 mm receivers and a bandwidth of 3.7 GHz per sideband. HATLAS12−00 is unresolved in these observations, with angular resolutions of 68 × 50 and 60 × 38 at 81.2 and 108.2 GHz, respectively (restored with natural baseline weighting).

Table 1 lists the photometry. We have included in the errors the absolute flux calibration uncertainties (3% for WISE, 3%–5% for PACS, 7% for SPIRE, 10% for SMA, and 15% for LABOCA, MAMBO, and CARMA).

We use lenstool (Kneib et al. 1996; Jullo et al. 2007) to find the best-fit parameters and their errors from the peak positions. lenstool implements a Bayesian Markov Chain Monte Carlo sampler to derive the posterior distribution of each parameter and an estimate of the evidence for the model.

The lensing system is mainly made of two red filaments that are ∼35 apart (Figure 1(a)). The outward curved shape of the northern arc can be explained if the source is intrinsically curved. Hence, we split each of the two arcs into three parts and build a simple lens model by putting two deflectors centered on G1 and G2. We find that the predicted counterimages can explain the additional features close to G1 and G2. Guided by the predicted counterimages, we define three systems of lensed images (Figure 1(a), inset). The 11 peak positions in three separate systems provide a total of 16 constraints (11 × 2 − 3 × 2), allowing us to include shear from nearby galaxies G3 and G4.

For the lensing galaxies, we find photometric redshifts of z_{G1 + G2} = 1.06 ± 0.16 and z_G4 = 0.80 ± 0.28 with the public photo-z code EAZY (Brammer et al. 2008). We obtain the nine-band photometry from the SDSS (u, g, r, i, and z) and the UKIDSS (Y, J, H, and K) surveys. At these wavelengths, the flux from the lensed galaxy is negligible (Section 4). The redshift of G1+G2 is measured from the total fluxes of G1 and G2, because they are blended in the seeing-limited data. G3 is undetected in SDSS but shows similar color as G1 and G2. Hence in the lens modeling we assume that all four galaxies are at z = 1.06. Note that although redshift errors of the lensing galaxies would lead to errors in the estimated lens masses, they would not change our conclusions on the lensed galaxy because the magnification factors would remain the same.

For the lens model, we assume that the dark matter plus baryonic mass profiles of the foreground lens galaxies G1 to G4 can be described as singular isothermal ellipsoids (SIEs; Kormann et al. 1994). The SIE profile is parameterized by the velocity dispersion (σ), the position (x, y), the axis ratio (q = b/a), and the P.A. (θ, E of N). We fix the positions to the centers of the galaxies. For G3 and G4, we further fix their q and θ to those from the light distribution, because they are not well constrained by the peak positions and there are significant correlations between the P.A. and ellipticity of the light and of the mass distribution (e.g., Sluse et al. 2012). Therefore, we have a total of eight free parameters. We find a best fit with χ² = 7.9 for dof = 8 (degrees of freedom) and an average positional error of 004 (∼1 pixel). The parameters and their errors are summarized in Table 2. We also list the mass enclosed by the critical curve for each SIE,

$$\begin{equation} M_E = \frac{4\pi ^2}{G} \frac{D_{\rm L} D_{{\rm LS}}}{D_{\rm S}} \frac{\sigma ^4}{c^2}, \end{equation} \tag{ 1 }$$

where D_L, D_S, and D_LS are the angular diameter distances to the lens, to the source, and between the lens and the source, respectively. The radius of the area enclosed by the critical curve can be approximated by the circularized Einstein radius:

$$\begin{equation} b = 4.5 \left(\frac{\sigma }{200 {\rm \,km\,s}^{-1}} \right)^2 \sqrt{\frac{2q}{1+q^2}} {\rm \,kpc}. \end{equation} \tag{ 2 }$$

In the errors of masses and velocity dispersions, we have included the 1σ uncertainty of the photometric redshift.

The nominal model described above is the most favorable description of the lensing system because of the following.

• 1.  
  Adding q's and θ's of G3 and G4 as free parameters does not substantially improve the fit: the Bayesian evidence⁴⁵ increases by only Δln (E) = 0.5 and the reduced χ² actually increases from 1.0 to 1.6 as a result of the decreased degree of freedom.
• 2.  
  Excluding the potentials of G3 and/or G4 does degrade the fit significantly. The Bayesian evidence decreases by Δln (E) = 2.5 and 30, and the reduced χ² increases from 1.0 to 1.8 and 6.5, when we exclude G4 and both G3 and G4, respectively.
• 3.  
  Including a group-scale potential with a Pseudo-Isothermal Elliptic Mass Distribution (PIEMD) profile (Kassiola & Kovner 1993) does not improve the fit. For the PIEMD profile, we adopt a cutoff radius of 500 kpc but allow the position, ellipticity, P.A., core radius, and velocity dispersion to vary. To limit the number of free parameters, we fix the q's and P.A.'s of the SIEs to those measured from the lensing galaxies but allow their velocity dispersions to vary. So we have a total of 10 free parameters. The reduced χ² of the best fit is 1.9, much higher than that of the nominal model. The Bayesian evidence also decreases by Δln (E) = 6.9 when compared with the nominal model.

Although the peak positions can constrain the deflectors through ray tracing, they cannot provide an accurate estimate of the magnification factor because the source plane light distribution is not taken into account. Because we want to estimate the intrinsic properties of the lensed galaxy, we are interested in the luminosity-weighted magnification factor, which depends on the source morphology because the magnification factor is different at each source plane position. In this section, we model the morphology of the lensed galaxy and refine the lensing potentials simultaneously with the K-band image. The PSF is derived from the most compact source in the field, which is 10'' SE of HATLAS12−00 (Figure 1(a)).

Following Section 3.1, we assume that the source consists of three clumps, each described as a Sérsic profile. Again we use SIE profiles for the lensing potentials. So we have a total of 29 parameters: 7 parameters for each Sérsic profile and 8 parameters for the SIE potentials. Our fitting procedure is as follows. For an initial set of parameters describing the source and the lenses from Section 3.1, we use lenstool to generate a lensed image of the source, which is then convolved with the PSF and compared with the observed image. We limit the comparison in a 33 × 45 (83 × 113 pixels) rectangular region that encloses the lensing features. This process is iterated with AMOEBA_SA to find the parameters that minimize the residual between the observation and the model. AMOEBA_SA is based on the IDL multidimensional minimization routine AMOEBA (Press et al. 1992) with simulated annealing added by E. Rosolowsky. We allow a maximum of 1000 iterations in each call of AMOEBA_SA. For the simulated annealing, we adopt an initial "temperature" of 100 and decrease it by 20% in each subsequent call to AMOEBA_SA. A good fit with a reduced χ² around unity is normally found after a few calls to AMOEBA_SA. For each iteration, we compute the total luminosity-weighted magnification factor (μ_K) by summing the pixel values in the image and the source planes with apertures matched by inverting the image plane aperture to the source plane. The 1σ confidence interval of μ_K is found with χ²(μ) − χ²_min ⩽ 1. Note that we compute the χ² values on the residual image binned by 4 pixel boxes (FWHM = 016 = 4 pixel), so that the noise becomes uncorrelated between pixels; or equivalently, one could divide the χ² values from the original residual images by a factor of 16. We find the luminosity-weighted magnification to be μ_K = 16.7 ± 0.8. The best-fit parameters for the deflectors are listed in Table 2. The results are very similar to those from fitting the peak positions, although the errors are smaller because the entire image provides more information than the peak positions alone.

Figure 2 shows the best-fit model. The lensed galaxy has a curved morphology, causing the northern arc bending in the opposite direction of the deflectors. In the source plane, the three clumps extend over only 021 or 1.6 kpc, and their effective radii are 021 ± 004 (1.5 ± 0.3 kpc), 0085 ± 0013 (0.6 ± 0.1 kpc), and 011 ± 005 (0.8 ± 0.4 kpc) from W to E.

Zoom In Zoom Out Reset image size

The nature of the feature NNE of G2 in the residual image is unclear, but it is unlikely to be at the same redshift as the lensed galaxy: tracing its position to the source plane and imaging it back predicts an unobserved equally bright counterimage 08 S of the southern arc. This feature could therefore be part of the galaxy G2.

Precise astrometry calibration is crucial for a joint analysis of images from different wavelengths. Because the only K source detected by the SMA is HATLAS12−00, we have to estimate the astrometry offset between the two images in a statistical way. Because the Keck image is tied to the SDSS astrometry and the SMA image is tied to the radio reference frame, we cross-correlate the VLA FIRST catalog (Becker et al. 1995) and the SDSS catalog within 1° of HATLAS12−00 and compute the optical–radio separation. Ninety-four radio sources have optical counterparts within 3''. We then fit an elliptical Gaussian to the two-dimensional distribution in ΔR.A. = −(α_FIRST − α_SDSS) and ΔDecl. = δ_FIRST − δ_SDSS. The systematic offset from the peak position of the Gaussian is consistent with zero (ΔR.A. = −007, ΔDecl. = +008). The best-fit Gaussian has σ's of 040 and 030, and a P.A. of 119° for the major axis. Therefore, the 1σ ellipse of the astrometry offset has major/minor semiaxes of 061/045. Our result is consistent with that of Ivezić et al. (2002), who found a ∼01 systematic offset and a 1σ error circle of 047 in radius between FIRST and SDSS astrometry.

We can constrain the astrometry offset further through lens modeling. As demonstrated by Kochanek & Narayan (1992) and Wucknitz (2004), interferometric data are most naturally modeled with the uv-plane visibilities, because this avoids beam deconvolution and naturally handles correlated noise. Here, however, we opt to model the cleaned map directly, because (1) the images are essentially unresolved in the SMA map and (2) we already have a good lens model from the K-band image (Section 3.2). Because of the limited spatial resolution of the SMA 880 μm image, the two centroid positions do not offer enough information to constrain the lens model. Hence, for the deflectors, we fix all the parameters to the best-fit values from Section 3.2; for the source, we assume a circular Gaussian profile with variable position and size. We shift the SMA image relative to the K model on a 2''×2'' grid with 01 steps. At each offset position, we find the best-fit model using the same fitting procedure as in Section 3.2. The modeling is performed on a 51 × 68 pixel (51 × 68) region enclosing the SMA sources. Figure 3(a) shows a map of the minimum χ² values at each offset position. The global best fit, with reduced χ² of unity, is reached when we shift the SMA image 06 E of the K image. The middle panels of Figure 2 show this global best-fit model.

Zoom In Zoom Out Reset image size

The noise of the SMA map is Gaussian but is highly correlated. We compute the rms noise of the SMA map after binning it by boxes of n² pixels. We find that the noise starts to decrease as 1/n for n ≳ 20 pixels (FWHM ≃ 2'' = 20 pixel), indicating that the noise becomes uncorrelated on 20 pixel scales. Therefore, we divide the χ² values from the residual images by a factor of 400, which is equivalent to computing χ² from residual images binned by 20 pixel boxes.

In combination with the 1σ error ellipse from FIRST–SDSS cross-correlation, we determine that the astrometry offset between 880 μm and K images is ΔR.A. = −05 ± 01 and ΔDecl. = 00 ± 02, i.e., the overlapping region between the ellipse and the 1σ contour of the χ² map. Collecting all of the solutions in this permitted offset region satisfying χ²(μ) − χ²_min ⩽ 1, we estimate a luminosity-weighted 880 μm magnification of μ₈₈₀ = 7.6 ± 1.5 and an 880 μm source size of FWHM = 015^+0.14 _{− 0.06} = 1.2^+1.0 _{− 0.5} kpc. Because we have fixed the deflectors with the best-fit parameters from the K band, the errors here do not include the uncertainties of the deflectors. Higher resolution far-IR images are required to constrain the deflectors and the source simultaneously.

Dust-emitting regions are often spatially offset from the UV/optical emitting regions in SMGs (Tacconi et al. 2008; Bothwell et al. 2010; Carilli et al. 2010; Ivison et al. 2010b; Riechers et al. 2010). This is clearly the case for HATLAS12−00, which shows distinctly different morphologies at the K band and 880 μm, even after convolving the K-band image with the SMA beam (compare Figures 2(c) and (f)). From the lens model, we estimate a source plane separation between the 880 μm source and the central K clump of 041 ± 007 or 3.1 ± 0.5 kpc (Figure 2(e)).

If we assume zero astrometry offset between SMA and Keck, then we obtain a model that poorly fits the observation (Δχ² ∼ 4; Figure 3(a)). The lens model gives a slightly larger magnification (μ₈₈₀ = 8.4 ± 1.6) and doubles the source size (FWHM = 2.5^+1.9 _{− 0.3} kpc). However, the source plane separation between the 880 μm source and the central K clump remains the same (3.2 ± 0.2 kpc).

We use the same technique to model the JVLA CO(1→0) map as in Section 3.3. The lensed images are better resolved in the JVLA image than in the SMA image, so we use an elliptical Gaussian instead of a circular Gaussian for the source profile. The model has a total of six free parameters (x, y, FWHM, q, P.A., and flux density). Again, we can constrain the astrometric offset between JVLA and Keck through lens modeling. Figure 3(b) shows the minimum χ² values at each offset position relative to the K image. To deal with the correlated noise, we scale the χ² values from the residual map by the product of the FWHMs of the major and minor axes of the beam. The global best fit, with reduced χ² of unity, is reached when we shift the JVLA image 02 E of the K image. The bottom panels of Figure 2 show the best-fit model. We estimate a CO magnification of μ_CO = 6.9 ± 1.6. Similar to the SMA image, we find a source plane separation of 4.7 ± 1.6 kpc between the cold molecular gas and the stellar emission (i.e., the central K clump). The CO(1→0) is emitted from a more extended region than the dust, but the two spatially overlap (Figure 2(e)). The CO source has FWHM = 09 ± 03 = 6.8 ± 2.3 kpc along the major axis, with an axis ratio of 0.8^+0.2 _{− 0.6}.

The molecular gas disk is massive. The velocity–area-integrated CO brightness temperature of L'_CO = (6.4 ± 1.0) × 10¹¹ K km s⁻¹ pc² indicates a molecular gas reservoir of M_gas = (7.4 ± 2.1) × 10¹⁰ M_☉ after lensing correction, assuming a conversion factor of α_CO = M_gas/L'_CO = 0.8 M_☉ (K km s⁻¹ pc²)⁻¹, which is commonly assumed for starburst environments (Solomon & Vanden Bout 2005). Note, however, that α_CO is uncertain by at least a factor of a few and it may depend on the metallicity, the gas temperature, and the velocity dispersion of the galaxy (Narayanan et al. 2012).

The CO magnification factor determined from the lens model (6.9 ± 1.6) is in excellent agreement with that estimated from the CO luminosity−FWHM correlation. The observed L'_CO and line width indicate a magnification factor of 7 ± 2, based on its deviation from the correlation established by unlensed SMGs (Harris et al. 2012; M. S. Bothwell et al., in preparation). This agreement demonstrates that strongly lensed SMGs may be effectively selected with CO spectroscopy in the future.

Adopting the photometric redshift of 1.06, we model the nine-band photometry (u, g, r, i, z, Y, J, H, and K) of G1+G2 with the stellar population synthesis models of Bruzual & Charlot (2003, hereafter BC03). We assume a Chabrier (2003) initial mass function (IMF), Calzetti et al. (1994) extinction law, and exponentially declining star formation history, with a range of e-folding times (τ = 0.1–30 Gyr) and ages (0.01–12.5 Gyr). For each template, we fit for the stellar mass (M_stellar) and extinction (E(B − V)). The best-fit model gives χ² = 6.4 for dof = 7 (Figure 4(a)). The derived properties of G1+G2 are listed in Table 3. The intrinsic extinction is small (E(B − V) = 0.04^+0.11 _{− 0.04}) and there is very little current star formation (SFR = 0.1^+0.4 _{− 0.1} M_☉ yr⁻¹). The dust-absorbed UV/optical luminosity ((3⁺¹⁴ _{− 3}) × 10¹⁰ L_☉) is less than 0.15% of the total integrated IR luminosity before lensing correction (L_8–1000 = 1.2 × 10¹⁴ L_☉). Therefore, G1 and G2 do not contribute significantly to the far-IR fluxes, in agreement with their absence in the SMA image. The stellar mass from SED modeling is ∼80% of the total mass within the critical curves from lens modeling (Table 3), implying that the galaxies are dominated by stellar mass within ∼7 kpc.

Zoom In Zoom Out Reset image size

In the J and K bands, we obtain the photometry of the SMG with an aperture contoured around the multiply imaged features after subtracting the foreground lenses. We measure K = 21.2  ±  0.1 and (J − K) = 2.1  ±  0.2 in AB magnitudes, consistent with the red J − K colors of unlensed SMGs (Frayer et al. 2004; Dannerbauer et al. 2004). Careful modeling is required to extract SMG photometry from the WISE data because the SMG is blended with the foreground galaxies G1 through G4 (FWHM = 6''–12''). We model the WISE 3.4 μm source with four elliptical Gaussians of the same shape. With galfit (Peng et al. 2010), we fix their positions to those determined from the Keck image, but we allow the Gaussian shape to vary. Then we measure the flux density of G1+G2 (∼205 μJy) decomposed from that of G3 and G4 (∼117 μJy). Finally, the flux density of the SMG (∼37 μJy) is estimated from the excess of G1+G2 relative to their best-fit stellar population synthesis model (Section 4.1). Unfortunately, we cannot separate G3+G4 from G1+G2 in the longer wavelength channels of WISE because of the inferior image quality, so we treat the excesses over the best-fit model of G1+G2 as upper limits for the SMG. The far-IR-to-millimeter SED is dominated by the lensed SMG; therefore no foreground subtraction is necessary. Table 1 summarizes the photometry for HATLAS12−00.

We opt to model the rest-frame optical and far-IR emission separately, instead of fitting them together in a self-consistent way with magphys (da Cunha et al. 2008), because our lens model shows that they are emitted from physically distinctive regions, i.e., the dust that attenuates the optical emission has little to do with the starburst-heated dust that emits in the far-IR.

It is difficult to constrain the stellar population with only three photometric detections (dof = 1) in rest-frame optical. However, we can limit the parameter space by excluding unphysical models, such as those that require negative extinction corrections and those whose ages exceed the cosmic age at z = 3.2592 (1.9 Gyr). We use the same BC03 templates as in the previous section. The shaded region in Figure 4(a) shows all of the permitted models with χ² < dof + 1. These models give a range of extinctions, stellar masses, and SFRs: E(B − V) < 0.94, M_stellar = (3.5 ± 2.4) × 10¹⁰ M_☉, and SFR <2000 M_☉ yr⁻¹ (Table 3). The dust-absorbed UV/optical luminosity ranges from 0% to 140% of the observed L_8–1000, but 90% of the models have dust-absorbed luminosity less than 50% of the observed L_8–1000.

We fit the far-IR SED with a single-temperature-modified blackbody,

$$\begin{equation} S(\nu _{\rm obs}) = \sigma (1-e^{-\tau }) B(\nu _{\rm rest},T) (1+z) \mu /d_L^2, \end{equation} \tag{ 3 }$$

where σ is the total absorption cross-section of dust particles at the optically thick limit (i.e., the size of the dust-obscured region), B(ν, T) is the Planck function, τ = (ν_rest/ν₀)^β = (λ₀/λ_rest)^β the optical depth, μ is the lensing magnification factor, and d_L is the luminosity distance. In the optically thin limit (λ ≫ λ₀), dust mass can be derived based on the knowledge of the opacity κ_d (absorption cross-section per unit mass):

$$\begin{equation} M_{\rm dust} = \frac{S(\nu _{\rm obs}) d_L^2}{\kappa _d(\nu _{\rm rest}) B(\nu _{\rm rest},T) (1+z) \mu }. \end{equation} \tag{ 4 }$$

It is generally assumed that the opacity follows a power law, κ_d(ν)∝ν^β, and has a normalization of κ_d = 0.07 ± 0.02 m² kg⁻¹ at 850 μm (Dunne et al. 2000; James et al. 2002). Both the general "optically thick" (S_ν∝(1 − e^−τ)B_ν(T)) and the optically thin (S_ν∝ν^βB_ν(T)) models provide good fits to the observed SED (Figure 4(b)). For the general model, we use all of the nine detections between 100 μm and 3 mm. The best-fit general model gives χ² = 1.6 for dof = 5, suggesting that the photometric errors have been overestimated. For the optically thin model, we exclude the PACS 100 μm point, which is clearly on the Wien tail where small grains tend to dominate the emission. The best-fit optically thin model gives χ² = 4.0 for dof = 5. The derived parameters are listed in Table 3.

The optically thick model yields dust properties similar to those of the local ultraluminous infrared galaxy (ULIRG) Arp 220 (Rangwala et al. 2011), with the optical depth exceeding unity below rest frame ∼250 μm. The intrinsic 8–1000 μm luminosity of L_8–1000 = (1.7 ± 0.3) × 10¹³ L_☉ classifies HATLAS12−00 as a hyper-luminous infrared galaxy (Hy-LIRG). The IR luminosity implies an SFR of 1900 ± 400 M_☉ yr⁻¹ for a Chabrier (2003) IMF (Kennicutt 1998). Using the values of L_8–1000 and T = 63 ± 2 K in the Stefan–Boltzmann law, we obtain a spherical source radius of 780 ± 100 pc, which is three times larger than that of Arp 220 (230 pc) because of the 10 times greater luminosity. The source radius is comparable to that we derive from the optically thick model (r = 500 ± 60 pc) and is consistent with the size we measure from modeling the SMA image (FWHM = 1.2^+1.0 _{− 0.5} kpc; Section 3.3). Therefore, the optically thick model is preferred.

The radio luminosity from the observed 1.4 GHz flux density is L_1.4 GHz = 4πd²_LS_1.4 GHz(1 + z)^{α − 1} = (7.3  ±  3.4) × 10²⁵ W Hz⁻¹ for a radio spectral index of α = 0.7. Assuming the radio emission is magnified by the same factor as the submillimeter emission, the IR-to-radio luminosity ratio of HATLAS12−00, q_L = log(L_IR/(4.52 THz L_1.4 GHz)) = 2.1 ± 0.2, is consistent with the radio–far-IR correlation of high-redshift starburst galaxies: e.g., Kovács et al. (2006) measured q_L = 2.14 ± 0.12 for 15 SMGs, while Ivison et al. (2010a) measured q_L = 2.40 ± 0.24 for 65 Herschel 250 μm selected galaxies. This suggests that the AGN contribution is insignificant in HATLAS12−00.

We also do not see significant AGN contribution in the mid-IR. In Figure 4(b), we fit the SEDs of the local ULIRGs Mrk 231, Arp 220, and the z = 2.3 SMG "Cosmic Eyelash" (Ivison et al. 2010c) to the far-IR SED. The WISE upper limits lie well below the AGN-dominated ULIRG Mrk 231 but are more consistent with Arp 220 and the Eyelash. Therefore, we conclude that HATLAS12−00 is predominantly a starburst system.