© ESO, 2011

1. Introduction

Understanding the chronology of bulge and disc formation at different cosmological epochs is a fundamental goal of observational cosmology. Different paths have been proposed for the formation of the spheroidal and disc galactic components including some simple scenarii like: bulges form first and discs are added later by accretion of material; bulges form out of discs by secular evolution; ellipticals form through mergers of spirals and it is likely that the reality is a complex mix of different processes. Recent simulations have pointed out that there might be a cycle in galaxy formation and evolution where spheroids, formed in gas rich mergers, can reform gas discs which could form a bar or spiral arms and reactivate star formation by gas funneled to the center (e.g. Governato et al. 2009).

The knowledge of the distribution of light in bulges and discs and its evolution would provide important constraints on competing scenarios of galaxy formation and evolution, as the ratio of spheroid luminosity to total luminosity measures e.g. the efficiency of the first burst of star formation relative to later slow accretion (Fall & Efstathiou 1980; Kauffmann et al. 1993). Measuring the distribution of light in bulges and discs is not an easy task. Observations show that galaxies are complex systems. Many observational and theoretical studies simply describe galaxies as constituted of two components with different photometric and dynamical properties: a thin, rotationally supported stellar disc which often also contains gas and dust, and a bulge, a spheroidal or ellipsoidal component made purely of stars. However, the existence of galaxies with a central bar–shaped component was already taken into consideration by Hubble in his first classification system and while a classification scheme matches the need to distinguish among different classes of objects, galaxies actually come in very different shapes. Many galaxies possess a disc, a bulge and a bar, but some, the elliptical galaxies, have no significant disc, while others, late-type spirals and irregulars, have little or no bulge. In the local universe there is a general consensus that the barred fraction, defined as the number of barred spirals over the whole spiral population, is  ~0.3 for strongly barred systems and  ~0.65 for all barred galaxies. High redshift estimates have instead proven difficult with measurements ranging from a dramatic paucity of barred galaxies at z > 0.5 to a constant bar fraction at 30% over the last 8 Gyr. Nonetheless, for L^∗ and brighter spirals the debate has been resolved by Sheth et al. (2008) who reconciled various studies, showing that the bar fraction declines.

Important progress in understanding how galaxies form and evolve has been achieved from both a theoretical and observational point of view. N-body simulations have reached the necessary resolution to show that bars not only form naturally in spiral galaxies but that the inflow of material to the center due to bar torques can grow bulges (Athanassoula & Misiriotis 2002; Athanassoula et al. 2005). Cosmological hydrodynamical simulations have shown that the role of the AGN feedback should be taken into consideration to match the observed number density of red sequence galaxies in the local universe (Khalatyan et al. 2008). Utilizing the HST/WFC3 imaging, van der Wel et al. (2011) have been able to investigate the internal structures and surface brightness profiles of massive, quiescent galaxies at z ~ 2 revealing that the great majority might be disc-dominated. Mid-infrared and sub-millimeter data (Viero et al. 2010) seem to indicate that these compact massive galaxies at z ~ 2 are actually actively star forming. New infrared and submillimiter observations as well as kinematical studies in the local universe (i.e. de Zeeuw et al. 2002; Sheth 2009; Cappellari et al. 2011) have also shown the complexity of what were once expected to be “red and dead” ellipticals, which instead turn out to have disky central substructures, star formation, molecular gas and dust. This complexity is difficult to take into consideration when modelling galaxy surface brightness profiles, specifically when comparing low and high redshift studies which do not present the same richness of details.

In the local universe galaxies of differing morphologies are segregated according to environmental density (Dressler 1980) with the low-density field composed largely of spirals and irregulars and the densest regions of clusters composed of lenticulars and ellipticals. Large spectroscopic redshift surveys, exploring the high redshift universe (i.e. VVDS, Le Fèvre et al. 2005; DEEP2, Davis et al. 2003; Faber et al. 2007; zCOSMOS, Lilly et al. 2007), have provided an accurate estimate of the environment from which the color-density (Cucciati et al. 2006; Cooper et al. 2006) and morphology-density (Tasca et al. 2009) relations have been shown to evolove significantly since z ~ 1. Many different interpretations have been proposed to explain this behavior, ranging from initial condition biases which imprint the differences at birth, to gravitational, gas dynamical or radiative processes through which galactic environment affects later evolution. Various studies using the large, complete samples provided by the Sloan Digital Sky Survey (SDSS, see below) show the distribution of galaxy mass to depend quite strongly on environment, apparently reflecting an initial condition bias, while for galaxies of a given mass, star formation properties depend much more strongly on environment than structural properties do, suggesting that external processes primarily affect the gas component from which stars form (Kauffmann et al. 2004).

The relative contributions of bulges and discs to the stellar content of galaxies of differing mass as well as to the overall stellar content of the Universe are clearly important quantities which should be reproduced within any viable picture of galaxy formation. It has long been known that the massive galaxy population is dominated by light from spheroids, while discs tend to dominate in lower mass populations (e.g. Efstathiou et al. 1982). The most commonly cited estimates of the fractions of all stars in bulges and discs trace back to the work of Simien & de Vaucouleurs (1986) and Schechter & Dressler (1987), finding that discs contribute roughly twice as much as bulges to the mean luminosity density of the universe. In these studies the use of a de Vaucouleurs bulge for every spiral galaxy has an impact on the bulge to disc ratio. As shown in Graham (2001) the B/D flux ratio is over-estimated if the bulge actually has a profile with n < 4, as we now know to be the case for a fraction of late-type spirals. The estimate of the fraction of baryons locked up in the spheroid and disc galactic components computed by Fukugita et al. (1998) using Gunn r-band CCD photometry agreed with previous findings. More recent studies (Benson et al. 2002, 2007; Driver et al. 2007; Durbala et al. 2008; Gadotti 2009; Oohama et al. 2009, and Weinzirl et al. 2009) also provide accurate measurements of the fraction of stellar mass in the various structural components of massive galaxies in the local universe. Some of these studies focused on the investigation of the fraction of stellar mass or luminosity density from specific morphological types (e.g. bright spirals in Weinzirl et al. 2009; Sb and Sc galaxies in Durbala et al. 2008) or environment (e.g. isolated galaxies in Durbala et al. 2008), therefore missing the global perspective from a complete galaxy sample.

The measurement of the light distribution in discs and bulges is best done using image decomposition techniques which fit specific light distribution models to imaging data. These techniques have been detailed by several authors who have made publicly available some popular decomposition codes like GIM2D (Simard et al. 2002), GALFIT (Peng et al. 2002) or BUDDA (Gadotti 2009). While some recent studies (i.e. Gadotti 2009; Weinzirl et al. 2009) have implemented increasingly sophisticated decompositions adding complex morphological components, a robust local reference of the distribution of light in the two major components of galaxies, bulges and discs, is still missing. Following the growth of these two main morphological components of galaxies along cosmic time remains constrained by the ability to identify those components at increasingly high redshifts. Sheth et al. (2008) show that HST/ACS images can be used to accurately identify the overall stellar structure (e.g., bars) in galaxies up to z ~ 0.85, but to quantify the stellar luminosity they contain at high redshift stays a not trivial exercise. Another limitation of this kind of photometric studies is that photometric components such as bulges and discs have distinct dynamical counterparts and the correspondence between the two does not always hold. The modelling of a photometric disc component does not imply the presence of an actual disc, since dynamically hot systems also have simple exponential profiles. On the other hand, a photometric bulge might not represent a genuine dynamically hot spheroid.

Here we present a simple and robust two components bulge and disc decomposition, and we aim to quantify the global fraction of light in each component from a representative galaxy sample and to produce a local reference for high redshift studies. We apply state of the art two-dimensional morphological decomposition to the images of a magnitude-limited sample of relatively nearby (median redshift of 0.05) and bright (r ≤ 15.9) galaxies with photometry available from SDSS. We use the results to study the overall contribution of discs and bulges to the light of the nearby galaxy population both as a function of intrinsic galaxy luminosity and for the population as a whole. In Sect. 2 we describe the sample. In Sect. 3 we describe the decomposition technique and the parametric functions it uses to fit the galaxy and we describe how we apply it to the data. images. In Sect. 4 results are compared for different photometric bands and for different parametric fitting functions in order to explore the robustness of the results. Deviations of the images from the best fit models are also quantified. In Sect. 5 we study the distribution of derived parameters for those galaxies which could be successfully fitted. We also identify and correct for a serious (and previously known) systematic which results in the assignment of significant disc light to bright galaxies where discs are in most cases actually absent. This section presents our main results, with estimates of the fraction of the total light of galaxies with a given luminosity which are in discs, in bulges, or in galaxies with no detectable disc or bulge. Combining these with previous measurements of the galaxy luminosity function allows us to obtain corresponding fractions for the galaxy population as a whole. A final section summarises our results and discusses them in the light of previous findings.

Throughout this paper, unless otherwise stated, we assume a Hubble constant of H₀ = 70 km s^-1 Mpc^-1 and an Ω_M = 0.3, Ω_Λ = 0.7 cosmology in calculating distances and luminosities.

2. Observational data

2.1. Galaxy sample

In June 2001 the SDSS (York et al. 2000) released its Early Data Release (EDR; Stoughton et al. 2002), roughly 462 square degrees of imaging data collected in drift scan mode. The imaging is conducted in the u-, g-, r-, i- and z-bands (Fukugita et al. 1996; Gunn et al. 1998; Hogg et al. 2001; Smith et al. 2002; Pier et al. 2003). The reader is referred to Ivezić et al. (2004) for details on the photometric quality assessment. The EDR contains around a million galaxies distributed within a narrow strip of 2.5 degrees across the equator. As the strip crosses the galactic plane, the data are divided into two separate sets in the North and South Galactic caps. The SDSS had the ambitious goal to image a quarter of the Celestial Sphere at high Galactic latitude as well as to obtain spectra uniformly for all the brighter galaxies. For the present project this has the advantage, in comparison to previous work, of having uniform photometry and spectroscopy over a much larger area, permitting a major improvement in sample size and homogeneity.

In the following analysis we are using a sample of galaxies defined by the Japanese Participation Group (JPG, Nakamura et al. 2003). This is an homogeneous sample obtained from the northern equatorial stripes of the SDSS EDR. The region of the sky covered is an area of 229.7 square degrees at 145.15° ≤ α ≤ 235.97° and  −1.27° ≤ δ ≤ 1.27°. The sample is complete to a magnitude limit of r ≤ 15.9 after Galactic reddening correction. The use of a magnitude-limited sample, selected from a large volume, is required in order to provide an unbiased description of the whole galaxy population, representative of different kinds of environments and of different galaxy morphological types. Our bright cut is additionally motivated by the fact that eye classifications cannot be made confidently beyond this magnitude, and for our current purpose this bright limit has the advantage that the galaxies are all large compared to the SDSS point-spread function (PSF). This is another fundamental aspect when performing morphological studies. Since the mean PSF full width at half maximum (FWHM) (commonly named “seeing”) in SDSS images is  ~1.4 arcsec and the pixel scale is 0.396 arcsec/pixel, we are Nyquist sampled. We are therefore assured that we can estimate dimensions about equal to or somewhat smaller than the FWHM. Beijersbergen et al. (1999) have also shown that a reliable decomposition of galaxies into a bulge and disc component requires the seeing to be smaller than the typical galaxy half-light radius. Since all galaxies in our sample have indeed half-light radii larger than 1.4 arcsec, we are therefore additionally reassured on the suitability of our sample for the study we want to perform.

Fig. 1 The distributions of galaxies in morphological type (upper panel) and in extinction-corrected Petrosian apparent magnitude in the r band (second panel) are shown for the photometric and spectroscopic samples, empty and dashed-line respectively. The redshift (third panel) and corresponding r band absolute magnitude (lower panel) distributions are shown for the galaxies in the spectroscopic sample.

All the 1862 galaxies in the sample were classified by eye on the system of the Hubble Atlas of Galaxies (Sandage 1961) by JPG scientists using the g band image of each galaxy. For each galaxy the final quoted type is the mean of 4 independent classifications by different scientists. The rms of these 4 classifications is also given; they typically agree within ΔT ≤ 1.5. The corresponding numerical classification as defined in the Third Reference Catalogue of Bright Galaxies (de Vaucouleurs et al. 1991) is also reported. The seven resulting subsamples separate galaxies according to morphology going from ellipticals (E) through lenticulars (S0) and early-type spirals (Sa, Sb), to late-type spirals (Sc, Sd) and irregulars (Irr). Those galaxies which do not enter any of the aforementioned classes are flagged as “unclassified”. Finally, for 1588 galaxies out of our sample of 1862 we have spectroscopic information. A summary of the distribution of our galaxies across the morphological classes is given in Table 1 and in the first panel of Fig. 1. The fractional morphological composition of our sample is E:(E/S0-S0):(S0/a-Sab):Sb:Sc:Sd:Im  ~0.13:0.27:0.26:0.16:0.12:0.03:0.01 in good agreement with the literature (e.g. Fukugita et al. 1998). The somewhat slightly higher fraction of E/S0’s is due to the primary r-band selection of our sample, a redder band if compared to the standard B-band selection used by previous authors. For the same reason, the fraction of Im galaxies is smaller by a factor of 2 − 3 than that from B-selected samples. We refer to Fukugita et al. (2007) and Nakamura et al. (2003) for further details. Further work exploring the properties of this sample can be found in Nakamura et al. (2003, 2004), Kelly & McKay (2004), Fukugita et al. (2004) and Oohama et al. (2009).

2.2. Photometric and spectroscopic data

Two important quantities used in this paper are taken directly from the SDSS database: the redshift and the Petrosian magnitude. The first is obtained by the spectroscopic pipelines idllspec2d (written by Schlegel & Burles) and spectro1d (written by SubbaRao, Bernardi and Frieman). A description of the tiling algorithm used to assign targets to each pointing is given in Blanton et al. (2003c). The SDSS spectroscopic galaxy samples consist of all galaxies brighter than r = 17.77 (Strauss et al. 2002) and of a sample of luminous red galaxies (Eisenstein et al. 2001) extending at r < 19.2. The distribution of galaxies with respect to z for our sample is shown in the third panel of Fig. 1. The second quantity is obtained by the Photo pipeline (see Lupton et al. 2001, 2002) and it is based on a modified form of the Petrosian system for galaxy photometry which is designed to measure a constant (and quite large) fraction of the total light of a galaxy independent of its characteristic surface brightness but it does depend on the profile shape of a galaxy. Stoughton et al. (2002) estimate that 80% of the flux is recovered for the de Vaucouleurs profile galaxies and close to 100% for exponential galaxies. Three related quantities also used in this paper are R₅₀ and R₉₀, defined as the radii which include respectively 50 and 90 percent of the Petrosian flux (see Stoughton et al. 2002), and the concentration index c ≡ R₉₀/R₅₀. The second panel of Fig. 1 shows the distribution of galaxies with respect to their Petrosian magnitude after correction for foreground Galactic extinction using the reddening map of Schlegel et al. (1998). Such extinction-corrected Petrosian magnitudes are used throughout this paper.

3. Image analysis

3.1. The fitting algorithm

We examine the structural properties of our galaxies using Gim2D (Simard et al. 2002), a two-dimensional photometric decomposition algorithm which fits each image to a superposition of an elliptical component with a Sérsic profile, representing the bulge, and a concentric elliptical component with an exponential profile, representing the disc. It is important to recognise that this separation into bulge and disc is based only on the observed image, and may not correspond to the “best” decomposition if additional information, for example from kinematics of the stars or the gas, is taken into account.

It is well known that the bulge component of most galaxies can be well represented by a surface brightness profile of Sérsic form: (1)where Σ(r) is the surface brightness at a distance r from the centre along the semi-major axis and Σ_e is its characteristic value, the effective surface brightness, defined as the value at the effective radius r_e. The parameter b is related to the Sérsic index n and is approximately set equal to 1.9992n − 0.3271 for 0.2 < n < 10 so that r_e is the projected radius enclosing half of the total light (Sersic 1968; Ciotti 1991). Many authors fit bulges and ellipticals with a more specific function, the de Vaucouleurs r^1/4 law, which is obtained by setting n = 4. When fitting Sérsic profiles in the following we will assume 0.2 < n < 10. This choice is driven by the knowledge that fits to low-luminosity ellipticals generally give n significantly smaller than 4 (Binggeli & Jerjen 1998; Graham & Guzmán 2003), while luminous elliptical galaxies do have n > 4 (Caon et al. 1993). Values of n in excess of 4 are also found for cD galaxies (Graham et al. 1996). In all these studies the whole galaxy surface brightness distribution has been modelled using a Sérsic profile and no attempt has been done to perform a bulge-disc decomposition since the galaxy samples used only included ellipticals and/or lenticular galaxies.

Discs are generally well described by an exponential profile (corresponding to n = 1), although non-axisymmetric features due to bars, spiral arms and dust lanes can be large. We use the standard parametrisation: (2)where Σ₀ is the central surface brightness and h the disc scalelength.

Despite being purely empirical fitting functions, with no strong theoretical justification, these laws have shown to provide robust modelling of galaxy profiles.

With this model there is a maximum of twelve parameters which are adjusted in fitting the galaxy image and that we retrieve as output from our decomposition: the total flux of the object; the bulge-to-total light ratio B/T, defined as the fraction of the total flux in the bulge component so that B/T = 1 corresponds to a pure bulge and B/T = 0 to a pure disc; the bulge effective radius r_e; the disc scalelength h; the disc inclination angle i defined so that i = 0 for face-on discs and i = 90 for edge-on ones – the disc axial ratio is then (b/a)_disc = cos(i); the bulge ellipticity e given in terms of the bulge axial ratio by e = 1 − (b/a); the bulge and disc position angles (hereafter PA) measured clockwise from north and allowed to be different; the Sérsic index n which we sometimes fix at n = 4; the x-y pixel shifts dx and dy of the galaxy centre position in the model and science thumbnail images; and the background intensity level.

Additional parameters could be introduced to model other features (e.g. bars, spiral arms, etc.) but such decompositions become somewhat arbitrary and may not converge to unique solutions. Our current choice is standard and we found it to be a good compromise between stability of results and flexibility of representation. Notice that for an axisymmetric galaxy the position angles of the bulge and the disc would be the same. Allowing them to differ makes it possible for the code to detect triaxial bulges or bars in suitably oriented galaxies. Even if to detect bars is out of the scope of our analysis, we have used this property to check that Gim2D does not increase the retrieved B/T ratio because of bars, as we have described in Sect. 4.1.

3.2. Image reduction

We perform our analysis starting from corrected frames of area 13′.52 × 8′.98 taken directly from the SDSS archive. In these large-scale images flat-field, bias, cosmic-ray, and pixel-defect corrections have already been applied. The pixel-size in these images is 0.396 arcsec. To proceed with our fitting, we begin with a list of source positions and apply the SExtractor galaxy photometry package version 2.2.2 (Bertin & Arnouts 1996) to each field to estimate the local sky background level at each point and to define the isophotal area where each object is above the detection threshold (we choose a threshold which is higher than background by 1.5 times the background noise). When SExtractor performs galaxy photometry, it constructs a segmentation (or mask) image in which pixels belonging to the same object all have the same value and sky background pixels are flagged by zeros. Our 2D image fit is carried out on all pixels belonging to the same SExtractor-defined object. In practice Gim2D uses thumbnail images, two for each galaxy, extracted around the object of interest. The area of these is chosen to be 20 times larger than the mask area defined by SExtractor. The first thumbnail is cut from the corrected frame and is corrected for the local background estimated by SExtractor, while the second contains the corresponding pixels from the mask image. As many other fitting algorithms, Gim2D needs initial first guess values to start the computation as well as an assigned specific range for each fitting parameter. We have determined these quantities with SExtractor: for each galaxy the photometric value computed by SExtractor (i.e. magnitude, position angle, half-light radius, etc.) is used as initial first guess quantity; the range instead is chosen to be the same for all the galaxies and to be large enough around the mean value of the distribution. Some fine tuning was additionally applied according to our experience. The fitting algorithm then produces values and uncertainty ranges for each model parameter. Subtraction of the best-fit model from the input thumbnail produces a residual image that can be used to characterise how well the model fits the galaxy. Examples of all these images for representative galaxies from each of our morphological classes are shown in Fig. 2.

Fig. 2 Examples of the science image, the mask, the model and the residual images for our two fits (de Vaucouleurs plus exponential or Sérsic plus exponential) in the r band (from left to right) for galaxies from our seven morphological classes: E, S0, Sa, Sb, Sc, Sd, Irr (from top to bottom).

When the fitting algorithm starts to sample the 12-dimensional parameter space, it considers not only the pixels assigned to the main galaxy by the mask but all pixels flagged as object or background in the SExtractor segmentation image. Important information about the galaxy could be contained in the pixels below the detection threshold. In the residual images (see Fig. 2) one can see that all pixels belonging to other objects in the vicinity of the one under consideration are masked out. The final flux is obtained by the integration of the best fit model over all pixels, assuring that we do not lose the flux in the masked regions.

The model image of each galaxy is convolved with a PSF before comparison with the real data. The PSF can be highly variable across a corrected frame (Stoughton et al. 2002) and for this reason it is important to interpolate the PSF parameters measured for individual stars to the position of each galaxy before convolving with the model. Our galaxy light model is thus the sum of an exponential disc and a Sérsic, or de Vaucouleurs, bulge, convolved with this “best” PSF.

Finally a photometric calibration and the redshift are required to retrieve physical quantities from the output of the B/D decomposition code.

4. Accuracy of the decompositions

We have carried out the above fitting procedure entirely independently for each of our galaxies in the r and the i bands. The two segmentation images differ slightly and in addition there are colour variations across many of our galaxies as a result of variations in the underlying stellar populations and in the dust distribution. It is thus reassuring that the structural parameters in the two bands are in good agreement in the great majority of cases (see below). This shows that the fitting procedure produces stable results. In addition Gim2D produces acceptable converged parameter sets for almost all the galaxies in both bands. In Table 2 we show the number of galaxies from our spectroscopic sample that are successfully modelled by the code in each band in the two cases when the Sérsic index is set equal to 4 and when it is allowed to float. We consider a fit to be successful when the code is converging (see Sect. 4.3 for details). As expected slightly more galaxies can be fit when n is kept free, and in this case only  ~50 of our photometric sample (~1900) galaxies cannot be fit acceptably in either of the two bands. These galaxies are almost all later type spirals (Sb-Sc-Sd) which are not modelled simply because the centroid position provided by the SDSS database does not match the one obtained with SExtractor within the defined 6 arcsec searching radius. It is interesting to notice that the modelling also fails occasionally for early-type objects when these are forced to follow a de Vaucolueurs law for the central photometric component. We now discuss aspects of these fits in more detail.

4.1. Comparison between fits

As already noted, we have fitted all our galaxies with a bulge model in which the Sérsic index n is free and also with a model in which it is fixed to the de Vaucouleurs value n = 4. The main reason for using a de Vaucouleurs profile is to allow comparison with the literature, since many previous studies made large use of this specific profile. In this subsection we show that the parameters of the decomposition of most interest to us are only weakly affected by this choice for most galaxies. In Fig. 3 we compare the values of the bulge fraction B/T, of the disc scale length h, the bulge effective radius r_e and of the galaxy half-light radius r_hl obtained for each galaxy in the two cases. We show results for the i band only (results for r are similar). It is evident from the scatter plots that the measurements agree with a dispersion of only a few percents, as better shown (at the bottom of each panel) by plotting the normalised differences between the two measurements as function of one of the parameters. It has been claimed (Balcells et al. 2003) that the innermost regions of galaxy profiles play a decisive role in determining the Sérsic index of the bulge and that extra central light results in higher n profile fits. Components such as inner discs, compact sources and star formation present in galaxy nuclei (Phillips et al. 1996; Carollo et al. 2002; Rest et al. 2001; Pizzella et al. 2002) might be responsible for light which cannot be separated from the spheroid (Ravindranath et al. 2001). Our results show that these innermost components are either modelled (inner discs) or found in the residual image after decomposition. We are therefore confident that our estimate of the Sérsic index of the bulge is not affected by this effect. Recently, the work by Laurikainen et al. (2004, 2005), claimed that the light attributed to the bulge model might be overestimated for galaxies with strong non-axisymmetric structures, such as bars, if the bar component is omitted in the decomposition. In line with Laurikainen et al. (2004, 2005), the fact that our decomposition method does not take bars into account, but only bulges and discs, might erroneously led to attribute the bar light to the bulge component. After extensive tests, we are now confident that Gim2D does an excellent job in disentangling bars light from the bulge contribution. In particular, we observe that in our sample the sérsic index of bulges in barred and unbarred galaxies follows the same distribution and that there is no evidence that barred galaxies have a B/T distribution biased towards high values. We therefore conclude that the B/T of barred galaxies is not artificially increased and biased towards higher B/T ratios, and that our results are not affected by the use of standard bulge-disc decomposition instead of bulge-disc-bar decomposition.

As already mentioned in Sect. 3.1, we do not include additional components in the photometric model to fit bar-like structures. One of the main reason is that, for poorly resolved galaxies, additional components might cause further systematics. We therefore tested whether Gim2D uses the bulge component, in default of more suitable solutions, to compensate this lack and erroneously attributes to the bulge the light coming from stars in the bars. This would increase the B/T and overestimate the bulge luminosity density. If Gim2D uses highly elliptical bulge component to fit bar features in galaxy discs, our measurements should include unrealistically elongated spheroid components. We consider as upper limit for ellipticity the one of the most elliptical observed galaxies (Lambas et al. 1992). Only 1.7% of galaxies in our sample have bulges with ellipticity higher than 0.83 and the mean B/T of these galaxies is  ~0.14. This means that they are disc dominated and that therefore any contribution of bars to the bulge light of these galaxies would be a very small fraction of the total luminosity budget. We conclude that we can confidently discard a possible bias against disc light contribution caused by an erroneous fit of bar structures by the spheroid component.

Fig. 3 The bulge fraction (upper-left), the disc scalelength (lower-left), the effective radius of the bulge (upper-right) and the half-light radius of the galaxy (lower-right), obtained using different parametric functions to perform the decomposition (de Vaucouleurs plus exponential on the x-axis and Sérsic plus exponential on the y-axis) are plotted against each other for the 1680 galaxies (in the photometric sample) successfully modelled by the two fitting functions in the i band.

4.2. Comparison between different bands

In this section we compare the parameters estimated for each galaxy when the same model is fit independently to images in each of the two SDSS bands analysed here. In Figs. 4 and 5 we plot the retrieved parameters in the Sloan i and r bands against each other for the 1636 galaxies modelled successfully in both bands assuming a de Vaucouleurs profile for the bulge and an exponential for the disc. The corresponding results are shown in Figs. 6 and 7 for the 1742 galaxies modelled successfully in both bands assuming a Sérsic profile for the bulge and an exponential for the disc.

Among various parameters, Gim2D also computes the flux in the pixels belonging to the galaxy. We observe that the derived magnitudes, in the two photometric bands considered, are in very good agreement with the small scatter being consistent with the expected variation in mean galaxy colour. The scatter in bulge-to-total light ratio is gratifyingly small (<1%) for the vast majority of the galaxies, particularly when n is fixed to the de Vaucouleurs value. The scatter is slightly larger (≈2%) when n floats because different best fit values of the index in the two bands lead to different splits of the luminosity between bulge and disc. In most cases, however, similar n values are found in the two bands. The “arrow” shape of the B/T plots reflects the fact that, in few extreme cases, the data in one band occasionally prefer a weak disc while no disc is present in the best fit in the other band.

There is also quite good agreement, again a scatter of the order of a few percents, between the values of the scalelengths for the bulge and disc components measured in the two bands. The agreement becomes slightly worse for bulges than for discs and for components of small angular size compared to larger ones. This is presumably a reflection of resolution problems due to the finite pixel size and to difficulties with the PSF deconvolution. Nevertheless the apparent axial ratios of both bulges and discs agree well in the two bands with the scatter increasing for rounder systems. The agreement between the position angles of the disc and of the bulge in the considered photometric bands is also good with a scatter of only few percents. In addition we check that the distribution of the position angle is consistent with a random orientation on the sky.

Fig. 4 The bulge fraction (upper-left), the disc scalelength (lower-left), the effective radius of the bulge (upper-right) and the half-light radius of the galaxy (lower-right) in the i and r bands are plotted against each other for the 1638 galaxies (in the photometric sample) modelled in both bands by fitting a de Vaucouleurs profile to the bulge and an exponential to the disc.

Fig. 5 The ellipticity of the bulge (upper-left), the inclination angle of the disc (lower-left), the position angle of the bulge (upper-right) and the position angle of the disc (lower-right) in the i and r bands are plotted against each other for the 1638 galaxies (in the photometric sample) modelled in both bands by fitting a de Vaucouleurs profile to the bulge and an exponential to the disc.

Fig. 6 The same as Fig. 4 but for the 1742 galaxies (in the photometric sample) modelled in both bands by fitting a Sérsic profile to the bulge and an exponential to the disc.

4.3. Error estimates and goodness of fit

When carrying out its fitting Gim2D constructs a χ² value for each PSF-convolved model by summing over all pixels within the mask the square of the difference between model and data divided by the variance of the pixel noise, assumed to be due entirely to photon statistics. This measure of goodness of fit is then minimised over all parameters (each required to lie within a prespecified “allowed” range) to locate the maximum likelihood model. Once the algorithm has converged, the region of parameter space surrounding the likelihood maximum is sampled in order to compute marginalised a posteriori one-dimensional probability distributions for each model parameter. These are then used to define best parameter estimates, taken to be the medians of these distributions, and 99% confidence ranges defined by their upper and lower 0.5% points. Nevertheless the computed χ² turns out to be not too sensitive to whether problems occur in the decomposition (i.e. a wrong PSF) or the decomposition reliably describes the light distribution in the galaxy.

With the aim of better understanding the goodness of the fit of our models we introduce two additional measures, G₁ and G₂, which characterise the size of the residuals without reference either to the overall luminosity and size scales of a galaxy or to the expected counting noise. They are defined using the region flagged as belonging to the galaxy in the segmentation image generated by SExtractor. G₁ is the difference between the model and observed luminosities in this region as a fraction of the model luminosity, while G₂ is the ratio of the sum of the absolute values of the residuals in all the pixels to the model luminosity. These two statistically meaningful parameters provide, by their definition, a different information: G₁ quantifies whether the modelling led to an over or an underestimate of the total galaxy luminosity while G₂ quantifies by how much our model differs from the real galaxy. Working from the individual pixels (ij), by definition the total counts in the science image (∑ O_ij) and in the model image (∑ M_ij) are due to the light of the galaxy plus a uniform sky (∑ S_ij), which is the same in the two images. The total counts in the residual image simply reflect the difference between the luminosity of the observed and modelled galaxy (∑ D_ij =  ∑ O_ij −  ∑ M_ij). Our definitions can consequently be formulated as: (3)Figure 8 shows the distribution of galaxies with respect to these quantities as a function of galaxy absolute magnitude. Results are presented for fits to the i band data both for floating n and for n fixed to the de Vaucouleurs value. Results for the r band are very similar. The G₁ parameter is narrowly distributed around zero, with a slight bias towards positive values. Thus Gim2D underestimates slightly the luminosities of these large galaxies, but typically by only a couple of percent. Misestimates by more than 5% are very rare. The distribution of G₂ peaks at 0.1 and is skewed with a longer tail towards higher values. There is a tendency for deviations from the models to be larger for intrinsically fainter galaxies, particularly below about M_i =  −20. Residuals are only slightly reduced by the extra freedom involved in allowing n to vary because, as can be seen from the examples in Fig. 2, the dominant residuals are often due to non-symmetric features such as spiral arms or dust lanes. Nonetheless small residuals do not necessarily mean good fits: there may be degenerate solutions with other profile shapes that fit the data equally well. However, without extra information (i.e. kinematics, stellar populations) it is impossible to know which decomposition is correct. We are encouraged that G₂ is less than 15% for about three quarters of our galaxies and almost never rises as high as 25%. We therefore believe that our model represents the images of the majority of galaxies adequately for our purpose, and that derived parameters can be used meaningfully to characterise physical properties of the galaxies themselves.

Fig. 7 Similar to Fig. 5 but for the 1742 galaxies (in the photometric sample) modelled in both bands by fitting a Sérsic profile to the bulge and an exponential to the disc.

Fig. 8 Goodness of fit measures for our two-component modelling of the i-band images. The two left-hand panels give results for the 1469 galaxies modelled successfully using a de Vaucouleurs profile for the bulge and an exponential profile for the disc. The two right-hand panels give results for the 1515 galaxies modelled successfully when we use the more general Sérsic profile for the bulge. The two upper panels show results for , the total light in the residual image in units of the total light of the model. The two lower panels give results for , the sum of the absolute values of the individual pixel deviations from the model again in units of the total model luminosity. G1 and G2 are plotted against absolute magnitude in the main part of each panel while a histogram of their marginal distribution is given in the inset.

5. Results

5.1. Comparison of Sérsic index distributions

Fig. 9 Distribution of the Sérsic index n. The upper-left panel shows the distribution for the total sample of 1766 galaxies successfully modelled in the i band. In the other panels the Sérsic index distribution is separated according to morphological type.

Previous studies have shown that the exponential fits and Sérsic formula can fit the light profiles of many nearby ellipticals and bulges extremely well. The first paper to suggest that exponential fits are preferable to r^1/4 fits for bulges in spiral galaxies was Andredakis & Sanders (1994), while Andredakis et al. (1995) was the first to use Sérsic’s model to describe bulges in disc galaxies. This work was refined in Balcells et al. (2003) using HST NICMOS data, and it was found that de Vaucouleurs bulges are also not so frequent in early-type disc galaxies in the near-infrared, in addition to late-type disc galaxies (e.g., de Jong 1996; Graham 2001; MacArthur et al. 2003). While massive systems usually require large values of n, similar to the de Vaucouleurs value, less massive ellipticals and bulges, particularly the bulges of late-type spirals usually demand smaller n values and indeed can even be well fit by an exponential law with n = 1. In particular, various studies which include early-type bulges (Andredakis et al. 1995; Khosroshahi et al. 2000; Graham 2001; Möllenhoff & Heidt 2001) find a continuous distribution of Sérsic shape indices n that scales with bulges luminosity from n < 1 to n > 4. For our JPG sample when we allow the Sérsic index to assume any value between 0.2 and 10 we obtain the distributions shown in Fig. 9. The upper-left histograms give results for the 1766 galaxies in the i band (out of the the total sample of 1862) which are successfully fit by Gim2D. These distributions show a peak at n = 4, confirming that the r^1/4 law provides an acceptable fit to the bulge component of a large fraction of the galaxies in a magnitude-limited sample. Nonetheless the gaussian shape clearly indicates that values higher than 4 are needed for a large fraction of objects. In the other panels of Fig. 9 we split the sample by morphological type and it becomes evident that most E/S0 galaxies are well described by a wide range of n value. The same is true for the bulges of most early-type spirals. For later-type spirals there is a clear shift to lower values of n. In particular, for Sd and irregular objects n ~ 1 is preferred, confirming the earlier studies referred to above. The distributions in Fig. 9 are almost independent of the band in which the decomposition is carried out, confirming the robustness of the results.

5.2. Disc and bulge luminosity

In this section we derive the fraction of the luminosity density in the local universe in bulges and in discs for the Sloan i and r bands, starting from our complete sample of r-selected galaxies with r < 15.9 (after correction for Galactic extinction). From our total sample of 1862 galaxies we here consider only the 1588 objects for which spectroscopic data are available and it is therefore possible to measure the absolute magnitudes needed by our estimation procedure. The absolute magnitudes used in this paper are k-corrected using the code of Blanton et al. (2003b), . In the SDSS main galaxy sample as a whole the median redshift is near z = 0.1, so Blanton chose to express results in the SDSS filter system shifted by 0.1. For our sample the median redshift is about 0.05. Nevertheless for consistency with other SDSS work (in particular, with the luminosity functions we use below) we follow Blanton’s convention and k-correct to z = 0.1. We denote absolute magnitudes in this system as and to distinguish them from those in the unshifted system. Since our primary results concern the ratios of luminosities in different components, this choice has no effect on our analysis. Finally, it is relevant to stress that the Petrosian magnitudes miss the galaxy’s flux outside the Petrosian aperture. The size of this flux deficit varies monotonically with the shape of a galaxy’s light profile. For this reason the absolute magnitude used in this paper are obtained after applying the correction to a total flux proposed in Graham et al. (2005).

There are 23 objects in the spectroscopic sample for which no morphological type was assigned by the JPG astronomers. We exclude these from further consideration here, leaving 1565 objects with redshift and a well defined “by eye” morphological type. From this sample we also excluded all galaxies for which Gim2D failed in the modelling, and in addition we removed three objects for which the k-correction is not reliable due to bad photometric data in the bluest and reddest bands. This reduced the number of objects used to estimate the luminosity densities in bulges and discs to the numbers in the last two columns of Table 2. To summarise, of the 1834 galaxies in our photometric sample for which the morphology is known, 1742 are successfully modelled by the code in both bands when using the Sérsic profile and 1638 if we adopt the de Vaucouleurs law. Galaxies in the spectroscopic sample are used to calculate the fraction of light in the local universe in discs and bulges. Out of the 1565 galaxies in our sample for which the redshift and the morphology are known, 1496 are successfully modelled by the code in both bands when using the Sérsic profile and 1407 if we adopt the de Vaucouleurs law.

Our strategy for computing the fraction of the local luminosity density which is in bulges and discs is as follows. We separate the galaxies in our sample into 11 bins according to their absolute luminosities in the r or i bands. The nine brightest bins each contain about 10% of the sample while the two faintest bins contain about 5% (we made this choice in order to get better luminosity coverage for faint galaxies). For each bin we then use our decompositions to estimate the fraction of the total light coming from discs for galaxies at that absolute magnitude. Assuming this fraction to be appropriate for all galaxies of similar intrinsic brightness, we combine it with luminosity functions determined from much larger SDSS samples by Blanton et al. (2003a) to obtain the fraction of the local luminosity density which is in discs. The complementary fraction is then the amount in bulges.

While this appears straightforward, a serious complication arises from the fact that 2D fitting codes like Gim2D tend to fit radial variations in axial ratio or position angle in ellipsoidal galaxies by assigning a fraction of their light to a disc, when in fact none is present. This systematic is well known and is commented on in Simard et al. (2002). We can demonstrate its presence in our purely luminosity-selected data by examining the distribution of disc apparent axial ratio returned by Gim2D. The distribution of b/a is expected to be uniform on  [0,1]  for randomly oriented thin discs. Figure 10 shows the distributions we actually obtain for the disc components in our sample, split into 10 equal bins by absolute total i-magnitude. While for faint galaxies these distributions are indeed consistent with being flat, in the brighter bins there is clearly a strong bias towards high b/a. Among the brightest galaxies almost no disc components are found with b/a < 0.5. Figure 11 shows that these bright bins are dominated by early-type galaxies according to the visual classifications of the JPG. The absence of galaxies with small b/a values demonstrates that only a small number of these systems actually have significant thin discs, despite the fact that Gim2D assigns most of the bright galaxies B/T ratios substantially smaller than unity (see Fig. 12).

Fig. 10 Distribution of apparent axial ratio for the disc components of the 1469 galaxies in the spectroscopic sample modelled successfully in the i band by a de Vaucouleurs plus an exponential. Galaxies are split by absolute magnitude into 10 bins containing approximately equal numbers of objects.

Fig. 11 As Fig. 10 but showing the distribution over morphological type.

Fig. 12 As Fig. 10 but showing the distribution over bulge-to-disc ratio.

In order to circumvent this problem in the following analysis we use only the galaxies in each absolute magnitude bin which have b/a < 0.5. This ensures that the great majority of these must be true discs since ellipticals with apparent axial ratios smaller than 0.5 are very rare. For a random orientation the total number of true discs expected in the bin is just twice the number with b/a < 0.5. The total light in discs in the bin is, however, more than twice the light in the discs with b/a < 0.5, since dust extinction is significantly stronger in edge-on than in face-on discs. This must be corrected if we wish to obtain an unbiased estimate of the luminosity density in discs. As we effectively use only galaxies with b/a < 0.5, we restrict our original sample by almost half. However using this procedure ensures that we properly eliminate fake disc detections, and are left with a population of galaxies from which the true distribution of light in discs and bulges can be estimated.

In practice, our procedure works as follows. For each absolute magnitude bin k we estimate the fraction of the light in the disc component as (4)where L_disc,k is the total luminosity of the discs of the galaxies in the bin and L_tot,k is the total luminosity from all components of these same galaxies. Since we assume that we can rely on our decomposition only for edge-on systems, we split the numerator into two parts (5)where L_{b/a < 0.5,disc,k} is the luminosity due to “edge-on” discs with b/a < 0.5 and is obtained directly from our decompositions. We estimate the luminosity L_{b/a > 0.5,disc,k} in “face-on” discs by assuming that true discs are randomly oriented and that their internal extinction A_λ depends on inclination according to the standard prescription (6)Here A_λ is the correction to exactly face-on orientation. Assuming this formula, it is a simple matter to relate the total luminosity density in discs to that in discs with b/a < 0.5. We find (7)hence (8)The bulge to total ratio B/T is one of the structural parameters returned by Gim2D for each galaxy. In combination with the corrected total absolute Petrosian magnitude it allows us to estimate the disc luminosity of each object contributing to L_{b/a < 0.5,disc,k}. We obtain L_tot,k simply by summing the Petrosian luminosities of all galaxies in the bin regardless of their b/a. We take values for γ_λ from the work of Tully et al. (1998). Since the numerical coefficient on the right hand side of equation 8 is only slightly smaller at i than at r band, we assume the value 2.56 for both photometric bands.

The results of these calculations are shown in Fig. 13 for the galaxies modelled successfully in the i and r bands. In each figure we show results separately for decompositions which force a de Vaucouleurs bulge and for decompositions in which we allow 0.2 < n < 10. The pattern is very similar in all cases. The average light fraction in discs varies smoothly from about 10% in the brightest galaxies to almost 100% in faint galaxies. For the brighter bins these fractions are much smaller than the values obtained from a direct naive average of the B/T histograms of Fig. 12 because of the systematic effect we have just been discussing.

The error bars on the points in Fig. 13 are important because they determine the precision of our final results. We neglect the formal errors on B/T returned by Gim2D because these are quite small, typically  ± 7%, and are well below the systematic error discussed above due to the assignment of isophote twists or axial ratio changes to spurious thin discs. We assume, however, that this systematic can be neglected for systems with b/a < 0.5. The uncertainty in our estimate of the disc light fraction in each bin is then dominated by sampling. As a statistical model for the population of a particular bin we assume that a randomly chosen galaxy has a detectable edge-on disc (b/a < 0.5) with probability p_eo where a priori we have 0 < p_eo < 0.5. We also assume that the B/T values of these edge-on systems are drawn at random from some unknown distribution with population mean and variance which we estimate using the sample mean and variance,  ⟨ B/T ⟩  and Var(B/T) respectively. If the bin contains N_t galaxies of which N_eo have discs with b/a < 0.5, then the maximum likelihood estimate of p_eo is provided N_eo/N_t < 0.5 (it is equal to 0.5 otherwise). To approximate the variance of we use the standard binomial formula , even though this is formally incorrect for . Our estimate of the mean light fraction in edge-on discs is then and we calculate the variance in this estimate as . These give the final results plotted when multiplied by the correction factor of Eq. (7) which accounts for the light in face-on discs. There is undoubtedly a systematic uncertainty associated with this last step, but our analysis shows that this should be small, and we ignore it here.

Fig. 13 The left-hand panel shows the luminosity fraction in discs for galaxies in the Sloan i band. The galaxies are modelled with two different sets of parametric functions: a de Vaucoulours profile for the bulge and an exponential for the disc (blue dots) or a Sérsic profile for the bulge plus an exponential for the disc (cyan triangles). The right-hand panel shows the same results but in the r band.

We can now average the disc light fractions of Fig. 13 over the galaxy population as whole in order to obtain the fractions of the total luminosity density in the local universe coming from discs and from bulges. In doing that we make the assumption that the fraction of light in discs in low luminosity galaxies (hence at lower redshifts of our sample) and the fraction of light in discs for bright galaxies (hence at the bright redshift end of our sample) is the same across the redshift range of the full sample. This assumption seems reasonable as the evolution of the fraction of light in discs is expected to be small since z ~ 0.1. The contribution of each of our absolute magnitude bins to the total luminosity density Φ_tot,k can be obtained by integrating the appropriate Blanton et al. (2003a) luminosity function across the bin. The final result for the fraction of the local luminosity density in discs is then, (9)Being our sample a complete hence representative sample of the galaxy population, and after performing a Kolmogorov-Smirnov test, we are reassured in using Blanton et al. (2003a) luminosity function parametrisation in the same bands. We prefer to use the values provided in Blanton et al. (2003a) because of their much better statistical precision. The final result we obtain for the total fraction of the light coming from discs in the local universe is (54    ±    2)%, with no detected dependence on observing band or decomposition parametrisation. The details are in Table 3. The error bars in this table are calculated directly from those in Fig. 13 assuming the uncertainties in the different absolute magnitude bins to be independent. Uncertainties coming from the luminosity function itself are negligible in comparison.

A slight variation of this analysis allows us to calculate a second interesting quantity: the fraction of galaxies in each of our absolute magnitude bins which contain no detectable thin disc and so may be considered “pure” bulge systems. Our hypothesis here is that Gim2D will detect any significant disc if it is sufficiently inclined to the line-of-sight that b/a < 0.5. Exactly one half of all disc galaxies should be at least this inclined. We can thus estimate the number of effectively discless galaxies in each bin by subtracting twice the number of objects with b/a estimates below 0.5 from the total number of galaxies in the bin. The result of this exercise is shown in Fig. 14 in similar format to Fig. 13. Again the results are very similar in the two pass-bands and for our two assumptions about bulge profiles. For the brightest bins we find that the great majority of galaxies (among  ~90% and 60%) are effectively discless, while in the faintest bins our statistics are consistent with at most a small fraction (less than 30%) of “pure bulge” systems. The fraction of “pure bulges” varies smoothly with absolute magnitude between these two extremes. The error bars reflect the binomial uncertainty in our estimate of the fraction of the bin population with detectable edge-on discs.

As before we can combine the results of Fig. 14 with the luminosity functions of Blanton et al. (2003a) to estimate the fraction of the luminosity density of the local Universe which is contributed by effectively discless systems. The result, given in detail in Table 4, is about (32 ± 2)% with no significant dependence on pass-band or bulge fitting function. Thus the breakdown of stellar luminosity in the local Universe is apparently 54% in discs, 32% in “pure bulge” systems with no photometrically detectable (by Gim2d) discs and at most 14% in the bulges of galaxies with detectable discs.

Fig. 14 The left-hand panel shows the fraction of pure bulge galaxies as a function of the i band absolute magnitude. Different parametric functions are considered: de Vaucouleurs plus exponential (blue dots) and Sérsic plus exponential (cyan triangles). The right hand panel shows the same results but in the r band.

6. Comparison with previous work

The knowledge of the amount of baryons present in different states and forms in the universe is fundamental to place constraints on large scale structure formation theories. This was the main motivation behind milestone studies in the late eighties: Schechter & Dressler (1987) were among the first to explore the fraction of the luminosity density produced by the spheroid or disc component of galaxies, since stars are the most prominent location of baryons. Using nearly complete magnitude-limited samples and an admittedly crude bulge/discs decompositions on photographic plate images they estimated that the relative contributions of bulge and disc to the mean mass density of the universe were very nearly equal. Comparing to more recent results, our results are in perfect agreement with Fukugita et al. (1998), who estimated that 55% of the luminosity density in the local universe is included in galactic disc components. A subsample of our sample is used by Oohama et al. (2009) to perform a complementary analysis. They carried out the bulge-disc decomposition using the growth curve fitting method, instead of the more commonly used surface brightness profile extracted from two dimensional surface brightness distribution, and they reached our same conclusions for what concerns the relative contributions to disc and spheroids (0.53:0.47). Benson et al. (2007) using a SDSS sample 1.4 mag fainter then our, already compared to an early version of our work and found an excellent agreement with our results. Various studies (Benson et al. 2002; Bell et al. 2003; Driver et al. 2007; Gadotti 2009 and Weinzirl et al. 2009) using different surveys, selections and tools estimated how stellar mass in the local universe is distributed among galactic components. While these studies present apparently discrepant results, the difference can be mainly assigned to the different selection criteria and different classifications used. In particular we stress that by construction, the use of mass-selected volume-limited samples (i.e. Gadotti 2009) selects against low mass galaxies thus the results obtained cannot be considered as representative of the whole population but biased in favour of more massive objects, that we know to be early type galaxies. It is therefore not surprising that some of these studies found a much higher stellar mass densities in the galactic spheroidal component. The parametric functions used in the models to fit the galaxy light distribution, and in particular the implementation of a bar component, might participate to the observed discrepancy. We pointed out that Gim2D does not allow to add an extra component to model bar-like structures. We have shown in Sect. 4.1 that Gim2D does not create a no realistically elongated spheroid component (e > 0.8) when a bar is present and therefore does not overestimate the bulge fraction. Nonetheless, it has been pointed out (Laurikainen et al. 2005; Gadotti 2008) that even when there is no bar in the model (as in our case) when fitting a barred galaxy one can get an acceptable (though wrong) fit. Gadotti (2008) proposed a formalism, to correct the luminosity density in the local Universe, which takes into account the light in stars which resides in bars. As a first approximation to the estimate of the total light in the local Universe coming from the bar component, we use the Gadotti (2008) results on the biases in the estimation of the bulge and disc luminosity fraction due to the non-inclusion of bars in the photometric models. We obtain that the luminosity density in bulges, discs and bars as extracted from disc galaxies is  ≈ 9%, 51% and 8% respectively, while the luminosity density in pure bulges (discless galaxies) remains unchanged at  ~ 32%. These results would be only slightly modify if we had included an additional correction, again proposed by Gadotti (2008), to correct for low spatial resolution. The aforementioned correction are based on the fact that: if bars are not modelled (i.e. there is no parametric function in the models to model the bar light) the bulge-to-total luminosity ratio is overestimated on average by  ~50% while the disc-to-total luminosity ratio by  ~10%; the fraction of barred galaxies, considering only disc galaxies, is  ≈ 70%. Our values are globally in agreement with Driver et al. (2007) who find that the stellar mass content in discs is 58 ± 6%.

7. Conclusions

We have used the two-dimensional photometric fitting programme Gim2D of Simard et al. (2002) on the images of an apparent magnitude-limited sample of 1862 galaxies selected from the SDSS (median redshift of 0.05), with r ≤ 15.9, and visually classified by Nakamura et al. (2003) and Fukugita et al. (2007). In almost all cases the code returns a well-defined decomposition of the galaxy into disc and bulge components with parameters which have small formal error bars and vary little either between the r and i band images analysed here or between decompositions in which the Sérsic index of the bulge component varies or is fixed to the de Vaucouleurs value n = 4. The total amount of light in the differences between the observed image and the best fit model is typically only about 10% of the galaxy luminosity. Despite this apparent success, we show that in most cases for intrinsically bright galaxies the “disc” component does not represent a true thin axisymmetric disc, since the sample distribution of axial ratios b/a deviates strongly from the uniform distribution expected for a population of such discs  −  near edge-on discs are grossly under-represented. This problem was already noted by Simard et al. (2002). Apparently the code is using the degrees of freedom provided by the assumed disc component to fit radial changes in isophote shape or position angle within galaxies which have no real thin disc.

We introduce a robust method to correct for this systematic problem by concentrating on galaxies for which Gim2D finds a disc axial ratio b/a < 0.5. We argue that ellipsoidal stellar systems with such extreme apparent axial ratios are quite rare so that the component isolated by the the decomposition programme is likely to be a highly inclined disc. Since the selection criteria for our sample depend at most weakly on inclination, we can use the assumption that our galaxy sample is randomly oriented to correct from the sample with highly inclined discs to the sample as a whole.

Using these multi-parameter fits, we estimate average disc light fractions for galaxies as a function of absolute magnitude. These range from about 10% for the brightest galaxies to almost 100% for the faintest ones. At each absolute magnitude we also estimate the fraction of “pure bulge” galaxies, defined as galaxies for which Gim2D would detect no disc with b/a < 0.5 even if the orientation were such that b/a < 0.5 would be expected for any true axisymmetric thin disc. We find that most of the galaxies in the brighter absolute magnitude bins are “pure bulge” by this definition, but that this fraction decreases steadily for fainter systems. We do not detect any population of “pure bulges” in our faintest bin. These numbers differ substantially from those inferred naively from the B/T distributions measured directly by Gim2D: for example, in the brightest absolute magnitude bin about 85% of galaxies are “pure bulge” and discs only contribute  ~10% of the total light, yet the median value of B/T returned by the code is 0.55 and very few galaxies are fitted with B/T > 0.8.

By combining these results for the absolute magnitude dependence of the mean light fractions in discs or in “pure bulges” with the luminosity functions of Blanton et al. (2003a), we have been able to estimate the fractions of all galaxy light in the local Universe coming from discs and from “pure bulges”. The results depend little on the waveband used or on the bulge luminosity profile we adopt. We find that 54 ± 2% of the local luminosity density is contributed by stars in discs and 32 ± 2% by stars in “pure bulges”. The remaining 14 ± 2% comes in half from the light in the bulges and from the other half from light in bars of systems with detectable discs. The mean bulge-to-total ratio of the latter systems is thus 14/68 ~ 0.2, substantially smaller than typical B/T values in the histograms of Fig. 12.

The fraction of light in bulges or discs is a quantity resulting from the interplay of processes and physical mechanisms which drive the formation and evolution of the galaxy population. Galaxy formation and evolution models must be able to reproduce these values and its evolution with redshift. There are two main mechanisms called upon to explain the formation of discs and spheroids: gas cooling and spinning is expected to form a disc, while mergers or disc instabilities would form spheroids (e.g. Benson et al. 2007). Semi-analytic models of galaxy formation, making assumptions on these physical mechanisms and introducing ad hoc prescriptions (i.e. Baugh et al. 2005; Bower et al. 2006), have proved succesful to match the SDSS luminosity functions of discs and spheroids. These mechanisms have been demonstrated to be very efficient to form discs and spheroids in recent high resolution simulations (e.g. Dekel et al. 2009). However, these models have yet to be applied to the study of galaxy morphology in a cosmological context, e.g. to reproduce the fraction of light in bulges and discs per unit volume, and its evolution with redshift, as one of the most fundamental galaxy properties. One of the main limitations of comparing simulations and observations of structural parameter in galaxies has been that the methods used to extract these parameters were different between models and observations. Only recently Governato et al. (2009, 2010) analysed simulated data to obtain “observables” more directly comparable to observational results. In doing that they used the Monte Carlo radiation transfer code SUNRISE (Jonsson 2006) to generate artificial optical images and spectral energy distributions of the outputs of N-body smoothed particle hydrodynamics (SPH) codes. The 2D light profiles are then fitted as the superposition of a Sérsic and an exponential component as done for real images. By resolving the inhomogeneous interstellar medium, resulting in strong outflows from supernovae explosions which remove low angular momentum gas, Governato et al. (2010) succeed for the fist time to form bulgeless galaxies, with an almost perfect exponential profile and a featureless disc, in a ΛCDM framework. Scannapieco et al. (2010) went even further and tried to reproduce the relative importance of the main stellar components, at redshift zero. They applied photometric decomposition methods, instead of the dynamical decomposition ones commonly used in simulations analysis. This work clearly shows that photometric decompositions match the ratio of components for spiral galaxies better than kinematic decompositions, still not being able to reproduce it, and suggest that results from simulations need to be analysed following observational techniques. These theoretical works mark a turning point to which our study provides new ingredients. While it is important to identify the physical mechanisms which allow to form bulgeless galaxies, the density of these sources should also be reproduced by simulation.

From this study, we have produced a quantitative measurement of the overall distribution of light in the two basic structural components (bulges and discs) which make up galaxies. We have also indirectly estimated the fraction of light contributed by bars in the local universe. We additionally provide the relative fraction of discless and bulgeless galaxies as a function of magnitude. These results offer a robust reference to theoretical studies and to future studies tracing the evolution of the fraction of light in discs and bulges as a clue of their formation process.

Acknowledgments

L.T. acknowledge support from Université de Provence. L.T. is grateful to the anonymous referee for his/her careful reading of the manuscript and suggestions. We thank Luc Simard for help in using his excellent two-dimensional fitting package Gim2D. L.T. wants to thank Alistar Graham for his careful reading of the manuscript and useful comments and Carlos Lopez San-Juan for helping with simulations. This work has been partially supported by a grant awarded for “EFIGI project” (grant # 45500) from the French research ministry. Funding for the creation and distribution of the SDSS Archive has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Aeronautics and Space Administration, the National Science Foundation, the US Department of Energy, the Japanese Monbukagakusho, and the Max Planck Society. The SDSS Web site is http://www.sdss.org/. The SDSS is managed by the Astrophysical Research Consortium (ARC) for the Participating Institutions. The Participating Institutions are The University of Chicago, Fermilab, the Institute for Advanced Study, the Japan Participation Group, The Johns Hopkins University, the Korean Scientist Group, Los Alamos National Laboratory, the Max-Planck-Institute for Astronomy (MPIA), the Max-Planck-Institute for Astrophysics (MPA), New Mexico State University, University of Pittsburgh, University of Portsmouth, Princeton University, the United States Naval Observatory, and the University of Washington.

References

Appendix A: Influence of bar-like structures on the measurement of B/T

In this appendix we test the ability of Gim2d to recover the correct value of the bulge to disc ratio (B/T) in the presence of a bar. Since Gim2d does not implement a parametric function to model a bar-like structure in the galaxy, the bar component should be found in the residual image. This is what we observe looking at the residual image of the fit but we additionally want check if the Gadotti (2008) correction we applied in Sect. 6 is justified. We look if part of the bar light might be erroneously attributed to the bulge or to the disc component, affecting the measurement of B/T.

Fig. A.1 Correlation between the input B/T for 100 model barred galaxies and the best-fitting B/T recovered using Gim2d. The black dashed line is the 1:1 correlation while the red line is the 1:1 line shifted according to the measured offset of 0.237.

Fig. A.2 Distribution of the offsets between the recovered and the input B/T ratio.

We create 100 simulated barred galaxies for which the structural parameters, and in particular the value of B/T, are known. Model barred galaxies are constructed assuming a Sérsic profile for the bulge, an exponential profile for the disc and a Sérsic profile, with initially a relatively flat inner and a steep outer profile (n < 1) and a boxy shape, for the bar structure. We adopted parameter values measured on galaxies and in particular we fixed the total magnitude to 15. The physical properties of bars have been taken from Gadotti (2010), the bar fraction has been fixed at 10% and the B/T ratio can assume values in the range [0–0.8], allowing for the fixed 10% of the light in the bar-like structure and for at least 10% of the light in the disc component. Poisson noise is added to the model galaxy after its image is convolved with an analytic PSF corresponding to a typical value of the seeing. The model galaxy is finally added to a real SDSS sky.

In a second step, model barred galaxies are analysed with Gim2d using exactly the same procedure used for real galaxies (i.e. a Sérsic parametric function for the central component and an exponential profile for the disc, without a parametric function to model bars). As expected, we observe a bar-like structure in the residual maps. We present the example of 100 simulations of one specific real galaxy with the comparison of the input and recovered B/T in Fig. A.1. Imposing a slope of one, the dispersion of the data is 0.01, with an offset of 0.237, as shown in Fig. A.2. Running an equivalent test, but on another set of 100 simulated galaxies for which we impose a larger bar size, we observe that the dimension of the bar has an impact on the offset value. The smallest the bar, the largest the offset. We are therefore reassured that the correction proposed by Gadotti (2008), that we applied to our results, is justified.

All Tables

All Figures

Fig. 1 The distributions of galaxies in morphological type (upper panel) and in extinction-corrected Petrosian apparent magnitude in the r band (second panel) are shown for the photometric and spectroscopic samples, empty and dashed-line respectively. The redshift (third panel) and corresponding r band absolute magnitude (lower panel) distributions are shown for the galaxies in the spectroscopic sample.

In the text

	Fig. 2 Examples of the science image, the mask, the model and the residual images for our two fits (de Vaucouleurs plus exponential or Sérsic plus exponential) in the r band (from left to right) for galaxies from our seven morphological classes: E, S0, Sa, Sb, Sc, Sd, Irr (from top to bottom).
In the text

Fig. 3 The bulge fraction (upper-left), the disc scalelength (lower-left), the effective radius of the bulge (upper-right) and the half-light radius of the galaxy (lower-right), obtained using different parametric functions to perform the decomposition (de Vaucouleurs plus exponential on the x-axis and Sérsic plus exponential on the y-axis) are plotted against each other for the 1680 galaxies (in the photometric sample) successfully modelled by the two fitting functions in the i band.

In the text

	Fig. 4 The bulge fraction (upper-left), the disc scalelength (lower-left), the effective radius of the bulge (upper-right) and the half-light radius of the galaxy (lower-right) in the i and r bands are plotted against each other for the 1638 galaxies (in the photometric sample) modelled in both bands by fitting a de Vaucouleurs profile to the bulge and an exponential to the disc.
In the text

	Fig. 5 The ellipticity of the bulge (upper-left), the inclination angle of the disc (lower-left), the position angle of the bulge (upper-right) and the position angle of the disc (lower-right) in the i and r bands are plotted against each other for the 1638 galaxies (in the photometric sample) modelled in both bands by fitting a de Vaucouleurs profile to the bulge and an exponential to the disc.
In the text

	Fig. 6 The same as Fig. 4 but for the 1742 galaxies (in the photometric sample) modelled in both bands by fitting a Sérsic profile to the bulge and an exponential to the disc.
In the text

	Fig. 7 Similar to Fig. 5 but for the 1742 galaxies (in the photometric sample) modelled in both bands by fitting a Sérsic profile to the bulge and an exponential to the disc.
In the text

Fig. 8 Goodness of fit measures for our two-component modelling of the i-band images. The two left-hand panels give results for the 1469 galaxies modelled successfully using a de Vaucouleurs profile for the bulge and an exponential profile for the disc. The two right-hand panels give results for the 1515 galaxies modelled successfully when we use the more general Sérsic profile for the bulge. The two upper panels show results for , the total light in the residual image in units of the total light of the model. The two lower panels give results for , the sum of the absolute values of the individual pixel deviations from the model again in units of the total model luminosity. G1 and G2 are plotted against absolute magnitude in the main part of each panel while a histogram of their marginal distribution is given in the inset.

In the text

	Fig. 9 Distribution of the Sérsic index n. The upper-left panel shows the distribution for the total sample of 1766 galaxies successfully modelled in the i band. In the other panels the Sérsic index distribution is separated according to morphological type.
In the text

	Fig. 10 Distribution of apparent axial ratio for the disc components of the 1469 galaxies in the spectroscopic sample modelled successfully in the i band by a de Vaucouleurs plus an exponential. Galaxies are split by absolute magnitude into 10 bins containing approximately equal numbers of objects.
In the text

	Fig. 11 As Fig. 10 but showing the distribution over morphological type.
In the text

	Fig. 12 As Fig. 10 but showing the distribution over bulge-to-disc ratio.
In the text

Fig. 13 The left-hand panel shows the luminosity fraction in discs for galaxies in the Sloan i band. The galaxies are modelled with two different sets of parametric functions: a de Vaucoulours profile for the bulge and an exponential for the disc (blue dots) or a Sérsic profile for the bulge plus an exponential for the disc (cyan triangles). The right-hand panel shows the same results but in the r band.

In the text

	Fig. 14 The left-hand panel shows the fraction of pure bulge galaxies as a function of the i band absolute magnitude. Different parametric functions are considered: de Vaucouleurs plus exponential (blue dots) and Sérsic plus exponential (cyan triangles). The right hand panel shows the same results but in the r band.
In the text

	Fig. A.1 Correlation between the input B/T for 100 model barred galaxies and the best-fitting B/T recovered using Gim2d. The black dashed line is the 1:1 correlation while the red line is the 1:1 line shifted according to the measured offset of 0.237.
In the text

	Fig. A.2 Distribution of the offsets between the recovered and the input B/T ratio.
In the text