A hyper luminous starburst at z=4.72 magnified by a lensing galaxy pair at z=1.48

[Abridged] We discovered in the Herschel Reference Survey an extremely bright IR source with $S_{500}$~120mJy (Red Virgo 4 - RV4). Based on IRAM/EMIR and IRAM/NOEMA detections of the CO(5-4), CO(4-3), and [CI] lines, RV4 is located at z=4.724, yielding a total observed L$_{IR}$ of 1.1+/-0.6x0$^{14}$L$_{\odot}$. At the position of the Herschel emission, three blobs are detected with the VLA at 10cm. The CO(5-4) line detection of each blob confirms that they are at the same redshift with the same line width, indicating that they are multiple images of the same source. In Spitzer and deep optical observations, two sources, High-z Lens 1 (HL1) West and HL1 East, are detected at the center of the three VLA/NOEMA blobs. These two sources are placed at z=1.48 with XSHOOTER spectra, suggesting that they could be merging and gravitationally lensing the emission of RV4. HL1 is the second most distant lens known to date in strong lensing systems. The Einstein radius of the lensing system is 2.2"+/-0.2 (20kpc). The high redshift of HL1 and the large Einstein radius are highly unusual for a strong lensing system. We present the ISM properties of the background source RV4. Different estimates of the gas depletion time yield low values suggesting that RV4 is a SB galaxy. Among all high-z SMGs, this source exhibits one of the lowest L$_{[CI]}$ to L$_{IR}$ ratios, 3.2+/-0.9x10$^{-6}$, suggesting an extremely short gas tdepl of only 14+/-5Myr. It also shows a relatively high L$_{[CI]}$ to L$_{CO(4-3)}$ ratio (0.7+/-0.2) and low L$_{CO(5-4)}$ to L$_{IR}$ ratio (only ~50% of the value expected for normal galaxies) hinting a low density of gas. Finally, we discuss that the short tdepl of RV4 can be explained by either a very high SFE, which is difficult to reconcile with major mergers simulations of high-z galaxies, or a rapid decrease of SF, which would bias the estimate of tdepl toward low value.


Introduction
In the local Universe, a third of the total bolometric luminosity of galaxies is emitted in the infrared (IR) and submillimetric (submm) domains by dust grains, which reprocess the energy absorbed from the stars and active galactic nuclei (AGN) (Soifer & Neugebauer 1991). In submm galaxies (SMGs, Smail et al. 1997; Barger et al. 1998), highly dust-obscured systems, the energy balance between optical and IR/submm domains is skewed even more towards long-wavelength emission. The dust is heated by numerous young stars causing the SMGs to be extremely bright in the submm regime and reach huge IR luminosities (L IR ) higher than ∼ 10 12 L . More and more of these extreme star-forming galaxies are being found (e.g., Daddi et al. 2009;Negrello et al. 2010;Frayer et al. 2011;Yun et al. 2012;Riechers et al. 2013;Weiß et al. 2013;Ivison et al. 2013;Vieira et al. 2013;Dowell et al. 2014;Cañameras et al. 2015;Díaz-Santos et al. 2016;Strandet et al. 2017;Marrone et al. 2018;Zavala et al. 2018), among which are sources at z > 4, indicating an extremely rapid assembly of these structures. Their high SFRs suggest the presence of large gas reservoirs (e.g., a gas fraction of ∼50%; Fu et al. 2013;Béthermin et al. 2015;Aravena et al. 2016) and their large dust content indicates that their interstellar medium (ISM) is metal enriched. One of the challenges posed by these sources is to understand how this massive and mature ISM can be in place so early in these galaxies. Indeed, these massive galaxies are still not well reproduced by cosmological simulations (Davé et al. 2010;Cousin et al. 2015;Sparre et al. 2015). The study of these starbursting systems, which usually lie above the main-sequence of star-forming galaxies (e.g., Noeske et al. 2007;Elbaz et al. 2007Elbaz et al. , 2011 and whose mid-IR properties resemble those of the very brightest, nearby IR-luminous galaxies (Díaz-Santos et al. 2011), is of paramount importance for providing constraints on models of the formation and evolution of massive galaxies. In fact, the massive ellipticals that we observe in the local Universe could have formed rapidly ∼ 10 Gyr ago and be the remnant of this population of starbursting galaxies at high redshift (McCarthy et al. 2004;Daddi et al. 2007b,a;Tacconi et al. 2008;Cimatti et al. 2008;Haan et al. 2013). Within this scenario, because of their tremendous SFRs, the most vigorously star-forming SMGs rapidly quench their star formation by exhausting their gas reservoir in only a few hundred Myr, making them very rare. For instance, Fu et al. (2013) estimate the space density of SMGs with SFR > 2000 M yr −1 to be 1.4× 10 −5 Mpc −3 .
Nevertheless, samples of SMGs have been built with facilities like the South Pole Telescope (SPT) and the Herschel Space Observatory (e.g., Eales et al. 2010;Negrello et al. 2010;Bussmann et al. 2013;Weiß et al. 2013;Vieira et al. 2013), and statistical studies of their properties can be now found in the literature. Studies of high-z bright SMGs (S 850 >50 mJy) show that the majority of these sources are gravitationally lensed by foreground, massive galaxies, thus amplifying the background source emission. This magnification allows us to detect highredshift galaxies that could have been otherwise missed or observed with a lower signal-to-noise ratio. The ISM of high-z sources can thus be studied through the easier detection of lines, thanks to flux boosting (e.g., Weiß et al. 2013;Alaghband-Zadeh et al. 2013;Béthermin et al. 2015;Bothwell et al. 2017;Yang et al. 2017;Cunningham et al. 2019). Another aspect of gravitational lensing is that the deflection of the light emitted by a background source allows us to probe the mass distribution of the foreground source acting as the lens, constraining dark matter (DM) sub-halo structures (Hezaveh et al. 2016) and/or the initial mass function (IMF; Cañameras et al. 2017), for instance. Finally, gravitational lensing boosts the angular resolution with which we can observe background sources allowing spatially resolved studies of the ISM of high redshift galaxies (e.g., Swinbank et al. 2015;Cañameras et al. 2018;Litke et al. 2019;Yang et al. 2019a;Apostolovski et al. 2019).
In this study, we report the discovery of a peculiar strong lensing system. The high-z galaxy, HRS188.6868+7.1357, hereafter Red Virgo 4 (RV4), was first detected in a survey of local galaxies carried out with Herschel/SPIRE (Griffin et al. 2010;Boselli et al. 2010). It is detected in the background of the SPIRE pointing of IC 3521 at RA: 188.6868 and Dec: +7.1357 (J2000). The high flux density (122 mJy at 500 µm) as well as the redward increasing spectrum led to interest in this source. We gathered and obtained ancillary data to constrain its infrared (IR) emission, measure its redshift, identify the foreground source acting as a lens (High-z Lens 1, HL1, which turned to be one of the most distant lenses found until now) and constrain the strong lensing system. In Sect. 2, we list the set of ancillary data spanning from optical to radio, with both spectroscopic and photometric observations. The analysis of the data resulting in the lensing configuration conclusion is presented in Sect. 3. The depletion time in particular is discussed in Sect. 5. The characterization of RV4 is developed in Sect. 4, and conclusions are given in Sect. 6. A companion paper focusing on a detailed lens model and on the nature of HL1 is in preparation.
Throughout the paper, we use the ΛCDM cosmology of Planck Collaboration et al. (2016) and a Salpeter (1955) IMF.

Data
In this section, we describe the set of data that we obtained for the characterization of the lensing system composed of RV4 and HL1.

Unresolved observations of RV4
2.1.1. Herschel RV4 is detected in the background of the SPIRE (Spectral and Photometric Imaging Receiver) images of IC3521, a Virgo cluster galaxy observed as part of the Herschel Reference Survey (HRS, Boselli et al. 2010) and the Herschel Virgo Cluster Survey (Davies et al. 2010). As shown in Fig. 1, RV4 appears as a particularly bright, red source in SPIRE imaging of IC3521, clearly contrasting with the foreground galaxy.
The FWHMs of the 250, 350, and 500 µm maps of IC3521 are 18.2 , 24.5 , and 36.0 , respectively (Ciesla et al. 2012). RV4 is not resolved in these SPIRE images, and flux densities are extracted using the timeline-based PSF fitting approach which is the most appropriate method for point-like Herschel sources (Bendo et al. 2013). RV4's SPIRE flux densities are provided in Table 2. We refer the reader to Smith et al. (2012) for a complete description of the SPIRE data reduction made as part of the HRS. IC3521 has also been observed with PACS (Photodetector Array Camera and Spectrometer; Poglitsch et al. 2010) at 100 and 160 µm (Cortese et al. 2014). RV4 is not detected in the PACS images, and the derived upper limits (see Table 2) do not provide any useful constraint on the IR spectral energy distribution (SED).

IRAM 30 m/NIKA
After its initial detection by Herschel, we used the NIKA camera (Neel-IRAM-KID-Array; Monfardini et al. 2010) on the 30 m telescope at Pico de Valeta (Spain) to follow up RV4 (234-14, PI: M. Béthermin) at 1.2 and 2 mm. Nine Lissajous scans of 5 minutes were performed to detect the source, and the data were reduced using the standard pipeline. At the IRAM 30m telescope's resolution at these wavelengths (∼12 and 18 at 1.2 and 2 mm, respectively), the RV4 system is unresolved; therefore a standard PSF-fitting extraction is used. Flux densities are provided in Table 2.

GBT/Zpectrometer
We obtained 6 h of GBT 1 /Zpectrometer observations of RV4 (GBT14A-162, PI: L. Ciesla) to measure its spectroscopic redshift. The observations were carried out on March 1 st , 2014. We aimed at the detection of the CO(1-0) line in the frequency range of the Zpectrometer (25.6-36.1 GHz), covering 2.1 < z < 3.5 with a spectral resolution of 32 MHz per channel. No line was detected in the observations. Nevertheless, these observations allowed us to narrow the range of redshift solutions for RV4, and thus pinpoint the true redshift of the source.

IRAM 30m/EMIR
We obtained 6 h of IRAM 30 m telescope DDT time (D07-15, PI: L. Ciesla) aiming at the detection of CO(5-4) and CO(4-3) to measure the redshift of RV4. The observations were made on February, 24 th and March, 1 st of 2016 with EMIR in band E090 (3 mm) over the 80-111 GHz frequency range with a spectral resolution ∆v of 50 km/s. The wide-band line multiple auto-correlator (WILMA) and the fast Fourier Transform Spectrometer (FTS-200) were used simultaneously as backends during the observations. Bright planet/quasar calibrators including Jupiter and J1226+023 were used for pointing and focusing. The weather conditions were excellent with τ 225 GHz 0.2, reaching a sensitivity of 0.6 mK per 50 km/s and a system temperature of 100 K. The data were calibrated using the standard dual method. Data were then reduced with the GILDAS 2 package CLASS. The baseline-removed spectral scans were co-added according to the weights derived from the noise levels of each. We also include ∼ 10% absolute flux calibration uncertainty in our overall uncertainty.

VLA
To spatially resolve the emission of RV4, we obtained 9.5 h (7.4 h on source) of Karl G. Jansky Very Large Array 3 (VLA) time (VLA/2014-06-044, PI: L. Ciesla) in A-configuration with 28 antennas at 3 GHz (S-band, 2.0 GHz-3.9 GHz), reaching a sensitivity of 2 µJy/beam and a beam size of 0.65 . The observations were carried out from August 3 rd to August 22 nd 2015. The data were calibrated by the observatory pipeline. We produced a continuum image using the CASA software (McMullin et al. 2007). We used all the channels to produce a continuum map and maximize the SNR. We used Briggs weighting with robust=0.5 to achieve a good compromise between sensitivity (2.6 µJy/beam) and angular resolution (0.68 ×0.55 ). The flux densities of the detections are provided in Table 2.

NOEMA
We obtained 6 h (3.1 h on source) of PolyFiX NOEMA data (W17EG002, PI: L. Ciesla) to map the CO(5-4) emission at observed-frame 100.6 GHz, in the 3 mm band (A-configuration, 9 antennas, 82.9-90.6 GHz and 98.4-106.1 GHz), with a native resolution of 2 MHz, that we later rebinned by a factor of 8 before imaging. Observations were carried out on February 6 th and 9 th , 2017. We reached a sensitivity of 13 µJy/beam with a spectral resolution of 167.8 MHz per channel, and a spatial resolution of 1.56 ×0.84 . The data were calibrated using the GILDAS/CLIC package. The data cubes and continuum maps were generated using GILDAS/MAPPING. Because of the very large band width of Polyfix, we imaged separately the continuum in the lower and upper side bands. Both CO(5-4) and [CI] (492.161 GHz rest frame) lines are clearly detected (see Sect. 3.1).

Spitzer
Mid-IR data from Spitzer/IRAC (InfraRed Array Camera) are available only at 3.6 and 4.5 µm, acquired as part of the The Spitzer Survey of Stellar Structure in Galaxies (S4G; Sheth et al. 2010), and were downloaded from the NASA/IPAC Infrared Science Archive 4 . We refer the reader to Sheth et al. (2010), Muñoz-Mateos et al. (2013), and Querejeta et al. (2015) for detailed descriptions of the data acquisition and reduction. No observation with IRAC3 and 4, nor with MIPS, is available at the coordinates of RV4.

CFHT data from NGVS
Deep ground-based optical images of the Virgo cluster are available as part of the Next Generation Virgo cluster Survey (NGVS; Ferrarese et al. 2012) obtained with MegaPrime (CFHT). The region around RV4 has been observed in u * , g , i , and z bands. The depths for a point source with SNR of 10 are 25.9, 25.7, 24.9, and 24.6 AB mag, in the u * , g , i , and z bands, respectively. We refer the reader to Ferrarese et al. (2012) and Licitra et al. (2016) for detailed information on the data acquisition and reduction.

Magellan/FOURSTAR
We observed HL1 in March 2018 with the near-infrared imager Fourstar on the 6.5m Magellan Baade telescope using a randomposition dither pattern. Three filters were used: J1 (corresponding to the Y band), J, and Ks, with integrations of 44.8, 38.4, and 15.1 min, respectively. These data were initially calibrated with the FSRED pipeline.

VLT/XSHOOTER
We obtained 6 h of VLT/XSHOOTER (Vernet et al. 2011) time as part of cycle 97A (097.A-0511, PI: T. Diaz-Santos) aiming to achieve an SNR of 5 in the 1.5-2.2 µm range (NIR arm). The observations were carried out on February 1 st to 7 th , 2017. The total effective integration time of 4.8 h and 4.5 h in NIR and VIS arms, respectively, was split between 6 OBs with 4 nodding positions (ABBA) each. Due to a telescope pointing issue, we needed to reject one of the 6 OBs. We reduced each of the remaining 20 nodding pairs (AB) separately using the XSHOOTER pipeline (Modigliani et al. 2010). We flux-calibrated the data using standard pipeline recipes applied to observations of flux standard stars taken during each night of the observations. In addition, we (1-0), CO(4-3), and CO(5-4) emission lines detected with 30m/EMIR at the RV4 position. They correspond to the integrated emission of the three VLA/NOEMA blobs. The red solid lines shows the best fit of Gaussian.
corrected the data for telluric absorption using a model atmospheric transmittance spectrum created with molecfit Kausch et al. 2015) from observations of telluric standard stars taken close in time and airmass to the science observations. Subsequently, we optimally combined the individual nodding pairs (2D spectra) with a weighted average using our scripts. Our scripts also corrected the wavelength scale to vacuum and removed the heliocentric velocity. Finally, we extracted 1D spectra from the resulting 2D master spectrum. The 1.2 slits provide a nominal spectral resolution of R=6500 and R=4300 in VIS and NIR, respectively.
3. An unusual strong lensing system  panel). RV4 is unresolved in the three bands, and the SPIRE flux densities are 78.4±8.9, 118.4±10.9, and 122.8±11.1 mJy at 250, 350, and 500 µm, respectively. The NIKA millimeter flux density ratio excludes the possibility of a low-z radio AGN and confirms the high-z nature of this emission with fluxes 21.1±1.2 and 6.0±0.2 mJy at 1 mm and 2 mm, respectively, placing the IR peak of the SED between 350 and 500 µm. The redshift of RV4 is provided by the detection of the CO(5-4), CO(4-3), and [CI]( 3 P 1 → 3 P 0 ) lines with 30 m/EMIR (Fig. 2), yielding a spectroscopic redshift of 4.72401±0.00042. Combining the SPIRE and NIKA observations and the spectroscopic redshift of RV4, we derive a total IR luminosity of 1.07±0.19×10 14 L from IR SED fitting (see Sect. 4.1). The high S 500 flux density is above the 100 mJy Negrello et al. (2010) lens selection threshold, above which the probability for a SMG to be lensed is very high. This strongly suggests the presence of multiple sources contributing to the submm fluxes (e.g., Hodge et al. 2013) or a magnification from lensing (e.g., Negrello et al. 2010). The SPIRE/NIKA emission is resolved by VLA A-configuration observations at 3 GHz/10 cm into three sources, hereafter named A, B, and C ( Fig. 1, right panel). The measured 10 cm continuum flux densities are 34.3±4.6, 31.2±4.9, and 19.8±6.3 µJy for sources A, B, and C, respectively.
We determine the spectroscopic redshift of each of the three blobs from NOEMA/PolyFix observations using the CO(5-4) emission line. The [CI](1-0) line is also detected but only for blob A and B; the [CI] to CO(5-4) line ratio is equivalent for both blobs, at 0.27±0.04 and 0.26±0.04 for blob A and B, respectively. For each blob, and for the CO(5-4) and [CI](1-0) lines, we measure the redshift using the slinefit 5 code (Schreiber et al. 2018). As shown in Fig. 3, the CO(5-4) lines of blob A and B are at nearly the exact same frequency (100.683 GHz and 100.680 GHz) with the same width (627±55 km/s and 622±67 km/s). The redshift of blob C is slightly offset with a CO(5-4) line at 100.66 GHz. These frequencies translate into z = 4.72359 +0.00022 −0.00011 , z = 4.72375 +0.00013 −0.00023 , and z = 4.72469 +0.00028 −0.00028 , for blobs A, B, and C, respectively. The errors quoted in the redshift measurements are underestimated, as the uncertainties from the line profiles are not included. Indeed, as discussed in different studies based on simulations (e.g., Hezaveh et al. 2012;Serjeant 2012) and observations (e.g., Riechers et al. 2008;Dye et al. 2015;Yang et al. 2017Yang et al. , 2019aApostolovski et al. 2019), the magnification can significantly vary along the velocity channels from blue to red and thus heavily distort the intrinsic line profile. This differential lensing effect may also cause different line profiles in different images. This effect could be related to the redshift offset of blob C, which might be caused by a higher magnification of the blue part of the line. Although the [CI](1-0) line is weaker, the line frequencies of blobs A and B are consistent with a single redshift. As for CO(5-4), the [CI](1-0) redshift of blob C seems to be offset, but in this case the line SNR is too low (≈2) for a meaningful comparison.
In Fig. 4, we show that the positions of the NOEMA blobs are consistent with the VLA continuum positions. The CO(5-4) positions of blobs A and B are also consistent with the continuum detections. However, there seems to be a shift between the position of the CO(5-4) emission from blob C and its continuum counterparts. The separation between the two centers of emission is 0.29 ±0.10 . In theory 6 , for the CO(5-4) line, a 4.5 sigma detection has a position uncertainty of 0.17×0.09 . The same uncertainty is expected for the continuum. Therefore, this spatial shift might not be significant.
Considering the extreme IR luminosity of the system, the spatial distribution of the three blobs,and the identical CO(5-4) line profile (redshift and width) for each blob, we conclude that RV4 is a lensed z = 4.72 submm galaxy, with A, B, and C, being multiple images of the same galaxy. 5 https://github.com/cschreib/slinefit 6 see https://www.iram.fr/IRAMFR/IS/IS2002/html 2/node131.html

Identification of the lens
Two sources, hereafter denoted as High-z Lens 1 West (HL1-W) and High-z Lens 1 East (HL1-E), are clearly detected in the u * , g , i , and z bands. The two sources, which are only 2" apart, are slightly offset from the three VLA and IRAM/NOEMA blobs, by approximately 2.2" to the north (Fig. 5, left panel). HL1-W and HL1-E are also clearly detected in IRAC imaging at 3.6 and 4.5 µm (Fig. 5, right panel).
The u * and g bands probe shorter wavelengths than the Lyman break at z = 4.72; thus we do not expect any emission from RV4 in these bands. Moreover, as shown by the VLA contours on the same figure, no i or z emission is seen from these three blobs either. The extended source detection limits of the NGVS are 26.3 and 25.8 AB mag arcsec −2 (2σ) for the i and z band, respectively (Ferrarese et al. 2012). In the Spitzer/IRAC images, although the VLA detections are close to the outskirts of the IRAC emission of the lens system, no strong emission from RV4 is detected (Fig. 5, right panel). Furthermore, the IRAC fluxes measured from PSF fitting for HL1-W and HL1-E are consistent with its SED, indicating no particular excess of flux that could be attributed to RV4.

Redshift determination for the lens
To determine the redshift of HL1-W and HL1-E, we obtained VLT/XSHOOTER observations of the system (Fig. 6). The spectroscopic redshifts are measured using the software slinefit. For HL1-E and HL1-W respectively, we obtain z spec =1.48379 +0 as well as a detection of the continuum emission. These values reveal HL1-E/W as one of the most distant lens known in strong lensing systems and correspond to a velocity difference of 343±38 km/s indicating that HL1-E and HL1-W are probably merging together. The projected distance between the two blobs of 2 (17.4 kpc) is also consistent with merger scenario. In the literature, we find only one lens at a higher redshift presented in Cañameras et al. (2017) at z spec = 1.525, which is gravitationally lensing an SMG detected by Planck.

Lensing properties
The data in hand allow us to obtain a better picture of the system in which RV4, at z spec = 4.72, is strongly lensed by a close pair of z spec = 1.48 galaxies. The redshifts of these sources as well as the positions of the multiple images of RV4 (A, B, and C) are used as constraints to provide a first order lens model. The redshifts of the lenses and sources as well as the centroids of images A, B and C are used to constrain a mass distribution assuming they arise from a multiply-imaged system. We use the lenstool 7 software (Kneib et al. 1996;Jullo et al. 2007) to optimize a parametric model of the mass distribution reproducing the locations of the lensed images. The model parameters are presented in Table 3. The modeling have been made following a Occam's razor principle in the fact that this the simplest model able to correctly reproduce the positions of the VLA/NOEMA 7 Publicly available at https://projets.lam.fr/projects/lenstool/wiki/ blobs as well as their flux ratios. A thorough model comparison will be presented in the second paper. Although the lens is composed of a pair of galaxies, we find that a simple model with a single dPIE (double Pseudo-Isothermal Elliptical) mass distribution provides a very good match, with an rms of 0.01 in the image plane. This simple model estimates the magnification to be 8.2±2.5 with an Einstein radius, constrained by the positions of blobs A, B, and C, of 2.2 ±0.2. The model provides magnifications of 3.4±1.9, 2.6±1.4, and 2.3±0.8 for blobs A, B, and C, respectively.
We compare the redshifts and Einstein radii of other strong lensing systems found in the literature in Fig. 7. In addition to being one of the most distant lenses known, HL1 has an Einstein radius of 2.2 ±0.2 (20 kpc at z =1.48) that places it in a new region of the Einstein radius versus redshift plane (Fig. 7).

Physical properties of the lens
To derive the physical properties of each component of HL1, we combine the photometric data (CFHT, FOURSTAR, Spitzer) with the XSHOOTER spectra. The XSHOOTER spectra are scaled using the i band flux density for each blob and then integrated into a set of 10 artificial filters. We use the SED modeling code CIGALE (Boquien et al. 2019) to perform the SED fit-ting. We assume a Salpeter (1955) IMF, a flexible delayed SFH (as described in Ciesla et al. 2017), Bruzual & Charlot (2003) models, and a Calzetti et al. (2000) attenuation law. We find that HL1-W and HL1-E have stellar masses of 2.4±0.4×10 11 and 3.3±0.6×10 11 M , respectively. The SED fitting performed on the whole system using the integrated photometry of the two components yields a total stellar mass of 6.2±1.0×10 11 M .
Although the stellar mass is well constrained thanks to FOURSTAR and Spitzer data, the star formation rate (SFR) estimates for HL1-W and HL1-E are more uncertain due to the faintness of the two galaxies. The SED fitting results in SFRs of 9.2 +46.0 −9.2 and 3.7 +5.1 −3.7 M yr −1 for HL1-W and HL1-E, respectively, where the errors reflect the difficulty of estimating the SFRs from the faint g and i emission. These values place the two counterparts below the MS of Schreiber et al. (2015) at z = 1.5, as shown in Fig. 8. However, we can probe the SF activity of these two sources using the Hα emission detected in the XSHOOTER spectra. As shown in Fig. 6, HL1-W has no Hα emission, confirming the fact that HL1-W seems to be passive. On the other hand, HL1-E seems to have a large and strong Hα line with a line flux of 1.79±0.33×10 −17 erg s −1 cm −2 associated with a strong [NII] line with a flux of 4.33±0.38×10 −17 erg s −1 cm −2 . Furthermore, the Hβ emission line is not detected whereas a strong [OIII] line with a flux of 2.40±0.74×10 −17 erg s −1 cm −2 is detected. These line fluxes, as well as the Hβ upper limit, place HL1-E in the AGN region of the BPT diagram, as shown in Fig. 9. Since HL1-E is likely hosting an AGN, we check if our measured stellar mass can be contaminated by AGN emission in the NIR. To do this, we use CIGALE to quantify the AGN fraction, defined as the contribution of the AGN to the total L IR of the galaxy, which is a scaling factor for the AGN emission model. The methodology is fully described in Ciesla et al. (2015). For HL1-E, the best χ 2 model yields f rac AGN = 0, meaning that an AGN emission component is not needed to reproduce the host galaxy UV-NIR emission. Therefore our stellar mass measurement is not contaminated by the AGN emission. We conclude that HL1-W and HL1-E seem to be passive, as no strong star formation activity is detected from either continuum or line analysis.
For the SED modelling, we adopt a Calzetti et al. (2000) attenuation law, which is better suited for star-forming galaxies than for passive galaxies. However, Calzetti et al. (2000) assume almost no attenuation in the NIR. Other attenuation laws such as Charlot & Fall (2000) and

Infrared properties
The combination of the three SPIRE detections with the NIKA 1 mm and 2 mm flux densities provides a good sampling of the  Myr 32.9±18.7 µL' CI 10 10 K km s −1 pc 2 9.6±2.6 µL' CO(4−3) 10 10 K km s −1 pc 2 17.8±2.3 µL' CO(5−4) 10 10 K km s −1 pc 2 17.1±1.9 Σ IR 10 11 L kpc −2 7.7±1.6 IR SED of RV4. Using the SED modeling code CIGALE, we fit the IR SED using two different dust emission models. It has been found that ∼20% of SMGs host an AGN (e.g., Coppin et al. 2010); therefore, we investigate if our data allow us to rule out the presence of an AGN contributing to the IR luminosity of RV4. The SMG is not dectected in the PACS images of IC3521 at 100 and 160 µm. The 3σ upper limits determined as in Ciesla et al. (2012) and Cortese et al. (2014) do not place any useful constraints on the IR SED of RV4, as shown in Fig. 10. However, we have retrieved WISE 12 and 22 µm maps of the region around RV4. The SMG is not detected, so we use the  3σ detection limits provided by the explanatory supplement 8 . Different models were assumed, and thus different values are quoted for the sensitivity limits. We assume RV4 to be unresolved in WISE bands, and take the less constrictive upper limit provided in the explanatory supplement, i.e., 0.52 and 3.24 mJy  (3σ) at 12 and 22 µm, respectively, as shown in Fig. 11. Using the best fit of the IR SED obtained in Sect. 4.1, we add different AGN components assuming a large range of AGN fraction, defined here as the AGN contribution to the total IR luminosity (Ciesla et al. 2015). As shown in Fig. 11, our data do not allow us to exclude the presence of a dust obscured, type 2 AGN that could contribute to the IR luminosity of RV4. However, the two WISE upper limits indicate that, were RV4 hosting a type 1 or intermediate type of AGN, we should expect a contribution to the total L IR less than 28%. We note that at these wavelengths, considering the redshift of RV4, we may need to take the stellar emission into account, therefore, the 28% contribution should be considered as an upper limit. As a final test, we determine the IR/radio coefficient, q IR , to understand if RV4 could be a radioloud AGN. Using the VLA continuum data point (Fig. 10), we obtain a q IR =2.70±0.09 using the CIGALE code. This is a typical value for star-forming galaxies (Seymour et al. 2009;Sargent et al. 2010). We therefore conclude that RV4 is not a radio-loud AGN, although the data in hand do not allow us to place a further constraint on the possible presence of a dust enshrouded AGN.

Gas mass and depletion time
Given that we detect high-J transitions of CO, i.e. CO (5→4) and CO(4→3), deriving a gas mass requires that uncertainties in the spectral line energy distribution (SLED) must be forded in with the uncertainty on the adopted α CO value. Therefore, we first use the [CI] flux measurement to derive the H 2 mass. Indeed, Bothwell et al. (2017) showed that [CI] and CO(2-1) have similar kinematic properties in a sample of dusty star-forming galaxies, suggesting that [CI] traces the same gas component as low-J CO emission. Furthermore [CI] as a gas mass tracer is far less affected by metallicity than CO (e.g., Bothwell et al. 2013), although it does require an assumption that the line is optically thin.
We follow the prescription of Papadopoulos & Greve (2004): and where h =0.75, the Einstein A-coefficient A 10 =7.93×10 −8 s −1 , X CI =3×10 −5 , and Q 10 =Q 10 (n, T k ) depending on n and T k the gas density and temperature (Weiß et al. 2003;. Assuming the abundance and excitation factors of Alaghband-Zadeh et al. (2013), we find µM CI H 2 =1.85±0.51×10 11 M , with µ the magnification factor due to lensing, which yields µM CI gas =2.65±0.73×10 11 M , assuming a 36% helium contribution and no significant atomic H component as in Yang et al. (2017). Converting the L IR to SFR we obtain µSFR=1.82×10 4 M yr −1 and derive a depletion time as: finding t dep =14.4±4.7 Myr. This is a short depletion time highly suggestive of a high-efficiency, bursty star formation, typical for SMGs (e.g., Aravena et al. 2016). However, as the gas mass measurement obtained from CI depends on the assumption that CI is optically thin, to provide an independent check, we also derive the gas mass using the CO(4-3) line (Table 4). Assuming L CO(4−3) /L CO(1−0) = 0.46 (Carilli & Walter 2013), we obtain µL CO(1−0) = 3.9 ± 0.5×10 11 K km s −1 pc 2 . Using α CO = 0.8 M /K km s −1 pc 2 , we obtain µM CO gas = 3.1 ± 0.4×10 11 M , which is consistent with the estimate obtained from CI. This gas mass estimate yields a depletion time of 17.2±2.2 Myr. Although the uncertainty of this estimate does not reflect the uncertainty in the L CO(4−3) /L CO(1−0) ratio and α CO assumptions, it is consistent with the value obtained from CI.
Finally, we use the dust mass as a third independent gas mass indicator: with δ GDR the gas-to-dust mass ratio. We assume here the average δ GDR of 56 (±28) obtained by Yang et al. (2017) for a sample of lensed SMG. We obtain µM dust gas =6.0±3×10 11 M and from this a t dep of 32.9±18.7 Myr. Thus, using the dust mass, we obtain a larger, yet still small value of the depletion time indicating a starbusting phase. Using the δ GDR obtained from the ALESS sample Swinbank et al. 2014), 75±10, would yield a depletion time of 44.0±5.9 Myr. However, the sample of lensed SMGs of Yang et al. (2017) has physical properties that are close to those of RV4, therefore we use the gas mass obtained using the gas-to-dust ratio of Yang et al. (2017) in the following. These derived values are short, comparable to those of other starburst galaxies, and pointing towards a rapid starburst episode. We note that when computing the depletion times, we assume that the magnification µ is the same for the gas component and the IR produced by star-formation since we cannot make a detailed investigation of differential magnification with the data in hand.
Given the good coverage of RV4's IR SED with 6 data points from 250 µm to 3 mm observed (44 to 524 µm, rest frame), we derived the dust mass of RV4 from SED modeling. We used an updated version of the model of Draine & Li (2007) (Draine et al. 2014), which is physically motivated and takes into account the different contributions from the dust heated in PDR and the dust heated by more evolved stars. As an alternative, Scoville et al. (2016) proposed a method useful when only one observations is available in the Rayleigh-Jeans part of the SED but that relies on assumptions such as the dust temperature that can lead to large uncertainties (Berta et al. 2016). Given our good coverage of both the peak of the IR emission and the Rayleigh-Jeans part of the SED, we prefer to rely on our measurement based on SED modeling. However, for comparison, we compute the gas mass of RV4 obtained by the Scoville et al. (2016) method and obtain 7.7±2.0×10 11 M , uncorrected from magnification. This leads to a depletion time of 42.8±11.1 Myr.

Gas density
In this section, we use the computed IR, CO(4-3), and [CI](1-0) luminosities of RV4 obtained from EMIR (summing the emission of all three blobs) to investigate its gas properties. First, we present in Fig. 12 (top panel) the relation between [CI](1-0) to IR ratio and the IR luminosity. For comparison, we add a selection of data from the literature (Walter et al. 2011;Alaghband-Zadeh et al. 2013;Béthermin et al. 2016;Bothwell et al. 2017;Lu et al. 2017;Cañameras et al. 2018;Dannerbauer et al. 2018;Valentino et al. 2018;Nesvadba et al. 2019). The position of RV4 is marked by the red star, taking into account the magnification from lens modeling. Even after this correction, RV4 is one of the most luminous IR sources of the compilation of SMGs shown in Fig. 12 (top panel). Its L IR is comparable to those of the recently published sources of Béthermin et al. (2016) (L IR = 1.1 ± 0.2 × 10 13 L ) and of Dannerbauer et al. (2018) (L IR = 1.1 ± 0.1 × 10 13 L ). Although the CO32-A source studied in Dannerbauer et al. (2018) has a [CI](1-0) to IR luminosity ratio consistent with sources with a less extreme L IR , RV4 and the source analyzed in Béthermin et al. (2016) exhibit a lower ratio.
The locations of RV4, CO32-A, and the source of Béthermin et al. (2016) in the upper panel of Fig. 12 seems to indicate a possible deviation from the main trend at high IR luminosity, especially if we also consider the z > 2.5 SMG source of Valentino et al. (2018) with the lowest L CI /L IR ratio. However, more statistics is needed to confirm this. The [CI] luminosity traces the total gas while the IR luminosity traces the SF activity over a typical scale of 100 Myr. Therefore, this ratio can be interpreted as a star formation efficiency indicator. The low L CI /L IR ratio is thus an indication that star formation is highly efficient in RV4.
To explore the properties of the ISM of submillimeter galaxies, Alaghband-Zadeh et al. (2013) investigated the position of galaxies in the L CI /L CO(4−3) versus L CI /L IR diagram. Although they acknowledged issues with simple PDR modeling, they compared a sample of quasars and SMGs to a set of PDR models from Kaufman (2009) and found that the SMGs have densities and radiation field strengths that are consistent with those of local starbursts, although ∼35% of their SMG sample has lower density and weaker radiation field. They interpret this as an evidence that star formation can be extended in some submillimeter galaxies. The quasar sample is found at higher densities and radiation fields than SMGs on average. In Fig. 12 (bottom panel), we plot RV4 with the two samples of Alaghband-Zadeh et al. (2013) on the L CI /L CO(4−3) versus L CI /L IR diagram, along with the SPT data of Bothwell et al. (2017). In terms of L CI /L CO(4−3) ratio, RV4 is consistent with what is found for SMGs, although the ratio is among the 25% highest for the two samples shown in Fig. 12 (lower panel). The position of RV4 in this diagram is close to the log(n H /cm −3 ) = 4 density line, indicating a rather low gas density compared to the other SMGs for the sample. However, in terms of L CI /L IR , RV4 has the lowest ratio of the SMG sample, close to the values obtained for quasars. This is an indication of a rather intense radiation field, which is consistent with RV4 being a starburst. This low gas density is also confirmed by the L CO(5−4) /L IR ratio for RV4, which corresponds to only 54±11% of the expected value for normal galaxies .

Size and IR surface brightness
In the VLA observations, with a beam size of 0.68 × 0.55 and a position angle of 30 degrees, all three blobs are compatible with being point sources. The emission from blob A and B can be deconvolved from the beam but with large uncertainty. Given the sensitivity (2.6 µJy/beam) of the VLA observations, measuring an accurate size is challenging, and we cannot exclude that the sources are marginally extended. However, we can measure the sizes of the three RV4 blobs in the uv plane of the NOEMA data, combining all continuum channels. Assuming a circular Gaussian profile, the resulting observed sizes are 0.59 ±0.05 , 0.50 ±0.05 , and 0.50 ±0.08 for blobs A, B, and C, respectively. All sources are clearly resolved and, within the uncertainties, have approximately the same size. We note that these sizes are not corrected for the shear that the strong lensing is causing. We also tried an elliptical Gaussian profile fit, but only the brightest source, blob A, has enough signal to allow for a meaningful fit. For blob A, we obtain an observed major axis of 0.64 ±0.23 , an observed minor axis of 0.31 ±0.6 , and a position angle of -75.8±14.9 degrees. The observed major axis is elongated approximately E-W, along the direction where the lens model predicts the greatest shear. The lens model applied in Sect. 3.4 provided the estimated magnification for each blob along the major and minor axis. Correcting the sizes obtained from CO(5-4) fits in the uv plane, we obtain a mean delensed size for RV4 of 0.40 ±0.09 × 0.26 ±0.06, corresponding to a physical size of 2.64±0.58 kpc × 1.71±0.38 kpc. RV4 is thus relatively compact.
Although lensing applies a shear to RV4's morphology, the IR surface brightness is conserved. Adopting the lens model that best reproduces the positions of RV4's blobs as well as their flux ratios, we distribute the IR luminosity among the three blobs taking into account their respective modeled magnifications. Using the continuum sizes derived from the NOEMA data, we obtain surface brightnesses of Σ IR = 9.5 ± 2.0 × 10 11 L kpc −2 , 7.1 ± 1.5 × 10 11 L kpc −2 , and 6.5 ± 1.4 × 10 11 L kpc −2 for blobs A, B, and C, respectively. These values are consistent within the error bars and provide an average Σ IR of 7.7±1.6×10 11 L kpc −2 for RV4. In Fig. 13, we compare RV4's Σ IR to those of other SMGs found in the literature.
RV4 has a notably low depletion time despite a Σ IR only slightly smaller than the sample median. Different assumptions could put RV4 closer to the median along both axis: the system's depletion times obtained from dust mass and gas-to-dust ratios are longer than the ones obtained from both [CI] and CO(5-4). Furthermore, if an obscured AGN is contributing to the L IR then the depletion time for RV4 would be higher. However, we do not expect these effects could erase all of the current shortfall relative to the Elbaz et al. (2018) MS relation, although the evolution of such relation with redshift is not yet understood. We  (Nesvadba et al. 2019), and the local LIRG sample of Lu et al. (2017). Lower panel: [CI] to CO(4-3) luminosity ratio versus [CI] to IR luminosity ratio (L ). The position of RV4 is marked by the red star. This figure is adapted from Alaghband-Zadeh et al. (2013). Their SMG sample is shown (green triangles) as well as its median position (dark green triangle). The blue squares indicate the quasar sample presented in Walter et al. (2011) with the dark blue square indicating the position of the median of their sample. The SPT sample of Bothwell et al. (2017), Planck's dusty GEMS (Cañameras et al. 2018;Nesvadba et al. 2019), and the local sample of LIRGs studied in Lu et al. (2017) are also shown in the figure. The grey lines indicate the contours of the gas density (n) and and the radiation field (G 0 ) for the corresponding L CI /L CO(4−3) and L CI /L IR ratios, as produced by the PDR models of Kaufman et al. (1999). also note that our method uses the size derived from 3 mm continuum emission which traces the relatively cold dust component (∼520 µm, rest frame) whereas the bulk of the IR emission comes from star-forming regions. Such an offset could have an impact on our derived Σ IR although it is difficult to quantify with the data we have in hand.

Interpreting the short depletion time
The depletion times estimated through different assumptions ([CI], CO, and dust) yield rather low values 14.4±4.7, 17.2±2.2, and 32.9.0±18.7 Myr, respectively, strongly suggesting that RV4 is experiencing a starburst event. As shown in Fig. 12, the L CI /L IR ratio is very low compared to a literature sample and seems to indicate a high intensity of the radiation field, close to those measured in quasars. In Fig.14, we place RV4 in a Kennicutt-Schmidt diagram (Kennicutt 1998) along with SPT sources from Bothwell et al. (2010) whose properties are close to those of RV4. The Σ gas determined from the CO and [CI] lines places RV4 on the relation of Daddi et al. (2010b) for ULIRGs and SMGs. The Σ gas determined from M dust is larger but still above the relation of Daddi et al. (2010b) for BzK galaxies but right on the universal relation determined by Bouché et al. (2007). The position of RV4 in this diagram confirms the starbursting phase suggested by the very short depletion time.
Even considering a starburst, the typical depletion timescale of a starburst being 100 Myr (e.g., Béthermin et al. 2015;Aravena et al. 2016), RV4's t depl is at least a factor 2 lower. A first explanation to this short depletion time could be that RV4 has an unusually high star formation efficiency, defined as the star formation rate divided by the gas mass. However, such high star formation efficiencies are hard to explain theoretically and it has been shown in simulations that major mergers are not that efficient in producing high bursts of star formation in massive gasrich galaxies . A second possibility would be that RV4 encounters a "recent" and rapid decrease of the SFR. In this case, the L IR would trace the SFR over a typical timescale of 100 Myr (Kennicutt & Evans 2012) while the gas mass is probed  by emission lines which trace the current gas content. If the SFR and the gas content are both rapidly decreasing, the L IR would overestimate the instantaneous SFR and thus bias the depletion time toward lower value. Therefore, if the decrease of star formation is very fast the actual depletion time would be closer to typical starburst values. If we assume that RV4 lies at least a factor of three above the MS, the stellar mass will be at most 1.7×10 11 M . Considering the three estimates of the gas mass derived in Sect.4.2, corrected for magnification, we find gas fractions larger than 16%, 18%, and 30% considering the [CI], CO, and dust based estimates, respectively. These values are not constraining enough for SMGs for which gas fractions are ∼50% (e.g. Daddi et al. 2010a;Tacconi et al. 2010); therefore we cannot disentangle between the two possible explanations for the short depletion time. Observations of RV4 with the JWST to probe its stellar content will be key to disentangle between the two scenarios.

Conclusion
We have serendipitously discovered a bright submillimeter galaxy in the field of view of the Herschel imaging of the nearby galaxy IC 3521. The rising SPIRE colors and high flux densities indicate a high-z galaxy with high IR luminosity. We have gathered ancillary data to shed light on the nature of this source and reached the following conclusions.
The redshift of RV4 is 4.724 as constrained from IRAM/EMIR and NOEMA data. The combination of Herschel and IRAM/NIKA data constrains the IR SED of RV4 resulting in a total IR luminosity of 1.06±0.6×10 14 L . However, the data in hand do not allow us to exclude the presence of a dust enshrouded AGN that could be contributing to this IR luminosity.
The Herschel emission splits into three blobs in VLA at 10 cm imaging. The IRAM/NOEMA detection of the CO(5-4) line of each blob confirms that the three sources are at the same redshift with CO(5-4) lines having the same width. Combined with the extreme L IR , we conclude that RV4 is lensed.
In Spitzer/IRAC and deep CFHT data, two sources acting as a gravitational lens are detected at the center of the virtual arc formed by the three VLA/NOEMA blobs. These two sources, High-z Lens 1 West and East (HL1-W and HL1-E), are located at z =1.48 by XSHOOTER data, indicating that they are probably merging. They constitute one of the most distant gravitational lenses found to date in the literature for strong lensing systems. The Einstein radius of the lensing system is 2.2 ±0.2, as determined from the positions of the three VLA/NOEMA blobs. HL1 is a very peculiar lens with a large Einstein radius combined with a high redshift. We derive a lens model and find that a single halo best reproduces both the positions and the flux ratios of the VLA/NOEMA blobs. This model yields a total magnification of 8.2±2.5 and thus an intrinsic L IR of 1.29×10 13 L for RV4. A detailed model of the system will be presented in a companion paper (Ciesla et al., in prep).
The SED modelling of HL1-W and HL1-E yields stellar masses of 2.4±0.4 and 3.3±0.6×10 11 M , respectively, and SFR compatible with the two galaxies being passive. The absence of Hα emission in HL1-W confirms the low SFR obtained from SED modeling. For HL1-E, a strong Hα line is detected together with a strong and broad [NII] line, indicating the possible presence of an AGN. A detailed analysis of the HL1 system will be provided in the companion paper (Ciesla et al., in prep) as well.
From the IRAM/EMIR and NOEMA data, we have measured the [CI], CO(4-3), and CO(5-4) apparent luminosities. On the one hand RV4 shows a relatively low L [CI] to L IR ratio, which can be interpreted as a high star formation efficiency compared to other high-z samples of galaxies, and/or with a hard radiation field. On the other hand, the L [CI] to L CO(4−3) ratio of RV4 is relatively high, indicating a lower gas density. Furthermore, the L CO(5−4) to L IR is only half of the value expected for normal galaxies, indicating a lack of gas.
We estimate the gas mass of RV4 from different tracers and obtain between 2.7±0.7×10 11 M and 6.0±3.0×10 11 M , not corrected for magnification. These numbers yield depletion times between 14.4±4.7 and 32.9±18.7 Myr. Such short values can be explained by either very high star formation efficiency or by a rapid and recent decrease of star formation. In the former scenario, RV4 would be an intense starburst, difficult to understand from the results of major mergers simulations of gas rich high-z galaxies. In the second scenario, the L IR probing the star formation activity on a timescale of ∼100 Myr and the gas content following the instantaneous SFR, if the star formation activity of RV4 is rapidly decreasing it would bias t depl toward low values. JWST will definitely be needed to determine the M gas to M * ratio of RV4 and investigate the possible low gas content of this SMG.