© ESO, 2016

1. Introduction

The measurement of the stellar mass function (SMF) is a powerful statistical tool for tracing the stellar mass assembly or galaxy growth over cosmic time. Galaxy formation models rely on the well established ΛCDM cosmological framework that governs the growth of the dark matter structures and the less well understood baryonic physics at play inside the dark matter haloes (gas accretion, minor or major merging, star formation activity, feedback mechanisms, etc.). The shape of the galaxy SMF compared to the expected halo mass function provides valuable information about the physical processes acting at the low- and high-mass ends of the mass function (Silk & Mamon 2012).

A decade ago, early deep extragalactic surveys have revealed that the average stellar mass density decreased gradually (the integrated form of the SMF) from z ~ 3 to z ~ 0 (e.g. Dickinson et al. 2003; Fontana et al. 2003). This trend is now confirmed up to redshift z ~ 8 (Song et al. 2015) and is consistent with a hierarchical build-up of the cosmic structures. Later on, larger surveys have measured the evolution at high redshift of the galaxy bimodality, the well-known separation between star-forming and quiescent galaxies observed in the local Universe (Baldry et al. 2004; Moustakas et al. 2013). They found that the bimodality was already in place at z ~ 1 with the quiescent galaxies dominating the massive end of the SMF and the star-forming galaxies dominating its low-mass end (Bundy et al. 2006; Borch et al. 2006). This quiescent population had its main build-up epoch between z = 2 and z = 1, where the stellar mass density increased by a factor 10 (Cirasuolo et al. 2007; Arnouts et al. 2007), while only a factor 2 increase is observed from z = 1 to z = 0 (Bell et al. 2004; Faber et al. 2007). According to the hierarchical scenario, such an early formation epoch of the quiescent population was not a problem as long as the stars formed before this in smaller units and galaxies continued to assemble their masses at a later time (through dry merging phases, e.g. De Lucia et al. 2006). This is a natural support of the star formation downsizing picture proposed by Cowie et al. (1996), where the onset of star formation begins earlier for most massive galaxies than for lower mass galaxies (see also Gavazzi & Scodeggio 1996). However, the models predict a continuous increase of stellar mass for these massive galaxies with cosmic time (e.g. De Lucia & Blaizot 2007), which is challenged by the last measurements of the SMF, where the massive end does not show significant evolution from z = 0 up to redshift z ~ 1 (e.g. Marchesini et al. 2009; Ilbert et al. 2013; Muzzin et al. 2013; Moustakas et al. 2013; Mortlock et al. 2015), suggesting a mass assembly downsizing.

The predominance of quiescent galaxies at the massive end (e.g. Baldry et al. 2012; Moustakas et al. 2013; Ilbert et al. 2013) supports the idea that the star formation activity is preferentially impeded in galaxies above a given stellar mass or a given dark matter halo mass, if we assume a stellar-to-halo mass relationship (e.g. Coupon et al. 2015). A wide variety of quenching mechanisms have been proposed to explain the star formation quenching in massive galaxies, such as major mergers (Barnes 1992), virial shock heating (Kereš et al. 2005), or radio-AGN feedback (Croton et al. 2006; Cattaneo et al. 2006) in massive haloes.

Several studies have emphasised the role played by the environment for the colour-bimodality and star-formation quenching in the local Universe (Hogg et al. 2003; Kauffmann et al. 2004; Baldry et al. 2006; Haines et al. 2007). Mechanisms such as ram-pressure stripping, in which the gas is expelled from the galaxy (Gunn & Gott 1972), or strangulation, in which the cold gas supply is heated and then halted (Larson et al. 1980), can be invoked as environmental quenching mechanisms. We emphasise that strangulation processes can either be linked to environment (e.g. when a galaxy enters the hot gas of a cluster) or to peculiar evolution (e.g. when a the radio-AGN feedback stops the cold gas infall). The latest measurements of the quiescent SMFs reveal an upturn at the low-mass end in the local Universe (Baldry et al. 2012; Moustakas et al. 2013), whose build-up is observed at higher redshift (Drory et al. 2009; Tomczak et al. 2014). This upturn in the low-mass end for quiescent galaxies could be associated to environmental quenching according to Peng et al. (2010, 2012), while Schawinski et al. (2014) suggested a fast process consistent with major merging. Constraining the quenching timescale at different masses might therefore help to highlight the quenching mechanisms.

Until recently, the above conclusions were mostly based on deep galaxy surveys such as GOODS (Giavalisco et al. 2004), VVDS (Le Fèvre et al. 2005), COSMOS (Scoville et al. 2007), and DEEP2 (Newman et al. 2013), which are perfectly suited to provide the global picture of the galaxy stellar mass assembly over a wide range of redshifts. However, given their small angular coverage (they explore a rather small volume below z < 1), they can be particularly sensitive to statistical variance (i.e. cosmic variance) at low redshift. This is particularly crucial for the very rare galaxies at the high-mass end of the exponentially declining SMF, and it has been claimed that its apparent lack of evolution may be dominated by observational uncertainties (Fontanot et al. 2009; Marchesini et al. 2009).

A first attempt to constrain the density evolution of the high-mass galaxies at z < 1 has been performed by Matsuoka & Kawara (2010). They combined the SDSS southern strip (York et al. 2000) and UKIDSS/LAS survey (Lawrence et al. 2007) over a total area of ~55 deg². They observed a mild-to-high increase of the number density of massive galaxies 10^11−11.5 M_⊙/10^11.5−12. M_⊙ with a corresponding drop of the fraction of star-forming galaxies in this stellar mass range from z ~ 1 to z ~ 0. While subject to large uncertainties in their photometric redshifts, stellar mass estimates, and reliability of the separation into quiescent and star-forming galaxies, this first result suggested that massive galaxies (M_∗ > 10¹¹ M_⊙) evolve since z ~ 1.

Moustakas et al. (2013) estimated the SMF between 0 < z < 1 over an area of ~5.5 deg² using PRIMUS, a low-resolution prism survey (for galaxies with i_AB ≤ 23; Coil et al. 2011). The wealth of multi-wavelength information from deep ultraviolet (GALEX satellite) to mid-infrared (Spitzer/IRAC) photometry allowed them to derive accurate stellar masses and a reliable separation between active and quiescent populations. Their SMF measurements confirmed the modest change in the number density of the massive star-forming galaxies (M_∗ ≥ 10¹¹ M_⊙), leaving little room for mergers, but observed a significant drop (50%) of the fraction of active star-forming galaxies since z ~ 1 that is in contrast with the classical picture, in which the star-forming population remains constant across cosmic time. Another major spectroscopic sample is provided by the VIMOS Public Extragalactic Redshift Survey (VIPERS; Guzzo et al. 2014), whose first ~50 000 galaxies down to i_AB = 22.5 over an area of 10.3 deg² have recently been released (PDR1, Garilli et al. 2014). Using the PDR1 combined with CFHTLS photometry and the same ultraviolet (UV) and near-infrared (NIR) data that we used here, Davidzon et al. (2013) produced the most reliable overall measurement of the high-mass end of the SMF in between 0.5 < z < 1.3 to date. The VIPERS SMF shows to high precision that the most massive galaxies had already assembled most of their stellar mass at z ~ 1, but that a residual evolution is still present. However, as discussed in Davidzon et al. (2013), although these two studies use spectroscopic redshifts, multi-wavelength information, and a large area, they disagree slightly concerning the general amplitude of the SMF. These discrepancies might be due to differences in the stellar mass estimates, for example, or to selection effects that are not fully accounted for. It highlights how subtle effects become crucial and can introduce significant systematic errors when statistical uncertainties are reduced so drastically.

In this paper we exploit the broad photometric coverage assembled over the footprint of VIPERS to build a unique multi-wavelength photometric sample covering more than 22 deg² down to K_s < 22, as part of the VIPERS-Multi Lambda Survey (VIPERS-MLS; see Moutard et al. 2016). We benefit of the synergy with the VIPERS spectroscopic survey by using the PDR-1 data to compute reliable photometric redshifts, and we derive stellar masses for 760 000 galaxies out to z = 1.5; This allows us to obtain a new estimate of the SMF that (a) has greater control over the low-mass slope because of the i < 23.7 /K_s < 22 depth of our sample for extended sources (more than 1 mag deeper in the i-band than VIPERS); (b) extends over a wider redshift range than VIPERS, from z = 0.2 out to z = 1.5; (c) is less affected by the cosmic variance because the effective area is doubled with respect to the VIPERS PDR-1 used in Davidzon et al. (2013) (we cover nearly the entire footprint of the final VIPERS survey and avoid the 30% area loss that is due to the detector gaps in VIPERS); (d) suffers a reduced Poisson error because our sample is ten times larger in the common redshift range; and (e) can be studied separately for star-forming and quiescent objects, which means that the quenching channels that characterise the massive galaxies up to z = 1.5 can be explored, as well as the low-mass galaxies at low redshift.

The paper is organised as follows. In Sect. 2 we describe our photometric and spectroscopic dataset. The photometric redshifts and galaxy classification are presented in Sect. 3, the stellar mass estimates in Sect. 4. We detail the measurements of the galaxy SMFs and the associated uncertainties in Sect. 5, where we also point out the effect of the photometric absolute calibration in the new generation of large surveys. We present the evolution of the stellar mass function and stellar mass density in Sect. 6. Finally, we discuss our results and their effects on the quenching channels in Sect. 7.

Throughout this paper, we use the standard cosmology (Ω_m = 0.3, Ω_Λ = 0.7 with H₀ = 70 km s^-1 Mpc^-1). Magnitudes are given in the AB system (Oke 1974). The galaxy stellar masses are given in units of solar masses (M_⊙) for a Chabrier (2003) initial mass function (Chabrier IMF).

2. Data description

The observations and the data reduction are described in detail in the companion paper (Moutard et al. 2016) and are briefly summarised below.

2.1. Optical CFHTLS photometry

The CFHTLS¹ is an imaging survey performed with the MegaCam² camera in five optical bands, u,g,r,i, and z. It covers ~155 deg² over four independent fields with sub-arcsecond seeing (median ~0.8”) and reaches a 80% completeness limit of u ~ 24.4, g ~ 24.7, r ~ 24.0, i/y ~ 23.7, and z ~ 22.9 for extended sources in AB system. We emphasise that the y-band refers to the new i-band filter, in accordance with the CFHTLS notation. We have used the y-band response curve in our analysis when appropriate, but we refer to the “i” filter term regardless of whether it was observed with the i or y filter.

In this work we use the W1 (+ 02^h18^m00^s−07°00^m00^s) and W4 (+ 22^h13^m18s+ 01°19m00^s) fields. Two independent photometric catalogues have been released to the community: the 7th and final release (noted T0007³) of the CFHTLS produced by Terapix⁴, and the release from the CFHT Lensing Survey team (CFHTLenS⁵). Both catalogues are based on the same raw images. The AstrOmatic software suite⁶ has been used to generate the mosaic images (SWARP, Bertin et al. 2002) and to extract the photometric catalogues (SExtractor, Bertin & Arnouts 1996). The two releases differ in several points, however.

• In T0007, detection is based on gri−χ² images, while the galaxies in CFHTLenS are i-detected.

• A point spread function (PSF) homogenisation is implemented in CFHTLenS (see Hildebrandt et al. 2012) to improve the colour estimates. In practice, the PSF is homogenised across the field of view for each filter and degraded in all filters to the one with the highest PSF.

• A new photometric calibration has been applied to the T0007 release. While the previous releases and the CFHTLenS release rely on Landolt standard stars (see Erben et al. 2013), T0007 is based on the spectrophotometric tertiary standards from the Super Novae Legacy Survey (SNLS; see the procedure described in Regnault et al. 2009). In brief, each tile of the CFHTLS-Wide is re-observed (with short exposures) during stable photometric conditions and bracketed by two observations of the CFHTLS-Deep field containing the SNLS tertiary standards.

The difference in the calibration scheme of the two releases affects the final photometry. A comparison of the magnitudes for point-like sources between the T0007 and CFHTLenS releases reveals systematic offsets that are significantly larger than the expected uncertainties. These offsets are reported in Table 1 (column Δmag). We emphasise that the differences listed in this table are entirely due to the new calibration scheme established for T0007. The procedure used by the T0007 release allows transferring the percent level accuracy of the SNLS photometric calibration to the entire CFHTLS-Wide survey. For this reason, we use the T0007 catalogue as reference in this paper. However, we also perform the complete analysis with the CFHTLenS catalogue to discuss the effect of such differences in the photometric absolute calibration.

2.2. WIRCam K_s photometry

We conducted a NIR K_s-band follow-up of the VIPERS fields with the WIRCam instrument at CFHT (Puget et al. 2004). The layout of the observations is shown in Fig. 1 (red background). We covered a total area of ~27 deg² with an integration time per pixel of 1050 seconds. The image quality is very homogeneous, with an average seeing of all the individual exposures of ⟨ IQ ⟩ = 0.6″ ± 0.09. The data have been reduced by the Terapix team⁷ and the individual images were stacked and resampled on the pixel grid of the CFHTLS-T0007 release (Hudelot et al. 2012). The photometry was performed with SExtractor in dual-image mode with a gri−χ² image as the detection image and the same settings as those adopted for the T0007 release. The images reach a depth of K_s = 22 at ~3σ. The completeness reaches 80% at K_s = 22; this was determined from a comparison with the deepest surveys UKIDSS Ultra-Deep Survey(UDS, Lawrence et al. 2007) and VIDEO (Jarvis et al. 2013) in overlapping regions. Because the primary optical detection is based on the gri−χ² image, we may miss the reddest high-redshift galaxies. To account for this possible bias, we measured our source incompleteness as a function of magnitude, K, and colour, (z−K), with all the sources detected in the deep VIDEO survey. We derived a colour-magnitude weight map that we show in Fig. 2 and use in the remaining paper as multiplicative weighting for our statistical analyses. We refer to the companion paper, Moutard et al. (2016), for a complete description of the method that was used to build this weight map.

2.3. GALEX photometry

When available, we made use of the UV deep-imaging photometry from the GALEX satellite (Martin & GALEX Team 2005). We only considered the observations from the Deep Imaging Survey (DIS), which are shown in Fig. 1 as blue circles (?~ 1.1°). All the GALEX pointings were observed with the near-ultraviolet (NUV) channel with exposure times of T_exp ≥ 30 ks. Far-ultraviolet (FUV) observations are available for ten pointings in the central part of W1.

The large PSF of GALEX (FWHM ~ 5′′) means that source confusion becomes a great problem in the deep survey. To extract the UV photometry, we used a dedicated photometric code, EMphot (Conseil et al. 2011) that will be described in a separate paper (Vibert et al., in prep.). In brief, EMphot uses the stamps in the u-band (here the T0007 release) as priors, and they are then convolved by the GALEX PSF to produce a simulated image. The scaling factors to be applied to each u-band prior is obtained by simultaneously maximising the likelihood between the observed and the predicted fluxes for all the sources in tiles of a few square arc-minutes. The uncertainties on the UV fluxes account for the residuals in the simulated or observed image. A typical depth of NUV ~ 24.5 at ~5σ is observed over the entire survey. The NUV observations cover part of the WIRCam area with ~10.8 and 1.9 deg² in the W1 and W4 fields, respectively.

2.4. Final photometric catalogue

The catalogue of sources comes from the T0007 release and is based on detection in the gri−χ² image. As mentioned above, the same procedure was applied to the K_s images. Following Erben et al. (2013) and Hildebrandt et al. (2012), we used the T0007 isophotal apertures for the photometry to estimate the colours. The apertures are smaller than the Kron-like apertures (Kron 1980), which provides less noisy colours and leads to an improved photometric redshift accuracy (Hildebrandt et al. 2012). We also confirmed this with our large spectroscopic dataset (see below), which is especially relevant for faint sources (i′ > 23.5).

To derive galaxy physical properties, we need to know the total flux in all wavelengths. Therefore, we rescaled the isophotal flux to the Kron-like flux, ${\mathit{m}}_{\mathrm{total}}^{\mathit{f}}\mathrm{=}{\mathit{m}}_{\mathrm{ISO}}^{\mathit{f}}\mathrm{+}{\mathit{\delta }}_{\mathrm{m}}$, by adopting a unique factor, δ_m, for each source to preserve the colours. δ_m is the weighted mean of the individual scaling factors, ${\mathit{\delta }}_{\mathrm{m}}^{\mathit{f}}$, and is defined as ${\mathit{\delta }}_{\mathrm{m}}\mathrm{=}{\sum }_{\mathit{f}}{\mathit{\delta }}_{\mathrm{m}}^{\mathit{f}}{\mathit{w}}^{\mathit{f}}\mathit{/}{\sum }_{\mathit{f}}{\mathit{w}}^{\mathit{f}}$, with f = g,r,i,K_s, and w^f its associated errors⁸. Finally, the GALEX photometry, which corresponds to the total flux measurement (i.e. model PSF photometry) was added in the same way as to the optical and NIR magnitudes.

We here limit the catalogue to galaxies brighter than K_s < 22. The catalogue includes a total of ~1.3 millions sources over an area of ~27.1 deg², which drops to one million sources over ~22.4 deg² after applying the masks provided by the CFHTLenS team.

2.5. Spectroscopic sample

Our WIRCam survey has been designed to cover the VIMOS Public Extragalactic Survey (VIPERS; Guzzo et al. 2014) that is carried out with the VIMOS spectrograph and therefore provides many high-quality spectroscopic redshifts. We also added a compilation of the best-quality spectra from the VVDS survey (I_AB ≤ 24, Le Fèvre et al. 2013), the K < 23 limited UKIDSS spectroscopic Ultra Deep Survey (UDSz, Bradshaw et al. 2013; McLure et al. 2013), the low-resolution spectra (λ/ Δλ ~ 40) from the PRIsm MUlti-object Survey (PRIMUS, i ~ 23, Coil et al. 2011), and the bright-limited (i < 19.9) spectroscopic survey BOSS from the SDSS (Baryon Oscillation Spectroscopic Survey, Dawson et al. 2013). The K_s < 22 spectroscopic sample we used is presented in detail in the companion paper (Moutard et al. 2016).

We selected only the most secure spectroscopic redshifts, which means confidence levels above 95% for high-resolution surveys and σ< 0.005 (8% of outliers with δz/ (1 + z) > 5σ) for PRIMUS best redshifts. When they are available, the redshift measurements from VIPERS were used. Otherwise, the measurements from the deepest high-resolution spectra were favoured. In total, we assembled a K_s < 22-limited sample of 45,951 high-quality spectroscopic redshifts to calibrate and measure the accuracy of our photometric redshifts over the unmasked area of the survey (we refer to the companion paper for more details).

3. Photometric redshifts

3.1. Photometric redshift measurement

The photometric redshifts were computed with the spectral energy distribution (SED) fitting code Le Phare   (Arnouts et al. 2002; Ilbert et al. 2006), using the templates of Coupon et al. (2015). The new templates are based on the Ilbert et al. (2009) library of 31 empirical templates from Polletta et al. (2007), complemented by 12 star-forming templates from the Bruzual and Charlot stellar population synthesis models of 2003 (Bruzual & Charlot 2003, hereafter BC03). These templates were optimised to be more representative of the VIPERS spectroscopic sample (for more details we refer to Coupon et al. 2015).

The extinction was added as a free parameter with a reddening excess E(B−V) < 0.3 following different laws: Prevot et al. (1984), Calzetti et al. (2000), and a dusty Calzetti curve including a bump at 2175Å. No extinction was allowed for SEDs redder than Sc. The extinction law of Prevot et al. (1984) was used for templates redder than SB3 templates (see Ilbert et al. 2009) and the law of Calzetti et al. (2000) for bluer templates.

Finally, any possible difference between the photometry and the template library was corrected for by Le Phare  according to the method described in Ilbert et al. (2006). In brief, in each band the code tracks a systematic shift between the predicted magnitudes at known redshift and the observed magnitudes. Since our observation area is divided into 47 tiles of ≲1 deg² with the relative calibration varying from tile to tile, we performed a tile-by-tile colour optimisation. We used the median offset over all the tiles when there were not enough galaxies with spectroscopic redshift in the tile (${\mathit{N}}_{\mathrm{gal}}^{\mathrm{spec}}\mathrm{\le }\mathrm{100}$) available, which was the case in 12 tiles. We stress that the corrections were computed to better fit the colours and are therefore relative. We normalised the median offset on the K_s-band because the NIR fluxes are the same (see Sect. 2). The median relative offsets thus calculated for each photometric band of the T0007 and CFHTLenS catalogues can be found in Table 1, with the associated tile-to-tile deviation estimates (namely, the normalised median absolute deviation, NMAD). The difference between the T0007 and CFHTLenS relative offsets is consistent with the difference Δmag. In other words, we retrieved the shift between the two absolute photometric calibrations through the relative offsets computed with Le Phare. This safety check confirms that the colour optimisation achieved with Le Phare  absorbs the uncertainties that are linked to the photometric calibration.

3.2. Accuracy and precision of photometric redshifts

The comparison between our photometric redshifts and the corresponding spectroscopic redshifts for our K_s< 22 -limited sample is shown in Fig. 3. Using the NMAD to define the scatter⁹, we find σ_{Δz/ (1 + z)} ~ 0.05 for faint (i> 22.5) galaxies, while the scatter reaches σ_{Δz/ (1 + z)} ~ 0.03 for the bright (i< 22.5) galaxies. Our photo-z outlier rate¹⁰ is η = 1.2% and η = 9% for corresponding bright and faint samples, respectively (see Fig. 3, top panels, lower right corners).

Although the spectroscopic sample has been assembled to be as representative as possible, it is not as deep as the photometric sample. Aiming to correct this effect, we computed estimators that are weighted with respect to the i-band distribution of the photometric sample. By using these weighted estimators (marked with an orange w in Fig. 3, top panels), we obtained an accuracy ${\mathit{\sigma }}_{\mathrm{\Delta }\mathit{z}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}\mathrm{)}}^{\mathit{w}}\mathrm{~}\mathrm{0.03}$ for bright (i< 22.5) galaxies with an outlier rate of η^w = 1.4%, and ${\mathit{\sigma }}_{\mathrm{\Delta }\mathit{z}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}\mathrm{)}}^{\mathit{w}}\mathrm{~}\mathrm{0.07}$ and η^w = 16.4% for faint (i> 22.5) galaxies in our K_s< 22 -limited sample.

Even though the T0007 and CFHTLenS calibrations differ, the photometric redshifts obtained in both cases agree well¹¹ and their accuracies are similar. This is expected from the colour corrections described in Sect. 3.1, which absorb the differences between the two calibrations.

Finally, based on Fig. 3, we can define a range of reliable redshifts up to z = 1.5, with σ_{Δz/ (1 + z)}(z) < 0.1. The highest redshift bin that we consider, namely between z = 1.1 and z = 1.5, is characterised by the weighted accuracy ${\mathit{\sigma }}_{\mathrm{\Delta }\mathit{z}\mathit{/}\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}\mathrm{)}}^{\mathit{w}}\mathrm{~}\mathrm{0.08}$ and weighted outlier rate η^w ~ 20%.

3.3. Star and galaxy classification

Being able to separate galaxies and stars is crucial in our sample, especially for the W4 field, which is close to the Galactic plane and therefore highly populated by stars. Garilli et al. (2008) have found that more than 32% of the objects in the VVDS-Wide survey are stars. This is a pure i< 22.5-selected spectroscopic survey lying in the CFHTLS W4 field. Aiming to better control the type of the objects that we select as galaxies without compromising the completeness of our sample, we performed a classification based on three different diagnostics. Our classification is presented in detail in the companion paper (Moutard et al. 2016) and is summarised below.

• First, we used the maximum surfacebrightness versus magnitude (hereafterμ_max−m_obs) plane where bright point-like sources are well separated from galaxies (see Bardeau et al. 2005; Leauthaud et al. 2007).

• Secondly, we compared the reduced χ² obtained with galaxy templates described in Sect. 3.1 and a representative stellar library (based on Pickles 1998). An object can be defined as a star when its photometry is better fitted by a stellar spectrum.

• Finally, we used the g−z/z−K_s plane (equivalent to the BzK plane of Daddi et al. 2004) to isolate the stellar sequence and imposed that a star belong to this colour region. This sine qua non condition enabled us to catch faint stars while preventing us from losing faint compact galaxies.

We also identified a sample of QSOs (Type-1 AGNs) as point-like sources lying on the galaxy side of the BzK diagram. Dominated by their nucleus, the emission of these AGN galaxies is currently poorly linked to their stellar mass. However, they represent less than 0.5% of the objects, and we removed them from our sample without compromising its completeness.

All the objects that were not defined as stars or QSOs were considered as galaxies. We verified on a sample of 1241 spectroscopically confirmed stars that we caught 97% of them in this way, while we kept more than 99% of our spectroscopic galaxy sample. With this selection we finally found and removed ~8% and ~19% of objects at K_s < 22 for W1 and W4, respectively, outside the masked area.

4. Stellar mass estimation

4.1. Method

Stellar mass, M_∗, and the other physical parameters were computed by using the stellar population synthesis models of Bruzual & Charlot (2003) with Le Phare. As in Ilbert et al. (2013), the stellar mass corresponds to the median of the stellar mass probability distribution (PDF_M∗) marginalised over all other fitted parameters. Two metallicities were considered (Z = 0.008 and Z = 0.02 i.e. Z_⊙) and the star formation history declines exponentially following τ^-1e^{− t/τ} with nine possible τ values between 0.1 Gyr and 30 Gyr as in Ilbert et al. (2013).

The importance of the assumed extinction laws for the physical parameter estimation has been stressed in several recent studies, for example, by Ilbert et al. (2010) and Mitchell et al. (2013) for the stellar masses, or by Arnouts et al. (2013) for the star formation rate (SFR). We considered three laws and a maximum dust reddening of E(B−V) ≤ 0.5: the Prevot et al. (1984), the Calzetti et al. (2000), and an intermediate-extinction curve (see Arnouts et al. 2013, for more details). As in Fontana et al. (2006), we imposed a low extinction for low-SFR galaxies (E(B−V) ≤ 0.15 if age /τ> 4). The emission-line contribution was taken into account following an empirical relation between UV and line fluxes (Ilbert et al. 2009).

Using a method similar to Pozzetti et al. (2010), we based our estimate of the stellar mass completeness limit, M_lim, on the distribution of the lowest stellar mass, M_min, at which a galaxy could have been detected given its redshift. For our sample, which is limited at K_s < 22, M_min is given by $\log \mathrm{\left(}{\mathit{M}}_{\min }\mathrm{\right)}\mathrm{=}\log \mathrm{\left(}{\mathit{M}}_{\mathrm{\ast }}\mathrm{\right)}\mathrm{+}\mathrm{0.4}\mathrm{(}{\mathit{K}}_{\mathrm{s}}\mathrm{-}\mathrm{22}\mathrm{)}\mathit{.}$(1)We then considered the upper envelope of the M_min distribution. In each redshift bin, M_lim is defined by the stellar mass at which 90% of the population have M_∗>M_min. We show the resulting stellar mass completeness limits (open circles) as a function of redshift in Fig. 4 over the M_∗ distribution for our K_s < 22-limited sample.

4.2. Stellar mass error budget

In this section, we quantify the uncertainties associated with the stellar masses, which will be propagated into the error budget of the stellar mass functions. The first to be considered is the uncertainty in the flux measurements. The photon noise is taken into account by the Le Phare  code during the χ² SED fitting procedure (rescaling and model selection), where it returns the 68% confidence interval enclosed in the probability distribution function marginalised on the stellar mass (PDF_M∗).

The second source of error is introduced by the photometric redshift uncertainty, which is not included in the PDF_M∗. One way to measure its effect is to compare the stellar masses derived with the photometric and spectroscopic redshifts. We emphasise that this analysis is probably limited by the completeness of our spectroscopic sample. By contrast, it is powerful when used to reflect all the photo-z error contributions (quality of the photometry and representativity of the templates). The difference between the two mass estimates is shown in Fig. 5 as a function of stellar mass ${\mathit{M}}_{\mathrm{\ast }}^{{\mathit{z}}_{\mathrm{phot}}}$. In the top panel, we show the difference in four redshift bins between z = 0.2 and z = 1.5. No dependence with redshift is observed. The linear regression of the whole sample, plotted as a black dashed line, also suggests that it is not mass dependent. The bottom panel shows the ${\mathit{M}}_{\mathrm{\ast }}^{{\mathit{z}}_{\mathrm{phot}}}\mathit{/}{\mathit{M}}_{\mathrm{\ast }}^{{\mathit{z}}_{\mathrm{spec}}}$ dispersion in five stellar mass bins and reveals a median dispersion of 0.06 dex, with a maximum of 0.19 dex at low mass.

We then defined the resulting mass uncertainty as the sum in quadrature of all contributions: ${\mathit{\sigma }}_{\mathrm{M}}\mathrm{=}\sqrt{{\mathit{\sigma }}_{\mathrm{fit}}^{\mathrm{2}}\mathrm{+}{\mathit{\sigma }}_{\mathit{z}}^{\mathrm{2}}}\mathit{.}$(2)However, we have to keep in mind that the stellar mass estimation relies on the numerous assumptions made when we generate our SED templates. For example, Maraston (2005) has shown that a different treatment of the thermally pulsing asymptotic giant branch (TP-AGB) phase in the SSP can lead to a global shift in the stellar mass estimation¹². Ilbert et al. (2010) showed that the use of the Salpeter (1955) IMF instead of the Chabrier IMF decreases the stellar masses by 0.24 dex. These systematic shifts are therefore not expected to affect the conclusions of our study. Mitchell et al. (2013) also pointed out the potential effect of the assumed dust attenuation on the stellar mass estimation¹³. As presented in the previous section, we considered three different extinction laws. This allows a higher diversity of possible values for dust attenuation, which is expected to limit the bias that may affect our stellar mass estimation.

4.3. Effect of the CFHTLS absolute calibration

As shown in Sect. 2.1, the absolute photometric calibrations of the T0007 and CFHTLenS magnitudes differ by more than 0.05 mag on average. Even if the T0007 were significantly improved in its calibration, we compare the photometric redshifts and the stellar masses computed with both catalogues blindly to quantify the effect of these offsets.

As seen in Sect. 3.2, the colour corrections applied during the photometric redshift computation allows us to obtain very similar photo-z despite the offset between their calibrations. However, these corrections are a combination of the photometry and the SEDs used to calculate photo-z. As described in Sect. 3.1, the templates used for photo-z are different from those used for the masses. Consequently, we did not apply the photo-z colour corrections with the BC03 templates. The differences in the T0007 and CFHTLenS photometries thus directly affect the stellar mass estimation.

Figure 6 presents these differences in the redshift bins 0.2 < z < 0.5, 0.5 < z < 0.8, 0.8 < z < 1.1 and 1.1 < z < 1.5. The difference between the stellar masses obtained with the T0007 (${\mathit{M}}_{\mathrm{\ast }}^{\mathrm{T}\mathrm{07}}$) and those obtained with the CFHTLenS (${\mathit{M}}_{\mathrm{\ast }}^{\mathrm{LenS}}$) is stellar mass dependent. On average, this systematic difference can reach ± 0.1 dex at low redshift (z< 0.8). At higher redshift, we do not observe a systematic difference between the two stellar mass catalogues since we used the same K_s-band calibration. Even if the object-by-object ${\mathit{M}}_{\mathrm{\ast }}^{\mathrm{LenS}}$ to ${\mathit{M}}_{\mathrm{\ast }}^{\mathrm{T}\mathrm{07}}$ ratio is characterised by a mean offset that never exceeds 0.2 dex, the comparison of the T0007 and CFHTLenS number counts above the mass completeness limit reveals different shapes around M_∗ ~ 10¹¹ M_⊙, notably at z< 0.8. This suggests a significant effect on the SMF massive end at low redshift, as we discussed below.

5. Measuring the stellar mass functions

To compute the galaxy stellar mass function, we selected a sample of ~760 000 galaxies at K_s ≤ 22 over an effective area of 22.38 deg². According to what we discussed in Sect. 3.2, we restricted our analysis to the range 0.2 ≤ z ≤ 1.5, where we combined reliable redshifts and large volumes.

The galaxy stellar mass function was derived with the tool ALF (Ilbert et al. 2005), which provides three non-parametric estimators: V_max (Schmidt 1968), SWML (the step-wise maximum likelihood; Efstathiou et al. 1988), and C⁺ (Zucca et al. 1997). The V_max estimator is most widely used because of its simplicity. The 1 /V_max is the inverse sum of the volume in which each galaxy was observed. The V_max is the only estimator that is directly normalised. The SWML (Efstathiou et al. 1988) determines the SMF by maximising the likelihood of observing a given stellar mass – redshift sample. The C⁺ method overcomes the assumption of a uniform galaxy distribution, as is the case when using the V_max¹⁴. As described in Ilbert et al. (2015), these estimators diverge below a stellar mass limit that should correspond to the limit calculated in Sect. 4. In Fig. 4 we verify that the V_max and SWML estimators (black dots) are consistent with our K_s-based stellar mass completeness limit (black open circles). We used the colour-magnitude weight map shown in Fig. 2 to correct the SMF for the potential incompleteness described in Sect. 2.2. In the remainder of this study, we work with stellar masses M_∗>M_min where all the non-parametric estimators agree.

5.1. Measurements by type and field

To separate quiescent and star-forming galaxies, we used the rest frame (NUV−r)° versus (r−K)° diagram (hereafter NUVrK) presented by Arnouts et al. (2013), which is based on the method introduced by Williams et al. (2009). Figure 8 presents the galaxy distribution in the NUVrK diagram for several redshift bins. This optical-NIR diagram allows us to properly separate red dusty star-forming galaxies from red quiescent ones. Edge-on spirals are clearly identifiable, as is illustrated by the morphological study of the NUVrK diagram at low redshift presented in the companion paper (Moutard et al. 2016).

When we computed the rest-frame colours, we adopted the procedure described in Appendix A.1 of Ilbert et al. (2005) to minimise the dependency of the absolute magnitudes to the template library. An absolute magnitude at λ⁰ was derived from the apparent magnitude in the filter passband that was the closest from λ⁰× (1 + z) to minimise the k-correction term, except when the apparent magnitude had an error above 0.3 mag, to avoid too noisy colours. The small break¹⁵ in the red clump is artificial and is an effect of the template discretisation, when our procedure used to limit the template dependency fails because of the low signal-to-noise ratio measurements (here due to the intrinsic low rest-frame NUV emission of quiescent galaxies)¹⁶.

As shown in Fig. 8, by following the low-density valley of the NUVrK diagram (the so-called green valley), the selection of quiescent galaxies can be defined with the general form $\mathrm{[}\mathrm{(}\mathrm{NUV}\mathrm{-}\mathit{r}{\mathrm{)}}^{\mathrm{◦}}\mathit{>}{\mathit{B}}_{\mathrm{2}}\mathrm{]}\mathrm{\cap }\mathrm{[}\mathrm{(}\mathrm{NUV}\mathrm{-}\mathit{r}{\mathrm{)}}^{\mathrm{◦}}\mathit{>}\mathit{A}\mathrm{(}\mathit{r}\mathrm{-}{\mathit{K}}_{\mathrm{s}}{\mathrm{)}}^{\mathrm{◦}}\mathrm{+}{\mathit{B}}_{\mathrm{1}}\mathrm{]}\mathit{.}$(3)A, B₁, and B₂ are three parameters to be adjusted in each redshift bin, as suggested by Ilbert et al. (2015) and Mortlock et al. (2015), because of the global ageing of the galaxy population.

In the four redshift bins, the slope A of Eq. (3) seems to be constant, with a typical value of A = 2.25. By projecting the galaxy distribution in a plane perpendicular to the axis of slope A¹⁷, we clearly distinguish the red and blue clouds as two normal distributions that we fitted by two Gaussians. We define B₁ as the position where the two Gaussians intersect.

In Fig. 7 we show the evolution of B₁ as a function of the look-back time (t_l). By assuming a linear relation between B₁ and cosmic time, we derive B₁(t_l) = −0.029t_l + 2.368 in our highest precision redshift range (0.2−1.0). Assuming that B₂ evolves as B₁, we empirically set B₂(t_l = 2.5 Gyr) = 3.3, and we find B₂(t_l) = B₁(t_l) + 1.004. We can write our selection of quiescent galaxies as $\begin{array}{ccc}& & \left[\mathrm{(}\mathrm{NUV}\mathrm{-}\mathit{r}{\mathrm{)}}^{\mathrm{◦}}\mathit{>}\mathrm{3.372}\mathrm{-}\mathrm{0.029}{\mathit{t}}_{{l}}\right]\\ & & \end{array}$(4)All the galaxies that are not selected as quiescent are considered to be star forming. In Fig. 8 the separations between quiescent and star-forming galaxies are shown as white solid line. We also define the green valley as the region around minimum B₁, reaching 10% of the peak of the red Gaussians, as shown by the white dotted lines. We consider these limits as possible systematic uncertainties when discussing the evolution of the quiescent and star-forming SMFs.

Figure 9 presents the global (black), star-forming (blue), and quiescent (red) galaxy SMFs for the two fields separately (W1: dot and W4: cross) in four redshift bins. The sample consists of 481 518 galaxies over 14.43 deg² in W1 and 268 010 galaxies over 7.96 deg² in W4. The error bars shown in the upper sub-panels reflect only the Poissonian contributions. The SMF comparison between the two fields agrees within the errors. In the lower sub-panels, we plot the stellar mass uncertainty by type, σ_M, defined in Sect. 4.2, as function of the stellar mass. First, σ_M decreases exponentially with stellar mass, as already noted in previous studies (e.g. Grazian et al. 2015). We can then fit the σ_M(M_∗) relation with a power law (Fig. 9 sub-panels, dashed lines). Secondly, the size of our galaxy sample allows for very small relative Poissonian errors down to densities of around ~10^-5−10^-6 Mpc^-3 even if we split by type and field. The cosmic variance contribution in the budget of the errors that affects our SMF measurement is therefore expected to be small, as discussed in the next section.

5.2. SMF uncertainties

In this section, we describe the error budget associated to our SMF measurements. All the contributions to the SMF uncertainties are expressed as a function of the stellar mass and redshift. In addition to the stellar mass and Poissonian errors already mentioned, the large-scale density inhomogeneities represent a source of uncertainty. This cosmic variance is known to represent a fractional error of 15−25% at the high-mass end (M ≥ 10¹¹ M_⊙) in the COSMOS survey and of around 20−50% in narrower pencil-beam surveys, generally dominating the error budget.

Following the procedure discussed by Coupon et al. (2015), we investigated the contribution of the cosmic variance in our sample by dividing our survey into N patches of equal areas. Since the effective surface can change from one patch to another, every patch was weighted according to its unmasked area. For a given observed area, we computed the number density dispersion N times over (N−1) patches by discarding a different patch every time. We then considered the mean number density dispersion over the N measurements as our internal estimate of the cosmic variance for a given effective area and the dispersion around the mean as an error estimate of the cosmic variance.

In Fig. 10 we plot our cosmic variance estimate σ_cv in the redshift bins [0.2, 0.5] and [1.1, 1.5], considering three stellar mass bins from M_∗ = 10¹⁰ M_⊙ up to M_∗ = 10^11.5 M_⊙ (with blue, purple and red dots, respectively) and for mean effective areas ranging from a ≃ 0.1 to a ≃ 2.8 deg². The relation of cosmic variance – area is well fitted by a power-law with σ_cv(a) = 10^βa^α (shown as dashed lines). To estimate the cosmic variance that affects our entire survey, we extrapolated the relations up to a = 22 deg², shown as squares.

For comparison, we also show the cosmic variance predicted for the same redshift and stellar mass bins (triangles) by using the code getcv (Moster et al. 2011). Our internal cosmic variance estimate (σ_cv) and the predicted one agree remarkably well up to our observed areas of a = 2 deg². For larger areas, the two estimates diverge slightly for high-mass (M_∗ > 10¹¹ M_⊙) galaxies at z < 0.5, where we slightly underestimate σ_cv with respect to the theoretical prediction. We have to stress that the Moster et al. (2011) procedure is optimised for pencil-beam surveys of areas a < 1 deg². At z > 0.5, the theoretical estimators always predict a cosmic variance lower than our own extrapolation. By using our internal estimate, we therefore adopt a conservative approach.

Finally, the last source of error that we need to consider is that of the stellar mass uncertainty defined in Sect. 4.2. To do so, we generated 200 mock catalogues with perturbed stellar masses according to the expected σ_M (which includes the photometric redshift uncertainties and the photon noise, Eq. (2)) and measured the 1σ dispersion in the density Φ of the reconstructed SMFs that we refer to as σ_Φ,M.

At the end, the error of the stellar mass function is the quadratic sum of all the contributions discussed above and is defined as ${\mathit{\sigma }}_{\mathrm{tot}}\mathrm{=}\sqrt{{\mathit{\sigma }}_{\mathrm{cv}}^{\mathrm{2}}\mathrm{+}{\mathit{\sigma }}_{\mathrm{poi}}^{\mathrm{2}}\mathrm{+}{\mathit{\sigma }}_{\mathrm{\Phi }\mathit{,}\mathrm{M}}^{\mathrm{2}}}\mathit{.}$(5)

5.3. Importance of photometric calibration in large surveys

As mentioned in Sect. 4.3, a mean offset of ~0.06 mag in the optical absolute photometric calibration (cf. Δmag in Table 1) can change the stellar mass estimate by 0.1 dex. In the top panel of Fig. 11 we show the difference between the two SMFs measured with the T07 and CFHTLens photometries, ΔΦ_calib = Φ_LenS−Φ_T07. This difference is normalised by the total statistical error discussed in the previous section, σ_tot. In general, ΔΦ_calib>σ_tot in our survey (solid lines), which means that the SMF variation induced by the calibration offsets is several times larger than the uncertainty of our SMF. It even reaches 5 σ_tot at low redshift (0.2 < z < 0.5; blue solid line), where the stellar mass is essentially driven by the optical photometry.

In contrast, by considering a subsample of 2 deg² (dashed lines), we find that $\mathrm{\left|}\mathrm{\Delta }{\mathrm{\Phi }}_{{calib}}\mathrm{\right|}\lesssim {\mathit{\sigma }}_{\mathrm{tot}}^{\mathrm{2}\mathrm{deg}}$ (green shaded area). This means that the SMF variation driven by the calibration offsets is smaller than the other uncertainties affecting the SMF in a 2 deg² survey (i.e. the variation is contained within the error bars). In other words, we cannot see the variation that is due to the calibration because it is hidden by other sources of uncertainties (Poissonian and cosmic variance). The systematic differences that are due to the T0007-CFHTLenS photometric offsets can therefore be neglected in a 2 deg² survey, while in surveys of 20 deg² and more, we reach a regime where the systematic uncertainty that is due to the photometric calibration dominates the error budget.

For comparison, we also investigated another source of systematic uncertainty: the photometric redshift bias. Using the photo-z bias (z_bias = z_phot−z_spec) presented in Sect. 3.1 (see Fig. 3, lower panels), we corrected our photometric redshifts. Instead of using a global correction, we applied a photo-z bias correction for different galaxy types¹⁸.

Similarly to ΔΦ_calib, the difference between the stellar mass functions computed with the corrected photometric redshifts and with the original ones (ΔΦ_zbias) is shown in the lower panel of Fig. 11. The effect of the photo-z bias on the SMF measurement is much weaker than the effect of the photometry. The SMF differences induced by the photo-z bias as measured in our sample are largely dominated by the statistical uncertainties in a 2 deg² subsample. Moreover, in the entire 22 deg² survey, the difference can only be detected at z ~ 0.65, while | ΔΦ_zbias | < 2 σ_tot.

Given the limited amplitude of its effect on our sample, the photo-z bias can be neglected in our study. By contrast, the SMF variations that are due to the difference in photometry stress the need of carefully controlling the absolute photometric calibration in large surveys. In the present study, the choice of the CFHTLS-T0007 photometry is supported by (1) the SNLS photometric calibration based on a new spectrophotometric standard for high-precision cosmology; and (2) the careful treatment by the Terapix team that enables homogeneously propagating the SNLS photometric calibration over the entire survey.

6. Evolution of the galaxy stellar mass function and density

As shown in Sect. 5.2, the large volume probed by our survey allows us to reduce both the cosmic variance and the Poisson uncertainties. We exploit this large volume to quantify the evolution of the galaxy SMF, especially at the high-mass end, where it is most relevant.

6.1. Evolution of the SMF

6.1.1. Comparison of the global SMF with the literature

In the upper sub-panels of the Fig. 12, we compare our global SMF measurements with the literature. Our results agree well overall with many previous studies, although some differences exist. We discuss these in this section. The error bars corresponding to our measurement (black) reflect the total error σ_tot defined in Eq. (5). In the lower sub-panels, we show the contribution of each error to the error budget. We show each contribution normalised by the total error σ_tot as a function of the stellar mass. First, we note that the Poissonian error (blue dash-dotted line) represents a minor contribution to the total SMF uncertainty up to the very high-mass end (i.e. above 10^11.5 M_⊙). Secondly, the contribution of the cosmic variance (σ_cv; red dashed line) is dominant up to stellar masses around the SMF knee (M_∗ < 10¹¹ M_⊙). We finally note that the contribution of the stellar mass uncertainty (σ_Φ,M; cyan solid line) drives the total uncertainty of the SMF high-mass end. We recall that while the Poissonian uncertainty is always taken into account in the literature, the error bars may reflect different contributions to the SMF uncertainty depending on the study considered in Fig. 12 (as specified in the caption).

The comparison with Davidzon et al. (2013) is straightforward since their observations were taken in the same two fields of the CFHTLS survey, covering an effective area of 5.34 and 4.97 deg² in W1 and W4, respectively. The authors derived the SMF between z = 0.5 and z = 1.3 using the VIPERS-PDR1 dataset (~50 000 galaxies), that is, the main spectroscopic sample used to calibrate our redshifts (Sect. 3.2). The work of Davidzon et al. (2013) clearly illustrates the advantages of using spectroscopic redshifts (e.g. the easier removal of stellar interlopers and QSO). However, to estimate the SMF through spectroscopic data, some difficulties need to be solved, such as the statistical weighting to account for the spectroscopic sampling rate (see Garilli et al. 2014, for more details about how these weights are computed in VIPERS). We observe a good agreement between the two SMF estimates, especially at M_∗> 10¹¹ M_⊙.

The statistical uncertainties are very low in both VIPERS and our analysis, and the two surveys are additionally collected almost in the same area. Any difference is likely due to some combination of the photometric redshift uncertainty of our sample, the spectroscopic incompleteness affecting VIPERS, the adopted SED fitting method, or the photometric calibration used in VIPERS (T0005) and in our survey (T0007). However, the only significant discrepancy is observed close to the stellar mass completeness limit of VIPERS, where the measurements of Davidzon et al. (2013) are ≲0.2 dex lower. This difference at low masses could be due to some incompleteness correction that is due to the i-band selection, while our sample is K_s-band selected.

Moustakas et al. (2013) also measured the SMF by relying on the spectroscopic redshift sample of PRIMUS (~40 000 galaxies between z = 0.2 and 1 and cover ~5.5 deg² over five fields). In general, we observe that the SMF measurements from PRIMUS form the upper limit of the literature. Their SMF estimate is significantly above the others at 0.5 < z < 0.8. In the range 10^10.5<M_∗ < 10¹¹ M_⊙, the difference reaches 0.2 dex, while the authors predict that the cosmic variance should not affect the measurement by more than 10%; a larger offset is observed at M_∗ > 10^11.5 M_⊙, which could be mainly explained by the cosmic variance affecting their measurement, which is estimated to be very strong at high mass¹⁹. In the next redshift bin (0.8 < z < 1.1), the SMF of Moustakas et al. (2013) is also significantly higher than ours. The reason for the discrepancy may be linked to the different recipe (dust models, template libraries, etc.) adopted by Moustakas et al. (2013) in their SED fitting procedure (see also Davidzon et al. 2013, for a discussion about the effect of different SED fitting methods on the SMF). We compared their stellar masses in the XMM-LSS field, which overlaps the W1 field. We found that the PRIMUS masses are higher than ours by 0.17 ± 0.09 dex at 0.2 < z < 0.5, 0.15 ± 0.08 dex at 0.5 < z < 0.8 and 0.12 ± 0.1 dex at 0.8 < z < 1²⁰. This could explain part of the observed shift in the SMFs. It is worth noting that the two largest spectroscopic surveys so far, VIPERS and PRIMUS, lead to the largest difference in the SMF measurements. This highlights the great effect of systematic uncertainties in the latest large surveys (see also Coupon et al. 2015).

The comparison of our measurements with deep photometric surveys shows that our results agree well with those of Tomczak et al. (2014) and Ilbert et al. (2013), down to the lowest stellar masses we can explore. Their analysis was based on much deeper data, which confirms the estimate of our lower mass limits. Only in the redshift bin 0.8 < z < 1.1, we note a significant difference with Ilbert et al. (2013) in the high-mass end of the SMF. This can be explained by the well-known over-density in the COSMOS field (Kovač et al. 2010; Bielby et al. 2012)²¹.

Finally, we show the local, z ~ 0.1, GALEX-SDSS SMF from Moustakas et al. (2013) in all panels of Fig. 12 (as a dash-dotted line). A small but clear and progressive deviation of the SMF with redshift is obvious, in comparison to the local SMF. Far from evident in the previous studies, the trend observed at high mass is confirmed and quantified in Sects. 6.1.3 and 6.2, and is discussed in Sect. 7.1.

6.1.2. Fitting the global, star-forming, and quiescent SMF

To quantify the evolution of the SMF, we adopted the parametrisation proposed by Schechter (1976). As already noted, the total stellar mass function is better fitted with a double Schechter function (Baldry et al. 2008; Pozzetti et al. 2010; Ilbert et al. 2013; Tomczak et al. 2014), defined as $\mathrm{\Phi }\mathrm{\left(}{\mathit{M}}_{\mathrm{\ast }}\mathrm{\right)}\mathrm{d}{\mathit{M}}_{\mathrm{\ast }}\mathrm{=}{\mathrm{e}}^{\mathrm{-}\frac{{\mathit{M}}_{\mathrm{\ast }}}{{\mathrm{ℳ}}^{\mathit{\star }}}}\left[{\mathrm{\Phi }}_{\mathrm{1}}^{\mathit{\star }}{\left(\frac{{\mathit{M}}_{\mathrm{\ast }}}{{\mathrm{ℳ}}^{\mathit{\star }}}\right)}^{{\mathit{\alpha }}_{\mathrm{1}}}\mathrm{+}{\mathrm{\Phi }}_{\mathrm{2}}^{\mathit{\star }}{\left(\frac{{\mathit{M}}_{\mathrm{\ast }}}{{\mathrm{ℳ}}^{\mathit{\star }}}\right)}^{{\mathit{\alpha }}_{\mathrm{2}}}\right]\frac{\mathrm{d}{\mathit{M}}_{\mathrm{\ast }}}{{\mathrm{ℳ}}^{\mathit{\star }}}\mathit{,}$(6)where ℳ^⋆ is the characteristic stellar mass, ${\mathrm{\Phi }}_{\mathrm{1}}^{\mathit{\star }}$ and ${\mathrm{\Phi }}_{\mathrm{2}}^{\mathit{\star }}$ are the normalisation factors, and α₁ and α₂ are the power-law slopes satisfying α₂<α₁.

It has been shown that the massive end of the stellar mass function can be significantly affected by the stellar mass uncertainty (Caputi et al. 2011) through the so-called Eddington bias (Eddington 1913). We corrected the SMF for the Eddington bias during the fitting process by convolving the SMF parametric form by the stellar mass uncertainty σ_M²² following the procedure described in Ilbert et al. (2013). The authors estimate dσ_M for each redshift bin, but Grazian et al. (2015) have pointed out the importance of using an estimate of σ_M that depends on the stellar mass in addition to the redshift dependence²³. We used the σ_M(M_∗,z) estimate described in Sect. 5.1 (cf. Fig. 9, sub-panels).

Figure 13 shows the SMF of the global (black stars), star-forming (blue stars), and quiescent (red stars) populations. We included the SMFs measured by Tomczak et al. (2014) and Ilbert et al. (2013), who probed the very low-mass populations for SF and Q galaxies. A simple Schechter function (i.e. ${\mathrm{\Phi }}_{\mathrm{2}}^{\mathit{\star }}\mathrm{=}\mathrm{0}$ in Eq. (6)) seems to be sufficient to fit the star-forming contribution above the stellar mass completeness limit (blue dashed lines). However, as already shown by several studies working with deeper surveys (see e.g. Drory et al. 2009; Ilbert et al. 2013; Tomczak et al. 2014), the SMF of star-forming galaxies reveals an upturn at low mass and is better fitted with a double-Schechter function (Tomczak et al. 2014, Sect. 3.2). Given our stellar mass completeness limit, we can only constrain the low-mass end of the star-forming SMF at z< 0.5. Therefore we set the low-mass components of the double-Schechter function to the values found at 0.2 <z< 0.5, α_{2 sf} = −1.49 and ${\mathrm{\Phi }}_{\mathrm{2}{sf}}^{\mathit{\star }}\mathrm{=}\mathrm{-}\mathrm{3.24}$. Our choice is supported the lack of evolution that is observed for α₂ and ${\mathrm{\Phi }}_{\mathrm{2}}^{\mathit{\star }}$ by Ilbert et al. (2013) and Tomczak et al. (2014). In addition, our values agree quite well with Tomczak et al. (2014), who probed the SMF at lower stellar mass. The resulting double-Schechter function is plotted in Fig. 13 (blue solid line).

For the quiescent galaxies, we clearly need a double-Schechter function to fit the SMF at low redshift (red stars in Fig. 13 upper left panel). The upturn at low mass is slightly more pronounced in our measurement than in the literature, regardless of the position of the quiescent galaxy selection in the NUVrK green valley (cf. Sect. 5.1). In other words, the low-mass slope that we measure does not depend significantly on our selection of quiescent galaxies. We also verified that we find the same shape when we select the quiescent galaxies based on their specific star formation rate (sSFR), using sSFR< 10^-11 Gyr^-1 (see Ilbert et al. 2010, for more details on this threshold). We find the upturn position around M_∗ ≃ 10^9.5 M_⊙, in good agreement with previous measurements, that is, M_∗ ≃ 10^9.2 M_⊙, M_∗ ≃ 10^9.4 M_⊙ and M_∗ ≃ 10^9.6 M_⊙ for Tomczak et al. (2014), Ilbert et al. (2013) and Drory et al. (2009), respectively. Even though several deep surveys show that the low-mass upturn of the quiescent SMF is still present at z> 0.5, using a single-Schechter function is sufficient given our survey stellar mass limit. The discrepancies between our star-forming and quiescent SMF and the literature are mainly explained by the different criteria used to separate quiescent and star-forming galaxies. If we include the galaxies lying in the green valley in the quiescent sample (i.e. by considering the upper or lower envelopes of the quiescent or star-forming SMF), our measurements of the SMF agree with those of Tomczak et al. (2014) and Ilbert et al. (2013) at z< 1.1. At higher redshift, including the green valley in the quiescent locus of the NUVrK diagram is not enough to reconcile the estimates.

We cannot exclude the possibility that we may have missed some fainter red galaxies as a result of the gri-detection described in Sect. 2.4. However, this effect should be limited since we corrected for this incompleteness according to the weight colour map shown in Fig. 2, as previously explained. To add an independent validation of our procedure of correcting for this incompleteness, we used the CFHTLS-Deep/WIRcam Deep Survey (WIRDS; Bielby et al. 2012), which overlaps our CFHTLS-Wide/K_s< 22 survey. We estimated the completeness of the gri selection as a function of redshift and stellar masses separately for quiescent and star-forming galaxies. Below z< 1.1, we did not find any completeness problems, regardless of galaxy type or stellar mass range. At 1.1 <z< 1.5, the quiescent sample is >85% complete after applying our weighting scheme, and applying a correction based on the WIRDS sample would shift the density by less than 0.1 dex, which is well inside our uncertainties. It appears that star-forming galaxies can also be affected by incompleteness around M_∗ ~ 10^10.9 M_⊙ (probably because of dust extinction in massive galaxies at high redshift)²⁴. However, comparison with the literature suggests that our SF sample does not significantly suffer from this incompleteness (Fig. 13).

Moreover, we have to highlight that our total SMF agrees with Tomczak et al. (2014) at M_∗>M_lim, while our SMF estimate for SF galaxies is continuously higher (by 0.02 dex at M_∗ = 10^10.5 M_⊙ and 0.08 dex at M_∗ = 10^10.75 M_⊙). This SMF difference for SF galaxies would allow a transfer (between the SF and Q populations) that is sufficient to reconcile the our SMF estimate for quiescent galaxies with the estimate of these authors. This stresses the sensitivity of the SMF to the Q/SF selection²⁵.

Since the low-mass end of the global SMF is strongly dominated by the star-forming population at z> 0.5, we assumed the same parametrisation of ${\mathit{\alpha }}_{\mathrm{2}}^{\mathrm{\ast }}$ and ${\mathrm{\Phi }}_{\mathrm{2}}^{\mathit{\star }}$ (i.e. ${\mathit{\alpha }}_{\mathrm{2}}^{\mathrm{\ast }}\mathrm{=}{\mathit{\alpha }}_{\mathrm{2}{sf}}^{\mathrm{\ast }}$ and ${\mathrm{\Phi }}_{\mathrm{2}}^{\mathit{\star }}\mathrm{=}{\mathrm{\Phi }}_{\mathrm{2}{sf}}^{\mathit{\star }}$). We derived two parametric forms of the global SMF, depending on whether the double or the simple Schechter form of the star-forming SMF is considered, as shown in Fig. 13 (with dashed and solid black lines, respectively). The corresponding best-fit parameters are reported in Tables 2 and 3, respectively²⁶.

6.1.3. Quantifying the SMF evolution

In Fig. 14 we plot the evolution of the SMF for all (left panel), star-forming (middle panel), and quiescent galaxies (right panel). Each redshift bin is coded with a different colour. As in Fig. 13, the shaded areas show the systematic uncertainty induced by the star-forming or quiescent classification in the NUVrK diagram, while the solid lines represent the parametric form of reference. The arrows show the position of the corresponding characteristic mass ℳ^⋆.

As mentioned above, the galaxy population at low masses is strongly dominated by its star-forming component, and the global SMF evolution is then mainly driven by the star-forming population. We note that the evolution of the global SMF that is characterised by a ~0.2 dex increase of the ℳ^⋆ (see the arrows in Fig. 14 left panel). However, there is almost no evolution of the star-forming population (Fig. 14 middle panel): the characteristic mass is nearly constant, with $\log \mathrm{(}{\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}\mathit{/}{\mathit{M}}_{\mathrm{\odot }}\mathrm{)}\mathrm{=}\mathrm{10.6}{\mathrm{6}}_{-0.03}^{\mathrm{+}\mathrm{0.02}}$ in the redshift range 0.2 <z< 1.5, while the evolution of the low-mass slope remains very stable, as discussed previously. This confirms that the probability of finding a star-forming galaxy declines exponentially above a certain stellar mass ${\mathit{M}}_{\mathrm{\ast }}\mathit{>}{\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$, which is constant with time. This stresses that the star formation seems to be impeded beyond this stellar mass independent of the redshift up to z = 1.5. This is one of the cornerstones of the empirical description proposed by Peng et al. (2010), in which the evolution of high-mass galaxy is dominated by internal quenching mechanisms (named mass quenching by the authors). Peng et al. (2010) suggested that the efficiency of mass quenching is proportional to SFR/ℳ^⋆ to keep the SMF of star-forming galaxies constant with redshift.

The right panel of Fig. 14 shows that the main contribution to the evolution of the total SMF is due to the quiescent population build-up. In addition to galaxies that are quenched by mass quenching (around ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$), the SMF evolution of quiescent galaxies reveals an increase of low-mass galaxies with time, as shown in Ilbert et al. (2010). In particular, the SMF upturn built at z< 0.5 suggests that the star formation of M_∗< 10^9−9.5 M_⊙ galaxies is efficiently quenched, at least at low redshift. Ascribed by Peng et al. (2010) to environmental quenching, the build-up of the low-mass quiescent population is discussed in Sect. 7.2. The increase of the very high-mass population that we observe in the quiescent sample (and consequently also in the total SMF) is discussed in Sect. 7.1.

6.2. Evolution of the number densities and stellar mass densities

We derived the galaxy number and stellar mass densities, n_∗ and ρ_∗, respectively, by integrating the stellar mass function ${\mathit{n}}_{\mathrm{\ast }}\mathrm{=}{\mathrm{\int }}_{{\mathit{M}}_{\mathrm{1}}}^{{\mathit{M}}_{\mathrm{2}}}\mathrm{\Phi }\mathrm{\left(}{\mathit{M}}_{\mathrm{\ast }}\mathrm{\right)}\mathrm{d}{\mathit{M}}_{\mathrm{\ast }}$(7)and ${\mathit{\rho }}_{\mathrm{\ast }}\mathrm{=}{\mathrm{\int }}_{{\mathit{M}}_{\mathrm{1}}}^{{\mathit{M}}_{\mathrm{2}}}\mathrm{\Phi }\mathrm{\left(}{\mathit{M}}_{\mathrm{\ast }}\mathrm{\right)}{\mathit{M}}_{\mathrm{\ast }}\mathrm{d}{\mathit{M}}_{\mathrm{\ast }}\mathit{.}$(8)We adopted the parametric form of the SMF corrected for the Eddington bias. We derived the number densities above the stellar mass completeness limit. The stellar mass density was calculated by integrating the SMF over the stellar mass range 9 < log (M_∗/M_⊙) < 13, as in Tomczak et al. (2014). We recall that at z> 0.5, the stellar mass density relies partially on the extrapolation of the SMF to the lower stellar mass limit.

In Fig. 15 we plot the cosmic evolution of the number densities, n_∗, in the stellar mass bins 10.5 < log (M_∗/M_⊙) < 11 (left), 11 < log (M_∗/M_⊙) < 11.5 (middle), and 11.5 < log (M_∗/M_⊙) < 12 (right), between redshifts z = 0.2 and z = 1.5. For every mass bin, we show the densities for the global, star-forming, and quiescent galaxy populations that we compare with the measurements from Moustakas et al. (2013, triangles) and Matsuoka & Kawara (2010, pentagons) when available. For the global population in our sample, we distinguish two types of evolution. In the two lowest stellar mass bins (10^10.5<M_∗/M_⊙< 10^11.5), we observe a two-phase evolution, with an increase of ~25−50% from z ~ 1.3 down to z ~ 1, followed by a plateau down to z ~ 0.2. For the most massive population (M_∗> 10^11.5 M_⊙), we observe a continuous increase by slightly less than a factor two from z ~ 1.5 to z ~ 0.2. A similar, but weaker, trend is seen in VIPERS because of the narrower redshift range. Our results are directly comparable with Matsuoka & Kawara (2010) for M_∗> 10¹¹ M_⊙. These authors also took the Eddington bias in their density estimates into account (the estimates are based on simulations). They also emphasised that their measurements at z< 0.5 are strongly biased because of their less reliable photo-zs. Within these limits, our n_∗ evolution measurements for the entire population agree well with their results. The trend observed with PRIMUS is also similar for the lowest mass bins, M_∗< 10^11.5 M_⊙, although they have systematically higher densities (~40% and ~25% for M_∗ ~ 10^10.75 M_⊙ and M_∗ ~ 10^11.25 M_⊙, respectively), as expected from the higher normalisation of their SMFs (cf. Fig. 12). In addition, it is important to recall that they did not take Eddington bias into account, which can enhance the differences, especially at M_∗> 10¹¹ M_⊙.

For the evolution by galaxy type, we observe a two-phase evolution of M_∗< 10^11.5 M_⊙ quiescent galaxies, while star-forming galaxies experience a constant evolution, if not a decreasing evolution. At low mass, M_∗< 10¹¹ M_⊙ (left panel), the density of quiescent galaxies increases with redshift and equals the star-forming density in the lowest redshift bin, at z ~ 0.3. For the intermediate masses, 10¹¹<M_∗/M_⊙< 10^11.5 (middle panel), the quiescent population becomes dominant at higher redshift, z ~ 0.9. In the highest stellar mass bin (M_∗> 10^11.5 M_⊙, right panel), the quiescent population always outnumbers the star-forming one by representing already 50−60% of the global population at z ~ 1.3 and more than 80% at z ~ 0.3 (i.e. n_∗ multiplied by 2.5). From z ~ 1 to z ~ 0.2, the number density of the massive star-forming galaxies has diminished by a factor of 1.5 and 2 in the two highest mass bins, respectively.

The number densities computed in VIPERS are not plotted since the stellar mass bins used by Davidzon et al. (2013) are different from ours. However, the authors observed the same general trends, though their uncertainties prevent them from distinguishing the two-phase evolutions observed in our survey (Davidzon et al. 2013, Fig. 6). We also generally agree with Moustakas et al. (2013) for star-forming galaxies, as our studies observe a decreasing n_∗ between z = 1 and z = 0.3 for M_∗> 10^10.5 M_⊙²⁷. The continuous increase of the corresponding quiescent population is also detected by Moustakas et al. (2013) between z = 1 and z = 0.1 when they measured the weighted linear fits of n_∗(z).

Figure 16 presents the cosmic evolution of the stellar mass density ρ_∗ for all (black), star-forming (blue), and quiescent (red) galaxies. We compare our results (filled stars) with previous studies. We also plot the stellar mass density obtained by assuming a different slope of the star-forming SMF low-mass end (open stars; cf. Sect. 6.1), but it does not change the results significantly. In good agreement with Ilbert et al. (2013, circles)²⁸ and Tomczak et al. (2014, squares), our measurement of the global evolution of ρ_∗ reveals two phases: a >50% increase from z ~ 1.3 down to z ~ 1, and a continuous 12–20% increase from z ~ 1 down to z ~ 0.3.

As mentioned in Sect. 6.1, our selection of quiescent galaxies is more compatible with the selections of Ilbert et al. (2013) and Tomczak et al. (2014) when we consider that galaxies lying in the green valley are classified as quiescent. This corresponds to the upper red and lower blue envelopes of ρ_∗ in Fig. 16. Still, our measurement for quiescent galaxies is smaller than previous measurements by up to 25%. We do not find this difference when we consider the global stellar mass density. The importance of the Q/SF selection is reinforced by the fact that the agreement is better with Ilbert et al. (2013), when they use the log sSFR = −11 selection²⁹ (Fig. 16, open red circles). Our measurement is also consistent with the ρ_∗ measured by Arnouts et al. (2007,red crosses) for quiescent galaxies, which are selected thanks to SED-fitting (we do not plot their star-forming ρ_∗³⁰).

As previously suggested, the evolution of the stellar mass density of star-forming galaxies seems to be quite stable at z< 1.5. At the same time, a rapid increase of the stellar mass contained in quiescent galaxies is observed, increased by a factor > 2.5 from z ~ 1.3 down to z ~ 1. At lower redshift, we detect a small and continuous ≳ 30% increase of ρ_∗ from z ~ 1 down to z ~ 0.3, which reflects the progressive quenching of less massive galaxies.

7. Discussion

7.1. High-mass end evolution

As highlighted above, our sample can be used to investigate the evolution of massive (M_∗> 10^10.5 M_⊙) and rare (M_∗> 10^11.5 M_⊙) galaxies, thanks to the large volume of our survey. Most importantly, we are interested in the evolution of these objects across cosmic time, in particular to understand which mechanisms determine their evolution from star-forming to quiescent galaxies. Several studies (e.g. Kauffmann et al. 2003; Bundy et al. 2006; Davidzon et al. 2013) have characterised galaxy quenching with the so-called transition mass, which is the stellar mass at which the quiescent and star-forming populations are equal in a given redshift bin. In the same spirit, we define the transition redshift, z_tr, at which the quiescent population becomes dominant. As shown in Fig. 15, the transition redshift is found to be z_tr ≳ 1.4, z_tr ~ 0.9 and z_tr ~ 0.2, for M_∗ ~ 10^11.75 M_⊙, M_∗ ~ 10^11.25 M_⊙, and M_∗ ~ 10^10.75 M_⊙ galaxies, respectively: globally, the more massive a galaxy, the earlier its star formation is stopped. This is qualitatively consistent with the redshift evolution of the transition mass (e.g. see Davidzon et al. 2013). As already mentioned, several physical mechanisms could explain this trend within a hierarchical context (e.g. De Lucia & Blaizot 2007; Neistein & Dekel 2008; Weinmann et al. 2012). For instance, based on the stellar-halo mass relation from Coupon et al. (2015), the star-forming galaxies with stellar masses of ${\mathit{M}}_{\mathrm{SF}}^{\mathrm{\ast }}\mathrm{(}\mathrm{~}$10^10.64 M_⊙) should reside in dark matter halos of masses of around M_h ~ 10^12.4 M_⊙. This value agrees well with the halo mass threshold invoked by Cattaneo et al. (2006), corresponding to halo’ quenching, but we cannot exclude that some radio-AGN quenching could also explain why massive galaxies cease forming stars and/or are not fuelled anymore by fresh infalling gas (Croton et al. 2006).

We find that the number density of the most massive (M_∗> 10^11.5 M_⊙) galaxies almost doubled from z ~ 1 to z ~ 0.3 (Fig. 15). This corresponds to the <0.25 dex increase of the SMF high-mass end that is seen between z ~ 1 and z ~ 0.3 (Fig. 14). Because the high-mass end is dominated by quiescent galaxies at z< 1, the increase of the M_∗> 10^11.5 M_⊙ population cannot be explained by incidental star formation (Arnouts et al. 2007). If we assume that, in general, these very high-mass galaxies do not experience significant star formation, they can still assemble stellar mass through mergers at z< 1, in particular through dry merging.

7.2. Taming of galaxies

In Sect. 6.1 we have shown that the characteristic stellar mass of the star-forming SMF does not vary significantly between redshifts z = 0.2 and z = 1.5. As described in Sect. 5.1, we performed three selections of the SF galaxies, and the values of ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ differed slightly from one selection to another. In Fig. 17 we plot ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ as a function of the redshift and the SF galaxy selection in the NUVrK diagram. First, we find that ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ is between 10^10.6 and 10^10.8 M_⊙ at 0.2 <z< 1.5, regardless of the SF selection in the NUVrK diagram. More precisely, we find

• log ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}\mathit{/}{\mathit{M}}_{\mathrm{\odot }}\mathrm{=}\mathrm{10.6}{\mathrm{9}}_{-0.05}^{\mathrm{+}\mathrm{0.04}}$ if the galaxies in transition are included in the selection of SF galaxies (upper dotted lines in Fig. 8);

• log ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}\mathit{/}{\mathit{M}}_{\mathrm{\odot }}\mathrm{=}\mathrm{10.6}{\mathrm{6}}_{-0.03}^{\mathrm{+}\mathrm{0.02}}$ for our intermediate selection; and

• log ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}\mathit{/}{\mathit{M}}_{\mathrm{\odot }}\mathrm{=}\mathrm{10.6}{\mathrm{4}}_{-0.01}^{\mathrm{+}\mathrm{0.01}}$ for the most conservative selection.

Therefore, the evolution of ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ is consistent with being constant if the galaxies transitioning in the green valley are excluded from the selection of SF galaxies. The invariance with respect to redshift of ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ for the most conservative selection strongly supports a mass-quenching process occurring around a constant stellar mass, which makes this selection suitable for investigating the galaxies that are about to quench.

7.2.1. Tracking galaxies in the green valley

To identify a potential quenching channel for ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ galaxies, we isolate and characterise the green valley galaxies in Fig. 18, where each panel shows a different redshift bin. The contours represent the density of quiescent and star-forming galaxies, when the galaxies in transition are excluded (i.e. using the strictest selection of Q/SF galaxies). The colour code expresses the stellar mass. In the lower panels, we show the rest-frame (r−K_s)° distribution of the transitioning galaxies (i.e. the galaxies lying in the NUVrK green valley). As explained in Sect. 5.1, the NUVrK diagram is very efficient in separating dusty star-forming galaxies from quiescent ones (see Fig. 16 of the companion paper), which allows us to properly define transitioning galaxies in the green valley. We observe that

• 1

  the (r−K_s)° distribution of galaxies in transition is narrow and does not evolve with redshift (>80% of these galaxies have 0.76 < (r−K_s)°< 1.23); and that

• 2

  the typical stellar mass of galaxies in transition is around ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ (>60% of these galaxies have 10^10.5<M_∗/M_⊙< 10¹¹).

Therefore, we isolated the quenching channel of the ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$-galaxies with the colour criterion (r−K_s)°> 0.76 in the NUVrK green valley (green dashed lines in Fig. 18, sub-panels).

We also detect a clear plume of young quiescent galaxies in Fig. 18, with (r−K_s)°< 0.76 (i.e. bluer than observed for galaxies following the ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ channel) at z< 0.5. It is well established that rest-frame optical-NIR colours are sensitive to both dust attenuation and age of the stellar populations (see e.g. Whitaker et al. 2012). Under the assumption that, on average, the (r−K_s)° colour of quiescent galaxies cannot become bluer with time, the young part of the quiescent population should have used another quenching channel. According to the limit that we defined to isolate the ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ quenching channel (green dashed line in Fig. 18), we separated the young quiescent (Q_yng) and old quiescent (Q_old) galaxies with (r−K_s)° = 0.76. Figure 18 also reveals that Q_yng galaxies are characterised by relative low masses (M_∗ ≲ 10^9.5 M_⊙), which seems to match the low-mass upturn of the quiescent SMF (see Fig. 13) at z< 0.5. In Fig. 19 we compute the SMF for Q_yng (magenta triangles) and Q_old (dark red circles) galaxies at 0.2 <z< 0.5. The Q_yng galaxies dominate at low mass, and they are responsible for the low-mass upturn in the quiescent SMF. At the same time, the SMF of Q_old galaxies peaks at ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$, which clearly supports the idea that the building of the quiescent SMF is led through two quenching channels that can be distinguished with a cut in the NUVrK diagram at (r−K_s)° = 0.76. The timescale might then be a key element for characterising the mechanisms that are involved in each channel.

7.2.2. Quenching timescales

In Sect. 7.2.1 we have identified two possible channels in which galaxies are transitioning to build the quiescent population. We now investigate the nature of these channels through their characteristic timescales.

The restframe UV is sensitive to timescales of 10^-2−10^-1 Gyr, and the scarcity of young/low-mass galaxies in the green valley allows us to expect that some quenching processes occur on timescales of the same order or shorter. To better constrain the timescale of the quenching that affects the star formation of low-mass and ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ galaxies, we explored the behaviour of simple scenarios of star formation history (SFHs) within the NUVrK diagram in a similar way as the approach adopted by Schawinski et al. (2014). We performed this analysis at 0.2 <z< 0.5, where both old and young quiescent galaxies are well identified. The use of simple e-folding SFHs implies that we assumed that galaxies can only become redder with time. This is motivated by the fact that the fraction of quiescent galaxies has continuously increased between z ~ 3 and z ~ 0.2 (e.g. Ilbert et al. 2010; Muzzin et al. 2013; Mortlock et al. 2015) and by assuming that most green valley galaxies are transitioning for the first time (Martin et al. 2007). Doing so, we neglect the green valley galaxies produced by rejuvenation processes, as observed in the local Universe (e.g Salim & Rich 2010; Thomas et al. 2010) and recently predicted at higher redshift in the eagle simulations (Trayford et al. 2016). However, in these simulations, the rejuvenation is responsible for a small fraction of the green valley galaxies.

Figure 20 presents the resulting tracks in the NUVrK diagram for SFHs constructed in the same way: a continuous star formation up to the quenching time at t_Q, followed by an exponentially declining star formation characterised by τ. To mimic the average properties of our K_s< 22 sample at 0.2 <z< 0.5, the example is plotted for one metallicity (Z = 0.008), one extinction law (Calzetti et al. 2000), one value of the dust attenuation (E(B−V) = 0.2), and with a stellar age of at least 1 Gyr. The stellar age is colour coded, and only the ages allowed by the given redshift bin are plotted. In the left panel of Fig. 20 the SFHs are characterised by t_Q = 1 Gyr, with τ = 0.1, 0.25, 1, 2, and 2.5 Gyr. The tracks are constructed in a very simple way, and the evolution assumes a constant dust attenuation based on its average value for the bluest SF galaxies. The arrows show the shift that is due to a 0.1 increase of E(B−V). It is expected that quiescent galaxies are less affected by dust, which would tend to make the tracks steeper in the NUVrK green valley. Keeping this effect in mind, we see as a first result that the presence of Q_yng galaxies is expected if any quenching process occurs early (t_Q ~ 1 Gyr) with a typical timescale of τ ≲ 0.25 Gyr (triangles and squares in the left panel of Fig. 20 ). As a second result, τ = 1 Gyr (inverted triangles in Fig. 20 left panel) seems to be a lower limit for the quenching timescale that is compatible with the channel drawn by ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ galaxies. The galaxies with a quenching τ> 2 Gyr do not reach the quiescent cloud.

In the middle and right panels of Fig. 20, we also investigate the effect of the quenching epoch. We fixed τ = 1 Gyr and τ = 0.1 Gyr for several values of t_Q between 1 and 9 Gyr. Any t_Q> 9 Gyr will produce the same result as t_Q = 9 Gyr since the NUVrK colours of SF galaxies saturate at ages >9 Gyr, as shown by the predicted track with a continuous star formation (circles). All the models with τ = 1 Gyr are able to explain the galaxy presence within the ℳ^⋆ channel. We could also imagine that a shorter timescale combined with a late quenching time can reproduce the observed ℳ^⋆ channel. However, if we consider an SFH with τ = 0.1 Gyr after 9 Gyr on the SF main sequence (middle panel, inverted triangles), the track seems to move away from the channel that is drawn by ℳ^⋆-galaxies in the NUVrK diagram. To produce a track that is compatible with this channel, we need to consider a quenching timescale τ ≳ 0.5 Gyr (red solid line), regardless of the considered SHF. We recall that we have considered the shortest timescales compatible with the ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$-quenching channel, and we could pick out SFHs that agree better. Namely, SFHs characterised by t_Q = 1 Gyr and τ = 1.5 Gyr, t_Q = 5 Gyr and τ = 1 Gyr, or t_Q = 9 Gyr and 0.5 <τ< 1 Gyr could also explain the presence of this channel. This suggests a quenching timescale range of 0.5 <τ< 2 Gyr for ℳ^⋆-galaxies, which corresponds to a quenching duration of between ~1 and 3.5 Gyr³¹. Therefore, the physical mechanism explaining the building of the quiescent SMF around ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ at z< 1 seems to be a slow process. Such a mass dependent mechanism is compatible with a strangulation picture where the star formation quenching occurs on several Gyr, moving slowly away from the SF main sequence in the NUVrK diagram, while the gas supply is progressively halted (Schawinski et al. 2014; Peng et al. 2015).

Figure 20 shows that the plume formed by Q_yng galaxies in the NUVrK plan is explained by a ~0.1 Gyr-quenching process occurring during the first ~5 Gyr of the galaxy life (squares, diamonds, and triangles in the middle panel of Fig. 20). The absence of these low-mass galaxies lying in the green valley can be first explained by the rapidity of their quenching. Indeed, a galaxy quenching with τ = 0.1 Gyr (triangles in the left and middle panels of Fig. 20) is expected to cross the green valley (delimited by the black solid lines) in ~0.4 Gyr, while a galaxy with τ ~ 0.5−2 Gyr spends ~1−3.5 Gyr there, on average. Nevertheless, the potential reservoir of SF M_∗< 9.5 galaxies that can quench is aboutten times larger than for galaxies around ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ (cf. Fig. 13). We could then expect to see more low-mass galaxies in transition. By adopting a conservative approach, we can assume that the ratio between the two quenching timescales is ~10 (0.1 Gyr for M_∗< 9.5 galaxies, 1 Gyr around ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$). The corresponding quenching rate should consequently be aboutten times lower for the low-mass galaxies that are the progenitors of the Q_yng galaxies than for the ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ galaxies. The resulting flux of quenching galaxies (i.e. quenching rate × SF reservoir) that cross the green valley is then expected to be of the same order of magnitude, both at low and high mass, except when only a fraction of the low-mass galaxies is likely to be affected by the quenching. The SF satellite galaxies, which are estimated to be ≳3 times less abundant than field galaxies (Yang et al. 2009; Peng et al. 2012), are therefore good candidates for this low-mass quenching mechanism. Moreover, its typical timescale is compatible with the scenario suggested by Schawinski et al. (2014) for the rapid formation of young early-type galaxies. In this picture, the quiescent low-mass galaxies are formed through dramatic events such as major mergers and not through ram-pressure stripping or strangulation, by explaining both the almost instantaneous star formation shutdown and the morphological transformation.

8. Summary

We analysed the evolution of the stellar mass function in an unprecedentedly large (>22 deg²) NIR selected (K_s< 22) survey. This allowed us to provide reliable constraints on the evolution of massive galaxies and to investigate quenching processes below redshift z ~ 1.2. Covering the VIPERS spectroscopic survey, we computed highly reliable photometric redshifts, with usual estimates of the precision σ_{Δz/ (1 + z)}< 0.03 and σ_{Δz/ (1 + z)}< 0.05 for bright (i< 22.5) and faint (i> 22.5) galaxies, respectively.

Paying particular attention to several sources of uncertainties (photometry, star-galaxy separation, photometric redshift, dust extinction treatment, and classification into quiescent and star-forming galaxies), we computed the SMF between redshifts z = 0.2 and z = 1.5. The unique size of our sample enabled us to drastically reduce the statistical uncertainties affecting the SMFs and stellar mass densities with respect to other current surveys over the stellar mass range we consider: the Poissonian error and cosmic variance are reduced by factors of ~3.3 and ~2, respectively, compared to a 2 deg²-survey. Combined with a careful treatment of the Eddington bias that is due to the stellar mass uncertainty, we produced an unprecedentedly precise measurement of the massive end of the SMF at z< 1.5. In particular, we stress the importance of constraining all sources of systematic uncertainties, which quickly become the dominant sources of error in large-scale surveys such as those planned with Euclid or LSST.

Using the (NUV−r) versus (r−K) rest-frame colour diagram to classify star-forming and quiescent galaxies in our sample, we measured the evolution of the SMFs of the two populations and investigated the possible quenching processes that could explain the build-up of the quiescent population. Our main conclusions are summarised below.

• 1)

  We provided clear evidence that thenumber density of the most massive(M_∗> 10^11.5 M_⊙) galaxies increases by a factor ~2 from z ~ 1 to z ~ 0.3, which was first highlighted by Matsuoka & Kawara (2010). This population is largely dominated by the quiescent population since z ~ 1, allowing for the possibility of galaxy mass assembly through dry-mergers in very massive galaxies.

• 2)

  The characteristic mass of the SF population was found to be very stable in the redshift range 0.2 <z< 1.5, with $\log \mathrm{(}{\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}\mathit{/}{\mathit{M}}_{\mathrm{\odot }}\mathrm{)}\mathrm{=}\mathrm{10.64}\mathrm{\pm }\mathrm{0.01}$. This confirms that the star formation is impeded above a certain stellar mass (Ilbert et al. 2010; Peng et al. 2010).

• 3)

  Using the NUVrK diagram as a tracer of the galaxy evolution, we identified one main quenching channel between the star-forming and quiescent sequences at 0.2 <z< 1.5, which is followed by galaxies with stellar masses around ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$. This channel is characterised by a colour (r−K_s)°> 0.76, typical of evolved massive star-forming galaxies, which should feed the majority of the quiescent population. We also identified a young quiescent population with (r−K_s)°< 0.76, whose galaxies likely followed another path to reach the quiescent sequence. We showed that this blue quiescent population, dominated by low-mass galaxies, is responsible for the upturn of the quiescent SMF at low redshift.

• 4)

  Assuming simple e-folding SFHs (galaxies can only become redder with time), we found that the ${\mathrm{ℳ}}_{{sf}}^{\mathit{\star }}$ channel is explained by long quenching timescales, with 0.5 <τ ≲ 2 Gyr. Galaxies in this channel are expected to turn quiescent after ~1−3.5 Gyr on average. This is compatible with strangulation processes occurring when the gas cooling or the cold gas inflows are impeded, allowing the galaxy to progressively consume its remaining gas reservoir (Peng et al. 2015). Conversely, the quenching of low-mass galaxies that is visible at low redshift is characterised by short timescales with τ ~ 0.1 Gyr. This quenching that halts star formation in ~0.4 Gyr can be consistent with major merging (Schawinski et al. 2014) and may preferentially affect satellite galaxies.

---

Acknowledgments

We gratefully acknowledge the anonymous referee, whose advice helped much in improving the clarity of the paper. We wish to thank J.-C. Cuillandre, H. Aussel, S. de la Torre and B. C. Lemaux for helpful discussions. We would also like to thank J. Moustakas for providing the PRIMUS stellar mass estimates. This research is in part supported by the Centre National d’Études Spatiales (CNES) and the Centre National de la Recherche Scientifique (CNRS) of France, and the ANR Spin(e) project (ANR-13-BS05-0005, http://cosmicorigin.org). L.G. acknowledges support of the European Research Council through the Darklight ERC Advanced Research Grant (# 291521). This paper is based on observations obtained with MegaPrime/MegaCam, a joint project of CFHT and CEA/DAPNIA, and with WIRCam, a joint project of CFHT, Taiwan, Korea, Canada and France, at the Canada-France-Hawaii Telescope (CFHT). The CFHT is operated by the National Research Council (NRC) of Canada, the Institut National des Sciences de l’Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. This work is based in part on data products produced at Terapix available at the Canadian Astronomy Data Centre as part of the Canada-France-Hawaii Telescope Legacy Survey, a collaborative project of NRC and CNRS. We thank the Terapix team for the reduction of all the WIRCAM images and the preparation of the catalogues matching with the T0007 CFHTLS data release. This paper is based on observations made with the Galaxy Evolution Explorer (GALEX). GALEX is a NASA Small Explorer, whose mission was developed in cooperation with the Centre National d’Études Spatiales (CNES) of France and the Korean Ministry of Science and Technology. GALEX is operated for NASA by the California Institute of Technology under NASA contract NAS5-98034. This paper uses data from the VIMOS Public Extragalactic Redshift Survey (VIPERS). VIPERS has been performed using the ESO Very Large Telescope, under the “Large Programme” 182.A-0886. The participating institutions and funding agencies are listed at http://vipers.inaf.it. This paper uses data from the VIMOS VLT Deep Survey (VVDS) obtained at the ESO Very Large Telescope under programs 070.A-9007 and 177.A-0837, and made available at the CESAM data center, Laboratoire d’Astrophysique de Marseille, France. Funding for PRIMUS is provided by NSF (AST-0607701, AST-0908246, AST-0908442, AST-0908354) and NASA (Spitzer-1356708, 08-ADP08-0019, NNX09AC95G). Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the US Department of Energy Office of Science. The Participating Institutions of the SDSS-III Collaboration are listed at http://www.sdss3.org/.

References

All Tables

All Figures