THE SPITZER SURVEY OF STELLAR STRUCTURE IN GALAXIES (S⁴G): MULTI-COMPONENT DECOMPOSITION STRATEGIES AND DATA RELEASE

Our decompositions use the GALFIT-software (Peng et al. 2002, 2010), which has become the de facto standard for detailed two-dimensional structural decompositions. It relies on parametric fitting, using the Levenberg–Marquadt algorithm to minimize the weighted residual ${\chi }^{2}{}_{\nu }$ between observed (OBS) and model (MODEL) images,

$$\begin{eqnarray}&&{\chi }^{2}{}_{\nu }=\displaystyle \frac{1}{N}\underset{x}{\displaystyle \sum }\underset{y}{\displaystyle \sum }\displaystyle \frac{{[\mathrm{OBS}(x,y)-\mathrm{MODEL}(x,y)]}^{2}}{\sigma {(x,y)}^{2}}.\end{eqnarray} \tag{ 1 }$$

The sum is taken over all the used (non-masked) image pixels, and $\sigma (x,y)$ indicates the statistical uncertainty of each pixel (sigma-image). The model image consists of a sum of model components, i.e., for bulge, disk, bar etc., convolved with the image point-spread function (PSF-image). Note that the reduced ${\chi }^{2}{}_{\nu }$ is used, with N denoting the degree of freedom, equal to the number of fitted pixels minus the number of free parameters in the fit.

GALFIT is extremely versatile in its selection of model components. Basically the user defines for each component its "radial" profile function, giving the surface brightness ${\rm{\Sigma }}(r)$ at each isophotal radial coordinate r. The isophotal coordinates are most commonly defined in terms of generalized ellipses (Athanassoula et al. 1990),

$$\begin{eqnarray}&&r(x\prime ,y\prime )={\left({| x\prime -{x}_{0}| }^{C+2}+{\left|\displaystyle \frac{y\prime -{y}_{0}}{q}\right|}^{C+2}\right)}^{\frac{1}{C+2\phantom{{}^{2}}}}.\end{eqnarray} \tag{ 2 }$$

Here ${x}_{0},{y}_{0}$ defines the center of the ellipse, q = b/a is the ratio between minor and major axis lengths. The $x\prime ,y\prime$ denote coordinates in a system aligned with the ellipse, with the major axis pointing at the position angle PA. For pure ellipses C = 0, while $C\gt 0$ indicates boxy and $C\lt 0$ disky isophotes.²⁰ For the pipeline decompositions, simple elliptical isophotes C = 0 are used for all components. Besides generalized ellipses, GALFIT provides several alternatives, such as definition of isophotal shape via azimuthal or bending modes, or via coordinate rotations, which would form a natural basis for detailed modeling of e.g., logarithmic spirals. To keep our models relatively simple (and uniform over the wide range of angular sizes and surface brightnesses spanned by the sample), we have not used these advanced GALFIT features. Keeping the models simple makes the interpretation of the observation minus model residuals more straightforward (see the NGC 1097 examples in Sheth et al. 2010).

The pipeline decompositions use five different choices for the model components/radial functions.

• 1.  
  The bulge component is described with a Sérsic profile ("sersic")

  $$\begin{eqnarray}&&{\rm{\Sigma }}(r)={{\rm{\Sigma }}}_{{\rm{e}}}\mathrm{exp}\left(-\kappa \left[{(r/{R}_{{\rm{e}}})}^{1/n}-1\right]\right),\end{eqnarray} \tag{ 3 }$$

  where ${{\rm{\Sigma }}}_{{\rm{e}}}$ is the surface brightness at the effective radius R_e (isophotal radius encompassing half of the total flux of the component). The Sérsic-index n describes the shape of the radial profile, which becomes steeper with increasing n. In particular, n = 1 corresponds to an exponential profile and n = 4 to a de Vaucouleurs profile. The factor κ is a normalization constant determined by n. In GALFIT the corresponding "sersic"-function is used, with the integrated magnitude m_bulge as a free parameter (instead of ${{\rm{\Sigma }}}_{{\rm{e}}}$).
• 2.  
  In the case of low or moderate inclination, the disk component is described with an infinitesimally thin exponential disk ("expdisk"),

  $$\begin{eqnarray}&&{\rm{\Sigma }}(r)={{\rm{\Sigma }}}_{{\rm{o}}}{q}^{-1}\mathrm{exp}(-r/{h}_{r}),\end{eqnarray} \tag{ 4 }$$

  where ${{\rm{\Sigma }}}_{{\rm{o}}}$ is the central surface brightness of the disk observed from the perpendicular direction and h_r denotes the exponential scale length. In this case the $q=\mathrm{cos}i$, where i is the disk inclination. Assuming no extinction, ${{\rm{\Sigma }}}_{{\rm{o}}}{q}^{-1}$ is the projected central surface brightness at the sky plane. The "expdisk"-function in GALFIT is used, with integrated ${m}_{\mathrm{disk}}=-2.5{\mathrm{log}}_{10}(2\pi {{\rm{\Sigma }}}_{0}\ {h}_{r}{}^{2})$ as a free parameter (instead of ${{\rm{\Sigma }}}_{0}$). Note that in cases that had more than one disk component, the inner disk was sometimes fit with a sersic or ferrer2 function, to allow the profile to drop faster than with expdisk.
• 3.  
  For a nearly edge-on disk (apparent axial ratio $q\lesssim 0.2$), the function ("edgedisk")

  $$\begin{eqnarray}&&{\rm{\Sigma }}({r}_{x},{r}_{z})={{\rm{\Sigma }}}_{{\rm{o}}}\displaystyle \frac{{r}_{x}}{{h}_{r}}{K}_{1}\left(\displaystyle \frac{{r}_{x}}{{h}_{r}}\right)\ {\mathrm{sech}}^{2}({r}_{z}/{h}_{z}),\end{eqnarray} \tag{ 5 }$$

  is adopted, where r_x and r_z are the (positive) distances along and perpendicular to the apparent major axis of the disk, and K₁ stands for a modified Bessel function. This function corresponds to the line of sight (viewing along the disk plane) integrated surface brightness of a 3D luminosity density distribution (van der Kruit & Searle 1981)

  $$\begin{eqnarray}&&L({r}_{x},{r}_{z})=\displaystyle \frac{{{\rm{\Sigma }}}_{0}}{2{h}_{r}}\mathrm{exp}(-{r}_{x}/{h}_{r})\ {\mathrm{sech}}^{2}({r}_{z}/{h}_{z}).\end{eqnarray} \tag{ 6 }$$
• 4.  
  For a bar component a modified Ferrers profile ("ferrer2") is assumed,

  $$\begin{eqnarray}{\rm{\Sigma }}(r)=\left\{\begin{array}{ll}{{\rm{\Sigma }}}_{{\rm{o}}}{\left[1-{(r/{r}_{\mathrm{out}})}^{2-\beta }\right]}^{\alpha } & r\lt {r}_{\mathrm{out}}\\ 0 & r\geqslant {r}_{\mathrm{out}}.\end{array}\right.\end{eqnarray} \tag{ 7 }$$

  Here r_out defines the outer cut of the profile, while α defines the sharpness of this cut. The parameter β defines the central slope of the profile, and ${{\rm{\Sigma }}}_{0}$ is the central surface brightness (in the plane of the sky).
• 5.  
  When the galaxy contains an unresolved central component it is fit with a PSF-convolved point source ("psf"). In this case the free parameter is the total magnitude m_psf. Typically, this component, if present, is not an active or starburst nucleus, but rather a small bulge with angular size so small that it cannot be resolved in the S⁴G images (${R}_{{\rm{e}}}\lesssim \mathrm{FWHM}=2\buildrel{\prime\prime}\over{.} 1$ of S⁴G images).

For the decomposition pipeline we chose to do three types of decompositions: (1) 1-component Sérsic-fits, (2) 2-component bulge-disk decompositions using Sérsic-bulges and exponential disks (or edge-on disk if appropriate), and (3) multi-component "final" decompositions, optionally with additional bar, disk and central components (the level of complexity of the models is discussed in more detail in Section 3). The first two types of models are made in an automatic manner, while the final models always include human judgment about what components should be included.

2.2.1. What is Needed?

2.2.2. Mask Images

2.2.3. Galaxy Centers, Sky Background, Isophotal Profiles

After the edited masks were completed, we run the galaxies through a semi-automatic IDL script which determines the galaxy centers, sky background levels and galaxy orientation parameters. The accurate galaxy center is measured with the cntrd-routine,²¹ after its approximate location is interactively defined. We also have the option to mark the center by force, in case the automatic center finding routine does not work satisfactorily even after repeated trials.

The regions used for estimation of the sky background are identified manually, by selecting several (typically 10–20) locations outside the visible galaxy, while avoiding the image edges or contaminated areas. The local sky values in these locations are obtained by taking medians of the non-masked pixels in 40 pix × 40 pix boxes. The global sky background (SKY) and its uncertainty (DSKY) are then estimated from the mean and standard deviation of these local values, respectively (see Figure 1). In Section 4, using the estimated DSKY, we show that the expected uncertainty of decomposition parameters caused by possible uncertainties in background subtraction is negligible. We also determine the average rms sky variation (RMS), by taking the median of standard deviations in different sky regions (after removing outliers by iterative 3σ clipping). We use the sky rms estimates in Section 2.2.5 for assessing the validity of theoretically calculated sigma-images.

Zoom In Zoom Out Reset image size

We calculate the isophotal profiles with a pyraf script called from IDL, using the standard IRAF ellipse algorithm (Jedrzejewski 1987). As inputs for the ellipse fitting the sky background subtracted data image and the edited mask image are used. We fix the ellipse center to the previously found galaxy center and use a logarithmic increment of 0.02 between isophote levels. As often happens with IRAF ellipse, the fit does not necessarily converge over the whole galaxy area: we have an option to re-try the fit with different starting locations until a successful fit is obtained over the whole galaxy region (see Figure 2). From the isophotal profiles, we choose a semimajor axis range from which outer orientations (${(b/a)}_{\mathrm{outer}}=1-{\epsilon }_{\mathrm{outer}},{\mathrm{PA}}_{\mathrm{outer}}$) are estimated. Also, a rough estimate of the galaxy outer radius, R_gal, is made to define the image region used in the GALFIT decomposition.

Zoom In Zoom Out Reset image size

Figures 1 and 2 give examples of typical plots produced during these preparatory steps, illustrating the sky background fitting and the elliptical isophote profiles. The estimated ${\epsilon }_{\mathrm{outer}}$ and PA_outer are marked. In all of our decompositions, we fix the orientation of the disk component to these outer values²² and interpret them to represent the galaxy viewing inclination. Therefore, extra care is taken to estimate the orientations reliably. For example, the corresponding inclination ${i}_{\mathrm{disk}}={\mathrm{cos}}^{-1}{(b/a)}_{\mathrm{outer}}$ is visually checked by de-projecting the galaxy images to face-on. Figure 3 shows an example of such a de-projection, also comparing the estimated inclinations with those calculated from axial ratios given in the HyperLeda database. Typically, our ${(b/a)}_{\mathrm{outer}}$ and PA_outer are determined at much lower surface brightness levels than those in HyperLeda (which are mainly from RC3 de Vaucouleurs et al. 1991). P4 values are thus less affected by bulges, bars, or prominent spirals, and should reflect better the orientation of the underlying extended disk, which appears with circular outer isophotes in face-on projection.²³ In some cases, S⁴G images are so deep that the outermost isophotes are dominated by an outer stellar halo rather than the disk. Good examples are NGC 681, NGC 1055, and NGC 4594. Possible misinterpretations of the outer isophotes were avoided by visually inspecting all the images: when the disk (identified with spiral arms, rings, and lenses) is clearly more inclined than suggested by the outermost isophotes of the image, the isophotes in the disk region were used for the estimate of galaxy orientation. For nearly face-on galaxies, the possible stellar halos are more difficult to distinguish, but in these cases the involved error in the orientation is less important. The final P4 axial ratios and position angles, center locations, and sky background values are listed in Table 1. For each galaxy we also include a flag indicating the inclination uncertainty: "ok" indicates that outer isophote axial ratio should give a reliable estimate of i_disk, "u" indicates that the inclination is uncertain, while "z" indicates that the galaxy is close to edge-on (the axial ratio is not used for an inclination estimate).

Zoom In Zoom Out Reset image size

Table 1.  Pipeline 4 Parameters: Galaxy Center, Outer Orientation, and Sky Background

IDE	xc	yc	PA ± dPA	ELL ± dELL	RANGE	FLAG	SKY	DSKY	RMS	${\sigma }_{\mathrm{conf}}$
ESO 011-005	780.84	299.02	42.7 ± 0.1	0.747 ± 0.002	13–16	z	0.0125	0.0039	0.0109	0.0092
ESO 012-010	760.68	294.95	146.2 ± 1.1	0.542 ± 0.012	75–90	ok	−0.0022	0.0025	0.0109	0.0096
ESO 012-014	782.57	459.99	31.0 ± 5.2	0.580 ± 0.040	67–75	u	0.0080	0.0033	0.0106	0.0091
ESO 013-016	478.49	282.78	−14.3 ± 3.1	0.343 ± 0.031	52–82	ok	0.0041	0.0027	0.0101	0.0084
ESO 015-001	291.43	294.80	125.7 ± 2.7	0.586 ± 0.018	56–75	u	0.0050	0.0018	0.0102	0.0086
ESO 026-001	795.50	495.47	19.3 ± 21.1	0.060 ± 0.028	52–60	ok	0.0125	0.0024	0.0105	0.0091
ESO 027-001	776.34	294.63	12.5 ± 10.1	0.216 ± 0.061	112–127	u	0.0109	0.0032	0.0112	0.0097
...
UGC 12791	295.88	289.63	82.6 ± 1.4	0.709 ± 0.037	45–60	ok	0.0406	0.0033	0.0110	0.0094
UGC 12843	290.20	287.90	17.5 ± 3.9	0.553 ± 0.058	52–67	u	0.0418	0.0037	0.0105	0.0086
UGC 12846	571.15	851.20	−4.7 ± 14.0	0.127 ± 0.031	60–67	u	0.0446	0.0007	0.0019	0.0000
UGC 12856	291.07	296.89	16.9 ± 1.5	0.615 ± 0.049	48–75	u	0.0401	0.0025	0.0111	0.0099
UGC 12857	284.91	280.12	33.5 ± 0.2	0.760 ± 0.009	15–30	z	0.0453	0.0018	0.0110	0.0096
UGC 12893	294.93	283.15	87.2 ± 3.8	0.131 ± 0.039	52–67	ok	0.0463	0.0022	0.0110	0.0093

Note. Galaxy center ${x}_{c},{y}_{c}$ is given in pixels, ELL ± dELL and PA ± dPA are the outer isophote ellipicity and position angle together with their standard deviations in the measurement range, given by RANGE (in arcsecs). FLAG indicates whether the inclination can be reliably estimated from the ellipticity (${i}_{\mathrm{disk}}={\mathrm{cos}}^{-1}(1-{\epsilon }_{\mathrm{outer}})$): ok = reliable, u = uncertain, z = nearly edge-on galaxy. SKY, DSKY, and RMS give the estimated sky level and its global and local variation (in MJy sr⁻¹). The last column ${\sigma }_{\mathrm{conf}}$ gives the estimated extra instrumental noise component during the Spitzer warm mission (See Section 2.2.5).

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

A scatter plot of P4 axial ratios versus HyperLeda values is presented in Figure 4 (upper left frame; only galaxies with flag = "ok" are shown). As anticipated, the P4 axial ratios are on the average closer to unity than those in HyperLeda, though the difference is not very large (median ${(b/a)}_{{\rm{P}}4}-{(b/a)}_{\mathrm{HyperLeda}}=0.024$). On the other hand the standard deviation of the difference is quite large (∼0.1). The upper right frame makes a similar comparison to P3 axial ratios (Muñoz-Mateos et al. 2015) which correspond to a fixed surface brightness level ${\mu }_{3.6}=25.5$ mag arcsec⁻². On average, P4 orientations are measured at about 0.9 times this distance. The scatter is now significantly reduced and no systematic difference is seen between P3 and P4. The lower frames in Figure 4 compares the position angles (now also "z" galaxies are included). In general the differences between P4 and HyperLeda are fairly small: for $b/a\lt 0.8$ the median absolute difference is 2°. The difference between P4 and P3 is even smaller (median absolute difference $0.9^\circ$ for $b/a\lt 0.8$). Nevertheless there are some exceptions, most notably NGC 4594 (The Sombrero Galaxy): in this case the P3 fixed isophote orientation corresponds to the extended halo, while the P4 orientation refers to the edge-on disk.

Zoom In Zoom Out Reset image size

Since we are fixing the disk orientations in the decompositions it is important to check the consistency of our inclinations. Figure 5(a) displays the histogram of the P4 axial ratios for Hubble types $-3\leqslant T\leqslant 10$. In case of a randomly oriented sample of thin disks, the distribution of $b/a$ should be flat. In case of finite vertical thickness a drop would be expected near a lower limit $b/a={q}_{i}$, where q_i is the intrinsic aspect ratio of the galaxies. According to Figure 5(a) such a drop is evident for $b/a\lesssim 0.15$. However, overall the sample contains an excess number of galaxies with small axial ratios $b/a\lesssim 0.5$ (see the dashed line in Figure 5(a)). Similar trend is seen also when using the HyperLeda axial ratios (Figure 5(b)) or P3 isophotal orientations (Figure 5(c)). A possible explanation for the excess of small $b/a$ ratios is that the S⁴G sample has been selected (Sheth et al. 2010) using an inclination-corrected blue magnitude limit (BT_corr = 15.5): if this dust correction were exaggerated, say for very late types, it would lead to an excess of faint, highly inclined galaxies. This explanation is supported by the solid curves in Figure 5 which display the histograms when limiting to galaxies with non-corrected $\mathrm{BT}\lt 15$: now the histogram of P4 values is quite flat. The histogram for P3 isophotal axial ratios is rather similar, though there are somewhat fewer small $b/a\lesssim 0.2$ values. This could be due to the above-mentioned faint stellar halos: in case of nearly edge-on galaxies a fixed surface brightness level could pick up the rounder faint outer envelopes, whereas in P4 we have in such cases tried to trace the disk isophotes. On the other hand, compared to both P3 and P4, the HyperLeda distribution has a clear deficit of large axial ratios, mostly likely due to the influence of inner non-axisymmetric structures.

Zoom In Zoom Out Reset image size

It is interesting to compare our sky background estimates to those in P3. In P3 an automatic sky measurement is made using 45 sky regions with 1000 pixels each. The regions are chosen close to the distance $2{R}_{25}$ from the galaxy center (R₂₅ is the blue band 25 mag isophotal radius from HyperLeda; if needed the distance of sky regions is modified manually). According to Figure 6 there is a very good agreement in the estimated sky backgrounds between P3 and P4 (see the right frame which takes into account that different P1 mosaics are used for some of the galaxies). This good agreement is remarkable as the measurements are made completely independently and with different methods. The median difference between the sky determinations (0.0006 MJy sr⁻¹) is only about 1% of the typical sky background value, and its standard deviation (0.003 MJy sr⁻¹) is comparable to the magnitude of global sky variations in both sets of estimates (see Figure 7). However, Figure 6 also reveals some cases where the difference between P4 and P3 is significant: inspection of the images indicates that this is due to a bright star (NGC 1055), a nearby interacting component (NGC 3327, NGC 4647), or a too small field of view (FOV; NGC 2655). In two cases (NGC 1300, UGC 10288) the final P1 mosaic used by P3 is much improved over the earlier version used in P4.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 7 compares our sky background variation estimates ("DSKY" denotes global variations between sky measurement regions and "RMS" the average of the locally determined rms-scatter) with the corresponding estimates in P3 (Muñoz-Mateos et al. 2015; their parameters ESKY1 and SSKY1, respectively). There is a good overall agreement in the level of estimated global variation (left frame): the somewhat larger values for P4 are likely to follow from the larger range of radii we used for the sky measurement regions compared to P3. Also the local sky rms values show good agreement (right frame).

2.2.4. Input Data Images

2.2.5. Sigma-images

The sigma-images quantify the statistical uncertainty of each image pixel and thereby determine the weights applied in GALFIT decompositions. This uncertainty contains two contributions: the noise contribution associated with the number of photons arriving at the instrument ("photon noise" or "shot noise"), and the noise originating from the instrument itself. The photon noise is assumed to follow a Poisson distribution, and it arises from two sources, the flux associated with the galaxy light and the flux coming from the sky background (zodiacal light). The main concern in the construction of the sigma-images is that the relative contributions of the photon noise and the instrumental noise are correctly estimated, so that correct relative weights are used in the decomposition for the bright central regions of the galaxies and for their faint outskirts.

The sigma-images are calculated using the pixel values and header information in the 3.6 μm data images and the pixel values of the weight images. The images provided by P1 are in flux units (MJy sr⁻¹), and for the calculation of the noise their pixel values F are converted to the number of electrons N_e,

$$\begin{eqnarray}&&{N}_{e}=\displaystyle \frac{F+{F}_{\mathrm{bg}}}{{F}_{\mathrm{conv}}}\;\times {T}_{\mathrm{frame}}\times {N}_{\mathrm{frames}}\times g,\end{eqnarray} \tag{ 8 }$$

where F_bg is the zodiacal light background which has been subtracted from the frame prior P1 by the automatic Spitzer pipeline (its value is given by the header keyword SKYDRKZB). Note that the flux F contains besides the galaxy light also the sky background which has been subtracted in P4, $F={F}_{\mathrm{gal}}+{F}_{\mathrm{sky}}$. The F_conv is the conversion factor between flux units and original digital units (header keyword FLUXCONV, in units of MJy sr⁻¹ per DN/s), T_frame is the integration time/frame in seconds, N_frames is the number of combined frames for each pixel, and g is the detector gain factor (GAIN in units of e/DN). The number of frames combined is coded to the pixel values W of the weight images, ${N}_{\mathrm{frames}}=W/10$. Note that ${T}_{\mathrm{frame}}=30$ s must be used instead of the original integration time/frame given by the header keyword FRAMTIME: this is because during the compilation of P1 mosaics the pixel values have been normalized to this value regardless of the original T_frame.²⁵ The statistical uncertainty of N_e in each pixel is then calculated as a combination of Poisson noise (photon noise) and the readout noise of the detector (RON),

$$\begin{eqnarray}&&{\sigma }^{2}({N}_{e})={N}_{e}+{N}_{\mathrm{frames}}\times {\mathrm{RON}}^{2}.\end{eqnarray} \tag{ 9 }$$

We use RON = 15.0, 14.6, and 21 electrons, for FRAMTIME = 12, 30, and 100 s, respectively. Note that these values, communicated by the Spitzer Science Center Helpdesk, deviate slightly from those given by the image header keyword RONOISE. The $\sigma ({N}_{e})$ is then converted to the estimated uncertainty of the image flux (note that $\sigma ({F}_{\mathrm{gal}})$ equals $\sigma (F)$ since F_sky is constant)

$$\begin{eqnarray}&&{\sigma }_{\mathrm{est}}(F)=\sigma ({N}_{e})\times {F}_{\mathrm{conv}}/({T}_{\mathrm{frame}}\times {N}_{\mathrm{frames}}\times g).\end{eqnarray} \tag{ 10 }$$

In order to assess the validity of this estimate we compare it to the actual noise measured directly from the image. In Figure 8 this is done for the sky measurement regions. In the left frame the measured sky rms (an average over all sky determination boxes) is plotted against the estimated ${\sigma }_{\mathrm{est}}$ from Equation (10). Colors distinguish between archival images from the cryogenic mission phase (original exposure time/frame either 12, 30, or 100 s) and the new observations during the warm Spitzer mission (time/frame 30 s, with the total exposure time of 240 s). For the archive images the overall agreement is quite good: there is a practically linear trend $\mathrm{RMS}\approx 0.9\ {\sigma }_{\mathrm{est}}$ holding for all three frame times, with the largest noise levels corresponding to the shortest frame times which have the largest contribution from the readout noise. The factor ∼0.9 is probably due to the P1 mosaicking process, during which the images have been combined and sampled to 075 pixel size from the native pixel size of 12. Because of this sampling the adjacent pixel values are strongly correlated, which is not taken into account in our theoretical estimate. Instead of trying to account in detail for the noise propagation during the mosaicking process we apply an empirical correction

$$\begin{eqnarray}&&{\sigma }_{\mathrm{decomp}}(F)=0.9\ {\sigma }_{\mathrm{est}}(F)\quad (\mathrm{cryogenic}\;\mathrm{mission})\end{eqnarray} \tag{ 11 }$$

to be used in decompositions of cryogenic phase archival images.

Zoom In Zoom Out Reset image size

In contrast, for the warm mission the observed rms is nearly 50% larger than the theoretical estimate (Figure 8), indicating the presence of an additional source of noise. Also, there is a noticeable drop in $\mathrm{RMS}/{\sigma }_{\mathrm{est}}$ ratio near the ecliptic plane (not present in the data from the cryogenic phase), indicating that the photon noise contribution to the ${\sigma }_{\mathrm{est}}$ (largest at $l\approx 0^\circ$) is overestimated compared to the instrumental contribution (constant with l). Following the advice of Spitzer Science Center Helpdesk, we include an additional instrumental noise component (${\sigma }_{\mathrm{conf}}$), which is added quadratically to the theoretical noise estimate. To account for the P1 mosaicking process, the multiplicative factor of 0.9 is again included. We thus adopt

$$\begin{eqnarray}&&{\sigma }_{\mathrm{decomp}}(F)=0.9\ \sqrt{{\sigma }^{2}{}_{\mathrm{est}}(F)+{\sigma }_{\mathrm{conf}}^{2}}\quad (\mathrm{warm}\;\mathrm{mission})\end{eqnarray} \tag{ 12 }$$

for the warm mission images. The value of the empirical correction term ${\sigma }_{\mathrm{conf}}^{2}$ is estimated by this formula when applied to the sky measurement regions; the same formula is then applied to all image pixels. The adopted values of ${\sigma }_{\mathrm{conf}}$ are listed in Table 1 above.

Figure 9 illustrates the magnitudes of different contributions to the sky background noise. For the archival images (cryogenic phase, left frame) the noise is dominated by the readout-noise, though the Poisson contribution due to zodiacal light (the F_sky we have subtracted in P4 + the SKYDRKZB subtracted during automatic Spitzer pipeline) still has a noticeable contribution. For the warm Spitzer mission (right frame) the extra noise term is even larger than the readout contribution.

Zoom In Zoom Out Reset image size

Nevertheless, at the central parts of the galaxies the photon noise associated with the galaxy light F_gal is the largest source of noise. This is illustrated in Figure 10 which compares the Poisson and background contributions as a function of surface brightness. The two horizontal lines indicate the typical background noise levels for the cryogenic (lower) and warm (upper) phases (includes both instrumental and noise due zodiacal light). The inset Figure illustrates how the σ-map looks for the galaxy NGC 3992 (observed during the warm mission). Near the center ($\mu \approx 17$), the photon noise due to F_gal completely dominates, though already in the bar region (μ = 20–21) both photon and instrumental contributions are important. Altogether the sigma-images and thus the applied relative weights between galaxy and background regions are intermediate between those typically encountered when decomposing ground-based optical and near-IR (NIR)-images. In the former case the photon noise due galaxy light usually dominates, while in the latter case the sigma-image is almost completely dominated by the background noise, so that the weight is almost constant for all pixels (this applies e.g., to Janz et al. 2014 GALFIT decompositions of Virgo dEs based on ground-based H-band images).

Zoom In Zoom Out Reset image size

In principle, the obtained ${\sigma }_{\mathrm{decomp}}$ is just a statistical estimate of the true underlying variance at each pixel. We did some experimentation by smoothing the sigma-images (median averaging with kernels amounting up to 20 pixels). Except in the case of a few galaxies with very centrally peaked light profiles, this smoothing had very little influence on the final decomposition parameters. For the galaxies where smoothing played a role, the derived parameters were in any case uncertain (for example, the bulge Sérsic index obtained unrealistically high values >10). In the end, we decided to apply no smoothing at all. Tests related to the sigma-images are presented in Section 4.3.

2.2.6. PSF-image

The IRAC data is not very well sampled: its native pixel resolution is 12, which is close to the Gaussian spread of a point source observed at channel 1. As discussed in detail in Peng et al. (2010), in such a case an oversampled PSF should be used. The IRAC PSF has also wide wings (see Figure 11), so that a relatively large convolution box size must be used in decompositions: we set this to 40'' × 40'' (in some cases with a very centrally peaked light profile this region was extended to 150'' × 150'' with considerable increase in CPU time). Note also that IRAC PSF depends slightly on the instrument orientation. Therefore, in principle a separate PSF should be used with each image, determined from point sources in the same frame, or a combination of appropriate PSFs, in case the final image is a combination of several images obtained at different times. Clearly, such a procedure would be very time consuming. Fortunately, such an accuracy is hardly needed in our decompositions. The common oversampled PSF provided by T. Jarrett was used for all images, made as a composite over several instrument rotation angles. Figure 11 displays the PSF, as well as a Gaussian profile (with $\mathrm{FWHM}=2\buildrel{\prime\prime}\over{.} 1$) approximately matching the core of the composite PSF. Also shown is an azimuthally averaged profile of the composite PSF. It will be shown in Section 4.1 that it is important to account for the central core, as well as for the nearly circular wings of the PSF, whereas the outermost spikes have less importance for the obtained decomposition parameters.

Zoom In Zoom Out Reset image size

The (ascii) input file for GALFIT specifies the galaxy data, mask, sigma, and PSF fits-files, and the region of the data image used in the decomposition. It also lists the components/functions used in the decomposition model, the initial guesses for the parameters, and specifies which of the parameters will be kept fixed, and which are iteratively varied in order to minimize the ${\chi }^{2}{}_{\nu }$. After convergence to a final solution, the final parameter values are written into an output file, with similar format as the input file. If needed, this output file can thus be used as an input for a new iteration (see Peng et al. 2002, 2010 for details).

The input file also specifies how to convert the image values to magnitudes. The data images from P1 are in flux units (MJy sr⁻¹). A conversion from pixel values F_i to (AB) surface brightness and integrated magnitudes is done with the formulas:

$$\begin{eqnarray}&&{\mu }_{3.6}=-2.5{\mathrm{log}}_{10}{F}_{i}+5{\mathrm{log}}_{10}\mathrm{pix}+\mathrm{zp}\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}&&{\mathrm{mag}}_{3.6}=-2.5{\mathrm{log}}_{10}\underset{i}{\displaystyle \sum }{F}_{i}+\mathrm{zp},\end{eqnarray} \tag{ 14 }$$

where pix = 075 and the zeropoint at 3.6 μm is zp = 21.097 (P3, Muñoz-Mateos et al. 2015). Values of pix and zp are inserted into GALFIT input file.

All P4 input files for 1-component (Sérsic) and 2-component bulge+disk (Sérsic+exponential) decompositions were generated automatically. Similarly, template files were created for the multi-component decompositions, which contained, in addition to bulge and disk components, entries for a Ferrers-bar, and a central unresolved PSF component. The user then manually choose which components are fit and which functions used in the final model (see Section 3 for more details). In all our decompositions we keep the centers of the components fixed to the galaxy center. The cases were this is clearly not appropriate (galaxies with off-center bulges and bars) are noted in the parameter files.

• 1.  
  In 1-component input files initial guesses are needed for five free parameters: the Sérsic index n, the effective radius R_e, the total magnitude m, the isophotal minor-to-major axial ratio q, and the position angle PA. The starting values of m and R_e were taken directly from the data (total galaxy magnitude and half-light radius), for the Sérsic index n = 2 is inserted as an initial guess, and q and PA were set to arbitrary values (0.9 and 10°, respectively). We thus avoided using the measured outer isophotes, to force GALFIT to search through a wider parameter space while minimizing the ${\chi }^{2}{}_{\nu }$. Typically 1-component fits converged after 10–20 iterations. When the fit did not converge, or if the final parameters were nonphysical (say, n > 10, $q\lt 0.05$, very large or small R_e), a new decomposition was started manually with new initial guesses. Usually this did not lead to any improvement, indicating that GALFIT is indeed very efficient in avoiding spurious local minima.
• 2.  
  The 2-component bulge-disk models apply a Sérsic-function for the bulge: they thus need guesses for the same Sérsic parameters as before, except that now these refer to the central component. Accordingly, we used the initial guess R_e (bulge) = 0.5 $\times {R}_{{\rm{e}}}$ (image) and ${m}_{\mathrm{bulge}}={m}_{\mathrm{image}}+1$. For the disk we use either "expdisk" or "edgedisk"-function, depending on the estimated galaxy inclination. In case of low or moderate inclination $b/a\gtrsim 0.2$ (corresponding to $i\lesssim 80^\circ$), we use the "expdisk" function, which needs two free parameters, the scale length (h_r) and the integrated magnitude of the disk, m_disk. We chose ${h}_{r}=$ 0.25 $\times \;{R}_{\mathrm{gal}}$ and ${m}_{\mathrm{disk}}={m}_{\mathrm{image}}+1$, thus starting with a model with fairly massive and extended bulge. The disk orientation was fixed to the shape of the outer isophotes determined from the image (see Section 2.2.3). In case of a nearly edge-on disk, $b/a\lesssim 0.2$, we use the "edgedisk" function with four free parameters: the central surface brightness ${\mu }_{0}$, radial scalelength h_r, vertical scalelength h_z, and the position angle of the disk. The first guesses are ${\mu }_{0}(\mathrm{disk})={\mu }_{0}(\mathrm{image})+3$, h_r as for the expdisk-model, while ${h}_{z}/{h}_{r}=0.1$. The position angle is left free, with PA_outer as an initial guess.
• 3.  
  In the template files for the multi-component fits the initial bulge and disk parameters are set as for the 2-component models. For the Ferrers-bar the free parameters are the surface brightness at the effective radius of the bar, ${\mu }_{{\rm{e}}}$, its outer truncation radius R_bar (denoted with r_out in Equation (7)), its axial ratio, and its position angle. As initial guesses we choose ${\mu }_{e}$ (bar) = ${\mu }_{{\rm{e}}}$ (image) + 3, ${R}_{\mathrm{bar}}=0.25\times {R}_{\mathrm{gal}}$, ${q}_{\mathrm{bar}}=0.5$, and ${\mathrm{PA}}_{\mathrm{bar}}={\mathrm{PA}}_{\mathrm{disk}}+90^\circ$. For the magnitude of the unresolved central component we used ${m}_{\mathrm{psf}}={m}_{\mathrm{image}}+5$. However, in practice we typically modified these pre-inserted template values even before starting the search of the final model, for example by adopting the output parameters from 2-component decompositions for the disk and bulge.

In its standard use, GALFIT is executed from the operating system command line, with an input file argument. This input file lists the input data files and the initial guesses for the parameters, as described above. The final decomposition parameters are written to an output file with a fixed name galfit.NN, where NN is a running number. Optionally, GALFIT makes a fits file containing the clipped data image (OBS; includes the region chosen for the fit), and total PSF-convolved model (MODEL), and the OBS-MODEL residual. Another GALFIT option is to write a FITS file containing model components in separate fits extensions.

In P4 we have used GALFIT via GALFIDL, which is a set of IDL routines designed for easy visualization of the output from GALFIT decompositions. In addition, GALFIDL includes wrapper routines for calling GALFIT from inside IDL, with the advantage that the GALFIT output files and the produced plots are automatically renamed in a systematic fashion, using the names of the input files. We have utilized this by coding the galaxy identification and decomposition model components to the name of each produced output file (see Appendix A)

The visualization options in GALFIDL follow those of the BDbar-decomposition program we developed earlier for the NIRS0S survey (Laurikainen et al. 2005), the most central of which is displaying a 2D plot of surface brightness versus distance from the galaxy center (see Figure 12). The advantage of this, compared to the more commonly used azimuthally averaged profile, is that the contributions of different model components, with different apparent ellipticities, are easily highlighted (Laurikainen et al. 2005; see also Gadotti 2008). The other visualization options include OBS-MODEL residual plots, profile cuts along a constant PA, comparison to observed profiles along isophotal major axis produced by IRAF ellipse, and plots showing schematically the different components included in the decomposition. The next section illustrates our decomposition strategies in more detail, concentrating on 2D-profiles. Additional plot types are illustrated in Appendix B, which describes the output released through a P4 web page for all S⁴G galaxies.

Zoom In Zoom Out Reset image size

The decompositions were made starting from simple 1 and 2-component models, and then adding as many components as necessary. For non-barred galaxies the process leading to the final model was the following.

• 1.  
  Accepting the automatic 1-component model (single Sérsic) as the final model. This was the case for elliptical galaxies (see NGC 3962 in Figure 13).
• 2.  
  Accepting the automatic 2-component bulge/disk decomposition as a final model. A typical example is NGC 3938 (Figure 14).
• 3.  
  Adopting a bulge/disk model, after interactively finding modified initial parameters that converged to an acceptable final fit.
• 4.  
  Adding a nucleus component or an inner disk (e.g., NGC 1357, Figure 15) to the bulge/disk model.
• 5.  
  When the galaxy had no obvious bulge we started from a single exponential disk, and if necessary, a second disk and/or nucleus was added (see NGC 723, Figure 16).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

When the outer profile was affected by a possible stellar halo, the outermost part of the profile was not fitted. The best model was vetted by looking at the original image, the residual image after subtracting the model, the 2D surface brightness profile, and the ellipticities of the structures. The value of final ${\chi }^{2}{}_{\nu }$ was not used as a criterion in assessing the relative merits of the models (often a simpler final model was preferred even if a more complicated model would have yielded slightly smaller reduced ${\chi }^{2}{}_{\nu }$).

Also, it is well known that many elliptical galaxies have small inner disks (Rest et al. 2001), and it has been shown by Huang et al. (2013) that many elliptical galaxies are better fitted with multiple Sérsic profiles. Nevertheless, such a detailed approach was not taken in this study, in which the emphasis is in the analysis of disk galaxies (paper 2). It is worth noticing that while using deep images like those in S⁴G, in an automatic fit the bulge profile even in late-type spirals is easily degenerate with the outer part of the disk. In automatic fits this may lead to an unrealistically large Sérsic n and R_e for the bulge, of which NGC 1357 is a good example (Figure 15).

For barred galaxies a similar step-wise approach was followed. NGC 936 (Figure 17) and NGC 5101 (Figure 18) are good examples demonstrating the importance of preventing the bar from mixing with the bulge flux. Adding a bar component to a simple bulge/disk model drastically changes the obtained properties of the bulge (for NGC 936 B/T drops from 0.46 to 0.19; for NGC 5101 from 0.99 to 0.22). NGC 5101 has also a type II profile in the disk break/truncation classification associated to a broad outer ring (Laine et al. 2014). Using the edge of the ring as a manifestation of a different flux distribution in the outer disk might be a bit misleading. Because of such ambiguities in the interpretation, we typically fit the type II disk profiles with a single exponential component. However, there are other barred galaxies in our sample, such as IC 4901 (Figure 19), in which two exponential components (+ Ferrers function for the bar) were used for fitting the disk. In this particular galaxy using two exponentials is necessary, and those clearly correspond to distinct surface brightness components. Note that the outer disk of NGC 5101 is clearly lopsided (see the residual plot in Figure 18). Such asymmetries are not taken into account in the pipeline decompositions (in case of strongly distorted galaxy no final model was made). See Zaritsky et al. (2013) for a detailed study of galaxy lopsidedness using S⁴G images.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

A third main group of galaxies in our sample are those having no obvious bulge. They may have a single exponential disk (NGC 3377A in Figure 20), or more than one disk component (NGC 723 in Figure 16). The structure fit as an inner disk in NGC 723 consists of broad, prominent, and tightly wound spiral arms. Bulgeless galaxies may have bars; NGC 3517 (Figure 21) is an example. To get a homogeneous estimate for the scale length of the disk for these galaxies, the outer disks were always fit with an exponential function, even in galaxies where the disk would have been somewhat better fit by a Sérsic function with n slightly less than unity. Generally, the assumption of an exponential disk is good, but there are also cases, like ESO 026-001 (Figure 22) in which a Sérsic function would actually be a better choice.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The GALFIT models for the nearly edge-on galaxies assume that the disk is viewed completely edge-on. A bulge, and in some cases also a bar or an additional thick disk component were included (Figure 23). In these models also the vertical thickness was an output parameter. However, these models are tentative, and are meant solely as starting points for better, scientifically oriented decompositions. There already exists detailed modeling of edge-on galaxies in S⁴G, based on fitting their vertical profiles to hydrodynamical thin-thick disk models (Comerón et al. 2011, 2012). Their radial luminosity profiles have been analyzed in Comerón et al. (2012, 2014) and Martín-Navarro et al. (2012).

Zoom In Zoom Out Reset image size

The P4 models for the spiral galaxies $T\gt 0$ are generally good, giving reliable estimates for parameters such as the B/T, the scale length of the disk (h_r), and its central surface brightness (${\mu }_{0}$). However, despite the fact that up to four components were fit, the pipeline decompositions for the early-type disk systems ($T\;\lesssim$ 1), because of their complex structures, are often insufficient. These systems may have nuclear bars, ovals, and lenses, which are not included in our models in any systematic fashion. Because of this, the pipeline B/T flux-ratios, particularly for S0 galaxies, can be over-estimated. Including all these structures will require even more complex decompositions, such as those done in the NIR by Laurikainen et al. (2005, 2006, 2009, 2010). Such time consuming modeling goes beyond the scope of our current P4 decompositions; nevertheless, the P4 decomposition output files provide good starting point for further fine-tuning.

As illustrated in Figure 11, the IRAC PSF has extended wings. Moreover, the PSF and the orientation of its asymmetric extensions vary from image-to-image, which has not been taken into account in our decompositions. To check the importance of the PSF wings, we compared differences in decomposition parameters obtained when the adopted composite PSF was replaced with a Gaussian PSF having the same $\mathrm{FWHM}=2\buildrel{\prime\prime}\over{.} 1$. Figure 24 compares the resulting effect on the Sérsic parameters in 1-component models. Clearly, decompositions with the Gaussian PSF yield n values that are systematically too small, differences reaching even tens of percents for some of the galaxies (though the median deviation is less than 5%). However, these rather large deviations are not representative of the true uncertainties, but rather give an idea of the magnitude of the error if the tails of the PSF were altogether ignored. A better measure of the actual uncertainty in P4 decompositions is obtained by comparing with an azimuthally symmetrized version of the composite PSF. Clearly, now the differences in n are much smaller (see the red symbols in Figure 24).

Zoom In Zoom Out Reset image size

We also checked the influence that the PSF has on the multi-component models. For that purpose we re-run all final decompositions that included both bulge and disk components (+ possible bar and center components; total of 524 models after excluding nearly edge-on galaxies), using both the Gaussian PSF and the symmetrized PSF. Table 3 lists the median of relative differences in bulge $B/T,n$, R_e, disk scale lengths h_r, and bar-to-total ratio Bar/T, when compared to the results obtained using the standard composite PSF. The largest differences are seen for the Sérsic parameters while using the Gaussian PSF, whereas h_r is barely affected. On the hand, the differences in parameters between those obtained using the composite PSF and using the symmetrized version are negligible. Based on these results we conclude that the spikes of the PSF have no significant effect as long as the nearly circular wings of the PSF are included. The use of single composite PSF for all S⁴G images should thus be acceptable.

Table 3.  The Effect of a Modified PSF on Final Decomposition Model Parameters

	GAUSSIAN PSF		SYMMETRIZED PSF
	Median(D)	Median($| D|$)	Median(D)	Median($| D|$)
B/T	−1.5%	4.5%	−0.1%	0.2%
n	−3.0%	7.8%	0.0%	0.5%
Re	8.9%	9.8%	0.1%	0.3%
hr	0.0%	0.3%	0.0%	0.0%
Bar/T	−1.5%	4.4%	−0.1%	0.2%

Note. D stands for the relative difference (e.g., $D=[n(\mathrm{mod})-n(\mathrm{ori})]/n(\mathrm{ori})$), where "ori" refers to the standard composite PSF. Medians are used to characterize the typical deviations and the scatter, to eliminate spurious cases where the decompositions with Gaussian PSF converged to a different type of solution.

Download table as:  ASCIITypeset image

In principle, poor sky subtraction can severely affect the decomposition results, in particular the parameters of the disk. To constrain the possible magnitude of such uncertainties, we re-run the multi-component decompositions that included both bulge and disk components (+ possible bar and center components; same 524 models as above). Two additional sets of sky values, SKY' = SKY ± DSKY, were used, where DSKY was the standard deviation of the different sky regions. Figure 25 shows the effect on the scalelength of the disk. Although individual changes can in few cases be large (${h}_{r}(\mathrm{mod})/{h}_{r}(\mathrm{ori})\gt 1.2$ in nine cases when too small a sky is subtracted), the median differences are less than 2% (and even smaller in the other parameters of interest, see Table 4). The sky subtraction is not a concern in the current decompositions.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 4.  The Effect of Sky Subtraction on Final Decomposition Model Parameters

	SKY + DSKY		SKY−DSKY
	Median(D)	Median($| D|$)	Median(D)	Median($| D|$)
B/T	0.2%	1.5%	−0.1%	1.6%
n	−1.4%	1.8%	1.8%	2.0%
Re	−1.2%	1.4%	1.7%	1.8%
hr	−2.5%	2.6%	2.6%	2.8%
Bar/T	0.2%	1.5%	−0.1%	1.6%

Note. D stands for the relative difference (e.g., $D=(n(\mathrm{mod})-n(\mathrm{ori}))/n(\mathrm{ori})$), where "ori" refers to the standard sky subtraction.

Download table as:  ASCIITypeset image

The weights applied to various pixels have an important role in decompositions, in particular when the galaxy structure is complicated, so that the differences between the applied model and the true structure are large. As mentioned in Section 2.2.5 the σ-image itself is a statistical estimate of the underlying σ in each pixel, so it might be reasonable to smooth it before applying it in the decompositions. In Figure 26 (left column) we examine the effect of sigma-image smoothing on the derived bulge parameters. A median filter is applied with a width of 5 pixels. Clearly the effect is quite small except for a few deviant cases marked on the plot. In these cases the bulge parameters are sensitive also to changes in the PSF or the sky background level.

For comparison, Figure 26 (right column) also illustrates the changes in bulge parameters if a constant sigma is assumed at all image pixels. A constant sigma exaggerates the relative weight of the central regions compared to the outskirts. Besides a large scatter, a systematic increase of the estimated n is also obvious: the median ${n}_{\mathrm{mod}}/{n}_{\mathrm{ori}}=1.2$ (the mean ratio is 1.4). What typically happens is that the fit tries to reproduce the central peak with an increased n, even if the outer disk then becomes too bulge dominated. Indeed, the bias (and the scatter) is particularly large for earlier-type disks (open circles in the plot indicate T ≤ 5; median ${n}_{\mathrm{mod}}/{n}_{\mathrm{ori}}=1.25$). This comparison reminds us that when decomposition parameters from different studies are compared, it is also important to make sure that similar weights have been applied.

Automatic 2-component Sérsic-exponential (or Sérsic–Sérsic) models are often applied to large data surveys. This is a natural approach as the data quality (depth/angular resolution) might be insufficient for more detailed modeling so that the large effort in multi-component decompositions does not seem justified. Moreover, it has been recently claimed (Tasca & White 2011) that 2-component decompositions (Sérsic + exponential) are sufficient also for fitting barred galaxies. Their argument was based on obtaining similar average B/T ratios for barred and non-barred galaxies in their 2-component bulge/disk models. They reasoned that if the omission of the bar were a problem it should manifest as a higher B/T for barred galaxies. However, to accurately address this matter one has to compare the different types of decompositions (2-component and multi-component) for well-defined samples of barred/non-barred galaxies.

Such a comparison between different decomposition models is shown in Figure 27. Again, those galaxies for which the final model contains both a bulge and a disk are studied. For the non-barred galaxies (those with no bar-component; leftmost column in the Figure) the bulge parameters (Sérsic n, B/T, ${R}_{{\rm{e}}}/{h}_{r}$) in automatic 2-component runs are almost identical to those in the final models. This agreement is expected because over 80% of the final non-bar models are just Sérsic-expdisk models (15% have two disk components, and 2% have an extra central component), and typically the automatically found BD models did not need any refinement. For barred galaxies (those with a bar-component in the final model; middle column), the obtained median values depend drastically on whether the bar is included. This result emphasizes that the examples of decompositions given in Section 3, highlighting the importance of modeling the bar (e.g., Figures 17 and 18) were not exceptional cases. Overall, ignoring the bar increases the estimated B/T ratios by a factor of 2–3 because of gross (even by as much as a factor of 5) overestimate of R_e and n. For example, for spirals in the range $1\leqslant T\leqslant 5$ the 2-component decompositions suggest $n\gtrsim 4$ whereas the multi-component runs indicate $n\approx 1-2$. Altogether, in the final models the difference in bulge parameters obtained in the multi-component decompositions for barred and non-barred galaxies is fairly small (right column in Figure 27).

Zoom In Zoom Out Reset image size

The conclusion that multi-component decomposition models are essential to measure realistic bulge parameters for barred galaxies is not new (Laurikainen et al. 2006, 2007; Gadotti 2008; Weinzirl et al. 2009). A similar conclusion, based on synthetic images, was reached also by Laurikainen et al. (2005).

In Figure 28 we compare the combined bar/non-barred sample of the previous Figure with the decompositions in Laurikainen et al. (2007). Because of the large fraction of barred galaxies, the difference in the obtained bulge properties between the 2-component and multi-component decompositions remains significant, even when barred and non-barred galaxies are considered together. We find an excellent agreement between the current multi-component results and those in Laurikainen et al. (2007), obtained with a different decomposition code (BDBAR; however, BDBAR uses IDL Curvefit and is thus based on the same Levenberg–Marquadrdt minimization as GALFIT), and based on different NIR image data. It is worth noting that in these decompositions the Sérsic n for Hubble types Sa-Sc is nearly $n\;\sim$ 1, whereas in the decompositions by Tasca & White (2011), for the same Hubble types, the Sérsic index is peaked at $n\;\sim$ 4. Small values of the Sérsic index, similar to ours for these Hubble types, are reported also by Graham & Worley (2008).

Zoom In Zoom Out Reset image size

One of the main goals of Pipeline 4 is to obtain measurements for the galaxy size–magnitude scaling relations. In order to be consistent with earlier analysis (e.g., Courteau et al. 2007) the P4 final models as a default use single exponentials for the disk. However, deep optical and NIR surveys (Erwin et al. 2005; Pohlen & Trujillo 2006; Gutiérrez et al. 2011; Muñoz-Mateos et al. 2013) have shown that only a fraction of galactic disks ($\sim 1/3$) are simple exponentials (=Type I in Pohlen & Trujillo 2006 classification). Instead the typical brightness profiles consist of two (sometimes three) exponential subsections with different radial slopes. When the outer disk has a steeper slope, the galaxy is classified as possessing a Type II break ("truncation"), and conversely if the outer slope is more shallow, it is classified as a Type III break ("antitruncation"). Kim et al. (2014) have recently made 2D decompositions for 144 barred S⁴G galaxies taking into account disk breaks in their decompositions with the BUDDA code (de Souza et al. 2004; Gadotti 2008, 2009). Their fitting function for the disk consists of two exponential sections, with different scale-lengths (h_in and h_out) inside and outside the break radius R_break. They also made decompositions where they fitted the disk with a single exponential component. Their result indicate that the inner scale lengths for two-component disks are typically about 40% longer than the scalelengths obtained in single disk fits; they thus conclude that "it is important to model breaks in Type II galaxies to derive proper disk scale lengths."

Nevertheless, it is not always obvious what is the "proper" disk scale length estimate to use in various scaling relations, in case the galaxy exhibits several exponential subsections. For example, it is well known (Pohlen & Trujillo 2006) that Type II breaks are often connected to outer rings associated with bar OLR resonances. Such breaks are indeed dominant for early type barred disks ($T\lt 3$; Laine et al. 2014). Since the bar torques are able to push material from the CR regions out toward OLR, this will promote a shallower distribution inside the break radius. However, beyond the OLR, the effect of the bar is insignificant, so that the underlying disk can remain more or less intact. In such a case, it might in fact be the outer, rather than the inner, scalelength that would better characterize the original overall mass distribution. On the other hand, for later Hubble types the Type II break is often connected with the end of prominent spirals (Laine et al. 2014) and could be due to suppressed star formation; for such a case, the inner scale length might indeed be more appropriate to characterize the disk as a whole. Laine et al. (2014) also find that for such spiral-related breaks the ratio ${h}_{\mathrm{inner}}/{h}_{\mathrm{outer}}$ is typically closer to unity than for OLR related breaks.

Figure 29 compares Pipeline 4 decompositions with several recent disk truncation studies, which use subsamples of the same S⁴G data set. In this plot the disk scalelengths are displayed against the stellar mass derived in P3 (Muñoz-Mateos et al. 2015). Besides the above-mentioned Kim et al. (2014) 2D decomposition study, we also compare with Muñoz-Mateos et al. (2013) and Laine et al. (2014), where fits to one-dimensional profiles were conducted. First of all, the figure (upper row) indicates a very good agreement for the scale lengths of Type I profiles between all four studies, conducted with independent methods. Second, it illustrates the significant difference between the inner and outer slopes for Type II (and III) profiles, amounting to roughly a factor of 2 (see Muñoz-Mateos et al. 2013; Kim et al. 2014). The P4 single disk scalelengths seem to fall quite close to being a geometric mean of h_inner and h_outer derived in earlier studies.

Zoom In Zoom Out Reset image size

To emphasize the possible "unperturbed" nature of Type II outer disks, we compare in Figure 30 the P4 scale lengths versus stellar mass for Type I galaxies with the Type II outer scale lengths derived in the above mentioned disk break studies. Indeed, the differences between the Type I single disk h and the Type II h_outer are quite small, much smaller than the differences compared to h_inner. The fits to the data also give the impression that the ${h}_{\mathrm{inner}}/{h}_{\mathrm{outer}}$ ratio approaches closer to unity for less massive galaxies: this is in accordance with the above-mentioned dominance of spiral-related less abrupt truncations for later, and thus on average less massive, spirals.

Zoom In Zoom Out Reset image size

The Pipeline 4 single exponential fits have a convenient feature of representing an effective average over inner and outer disks (when both present). They thus provide a homogeneous set of robust scale measurements, not sensitive to factors modifying the local slopes. Nevertheless, a possible caveat is that the fitted effective single h might become dominated by different degrees by the inner/outer parts, depending on the galaxy surface brightness. For example, the estimated h might be biased toward h_inner when the disk central surface brightness decreases toward less massive galaxies: this would be the case if the image depth was not sufficient to cover the galaxy regions beyond the break radius. Figure 31 addresses this potential problem by comparing the trends of the break radii with respect to galaxy mass, to that of the galaxies' visual outer extent (R_gal, see Section 2.2.3; a similar trend would result if ${R}_{25.5}$ were plotted instead of R_gal). The figure indicates that a break, if present, should be detectable through the whole range of S⁴G galaxy masses.

Zoom In Zoom Out Reset image size

In summary, we feel confident that the single disk fits provide a useful overall estimate of the disk original scale length (and its extrapolated surface brightness), though especially in case of barred massive galaxies secular evolution might have led to significant deviations from simple exponentials, important to include in detailed models for individual galaxies. Moreover it is likely that the slope differences associated with breaks are smaller for later types, which form a vast majority of S⁴G galaxies.

Nevertheless, as concluded by Kim et al. (2014), estimates of other decomposition parameters, such as the B/T for massive galaxies would become more accurate if the inner slopes are accounted for (say, leading to less disk light assigned to bulge). The situation is somewhat analogous to the benefit of including additional inner components like bars (Laurikainen et al. 2005; Gadotti 2008), lenses in S0s (Laurikainen et al. 2010), or barlenses (Laurikainen et al. 2014) into decompositions. However, for the goals of Pipeline decompositions, the expected magnitude of changes (about 10% relative change in B/T according to Kim et al. ) is quite small, compared to the uncertainties related to choice of the decomposition model components (say, including a bar versus ignoring it). The choice of the code might also sometimes have a bigger effect. For example, Kim et al. (2014) use NGC 936 as an example of Type II galaxy (see their Figure 4). For this galaxy they fit a break at 98'' and derive ${h}_{\mathrm{inner}}$= 53'' and ${h}_{\mathrm{outer}}$= 28'', all very close to the measurements in both Muñoz-Mateos et al. (2013) and Laine et al. (2014). On the other hand, the Pipeline 4 single disk fit (see Figure 17) gives h = 40''. We verified that truncating the disk in GALFIT decompositions at the break radius given in Kim et al. reproduces their inner slope quite well (we get 57''). At the same time, the B/T we obtain increases slightly (from 0.19 to 0.22), as anticipated by Kim et al.²⁷ Nevertheless, the B/T we obtain after accounting for the more shallow inner slope is still nearly 50% smaller than the value obtained by Kim et al. (B/T = 0.32), probably because of some model/code dependent factors, such as how the image pixels are weighted, or the PSF is treated.

The full S⁴G sample contains 2352 galaxies, chosen according to their internal extinction corrected B-magnitude (${M}_{\mathrm{Bcorr}}\lt 15.5$), apparent B-band 25 mag isophotal diameter (${D}_{25}\;\gt$ 60''), galactic latitude ($| b| \gt 30^\circ$), and H i recession velocity (${V}_{\mathrm{radio}}\lt 3000\;\mathrm{km}\;{{\rm{s}}}^{-1}$), obtained from the HyperLeda database. Due to its mag-limited character, it contains a large number of low surface brightness late-type spirals and irregulars. Also, galaxies with peculiar morphology were not specifically excluded. In some cases the FOV is not large compared to the galaxy size (the new Spitzer observations mapped regions covering at least 1.5 D₂₅, but this condition was not fulfilled by all the archival galaxies included in the sample). In such cases, the sky background is difficult to estimate reliably, and for galaxies near the ecliptic, the background may have larger gradients (see Muñoz-Mateos et al. 2015). Altogether, the sample contains a number of galaxies for which decompositions are less reliable, or not possible at all to carry out.

Because it is important to estimate the reliability of the derived structural parameters, we have assigned a quality flag to each galaxy, running from 1 (worst case) to 5 (most reliable). The judgment was done partly by visual inspection of the original data and partly by evaluating the decomposition models.

• 1.  
  Quality = 1 (31 cases)Reasons: Bad original data (very bright overlapping stars, strongly varying background, image defects).Action: Excluded from all analysis: galaxy identifications are listed in parameter tables but no parameter values are given. P4 web page illustrates the raw data + mask, but not any decomposition models.
• 2.  
  Quality = 2 (44 cases)Reasons: Original data is more or less fine, but the FOV is too small for reliable sky estimation. Alternatively, galaxies exhibit strongly distorted shapes which make even 1-component fits unreliable (mergers, interacting, peculiar, strong warp, very lopsided).Action: 1-component Sérsic fit was done and the parameters are given in Table 6 (and in the web pages), together with a comment indicating that they need to be taken with caution. Multi-component decomposition parameters are not given.
• 3.  
  Quality = 3 (61 cases)Reason: Original data is fine but the galaxies have complex structures that require detailed multi-component models beyond the scope of the pipeline decompositions (which have a maximum of four components).Action: 1-component Sérsic fit is considered reliable (Table 6). Multi-component decomposition was also made and parameters are listed in Table 7, with a cautionary comment. Web pages show both 1-component and multi-component decompositions.
• 4.  
  Quality = 4 (406 cases)Reason: Original data and decomposition are of good quality. However, the galaxy was either highly inclined (contained a "z" component; 333 cases), or it had complicated structure, so that there might be a degeneracy between model components (such as between inner and outer disk components; 73 cases).Action: All decomposition parameters given in Tables 6 and 7, and illustrated in the P4 web pages. However, these are omitted from the analysis of disk central brightness and scale length in paper 2.
• 5.  
  Quality = 5 (1810 cases)Reason: Both the original data and the decompositions are of good quality.Action: All decomposition parameters are given in Tables 6 and 7, and illustrated in the P4 web pages.

Table 5 summarizes the number of galaxies in different quality categories.

The output parameters of 1-component Sérsic fits are listed in Table 6. For 1-component fits the parameters are the Sérsic index n, effective radius R_e, integrated magnitude mag, axial ratio q, and major axis position angle PA (the centers are fixed to those given in Table 1; the isophotes are assumed to be elliptical and to have fixed a shape and orientation with radius). Additionally, there is a column indicating the reliability of the fit.

Single-component Sérsic-fits are routinely used in large data surveys. This gives objective, easily reproducible results, that provide useful characterization of the galaxy global characteristics. For example, Cappellari et al. (2013) argued that Sérsic n > 4 largely finds the most massive early-type galaxies ($M\gt 3\cdot {10}^{11}{M}_{\odot }$), which are also the slow rotators in their kinematic classification. Also, there are well-known correlations between galaxy color and Sérsic index. Therefore, in paper 2 we will present detailed analysis of the 1-component Sérsic parameters for the S⁴G sample. Here we report just the dependence of n on the morphological type.

Figure 32 displays the histogram of Sérsic index-values in the 1-component models. Galaxies are divided into four bins according to their mid-IR Hubble type (E with $T\leqslant -4$, S0 with $-3\leqslant T\leqslant 0$, spiral with $1\leqslant T\leqslant 9$, and irregulars T = 10). Clearly, the distribution of galaxies peaks at $n\;\approx$ 1–2, with a broad tail to larger values of n. There is also a secondary peak close to $n\approx 4$ (corresponds to a de Vaucouleurs profile), but this is not very prominent. The overall distribution reflects the nearly exponential profiles of many late type spirals, which dominate the S⁴G sample. For irregulars, the distribution peaks at $n\approx 1$. For S0's, the distribution is much broader. The distributions remain essentially the same for less inclined galaxies, for example if the sample is limited to those with apparent $b/a\gt 0.5$.

Zoom In Zoom Out Reset image size

For multi-component decompositions the tabulation is more complicated, as model components/functions vary from one galaxy to another. Also, the same function can be used to describe different structure components in different galaxies. We have decided to present the multi-component parameters in two different formats, one that is compact and easily human-readable, another more suited to automatic reading. In the first format (Table 7) the first entry for each galaxy indicates the used model and the number of components. The next lines, for each component included in the model, indicate the physical interpretation of the component (B, D (or Z), bar, N), and the GALFIT function used (sersic,expdisk, edgedisk, ferrer2, psf), followed by the component parameters. The Table caption specifies which parameters are listed for each function. The other table (available via IRSA and P4 web page) lists for each galaxy all possible components and their parameters: empty values indicate that this component was not included in the decomposition of this galaxy.

The S⁴G sample contains 358 galaxies which were considered to be close to an edge-on view and are excluded from further analysis in paper 2. Also, 26 are elliptical systems, modeled with a single Sérsic function. This leaves 1855 moderately inclined disk systems. As discussed in Section 3 there are over 20 different combinations of functions/components used in the final models, so that there is a need to group the decompositions in to major categories. A natural approach is to base this grouping on whether the decomposition model has a bulge component. Because the bulge can be modeled either with a "sersic" or "psf" component, depending on its apparent size, we have two categories, BD and ND models, respectively. When there is no trace of a bulge, the system can be either a single disk (D), possess a bar-like component (Dbar), or contain both inner and outer disk components (DD).

The numbers and relative fractions of galaxies in these categories as a function of Hubble type are shown in Figure 33. Here the mid-IR classification from Buta et al. (2014) is followed. Apparently the relative fraction of BD-models decreases gradually toward the late-type spirals. However, taking into account that most of the ND models describe small bulges (support for this claim is given in paper 2), indicates a much smoother distribution of galaxies with bulges, dropping rapidly only above $T\sim 5$. The fact that ND models cover a relatively large range of Hubble types, including S0s, is quite interesting because it indicates that S0s can possess very small bulges. This is in agreement with Laurikainen et al. (2010), where the same conclusion was made based on decompositions of NIRS0S data. The result is consistent with the idea that at least some S0s might be former late-type spirals with small bulges, devoid of gas, following quenching of star formation. A very small B/T in an early-type spiral has also been reported also by Kormendy & Barentine (2010). Another interesting feature in Figure 33 is that many galaxies that lack bulges, can still have bars. These Dbar galaxies peak at Hubble types T = 7. Beyond T = 9 they are replaced with single disks, becoming almost the sole type of models for irregulars (T = 10). The DD-models are most common (about 10%) for T = 9.

Zoom In Zoom Out Reset image size

Finally, Figure 34 gives examples of galaxies in these major categories, indicating the model components for four low mass and four large mass systems in each category. For BD and ND categories the barred and non-barred models are also distinguished. Note that "barred"/"non-barred" refers here to whether or not a bar component was included to the the final decomposition model, not to any morphological classification; for example a non-barred BD model has been adopted for NGC 5985, which has a SAB family classification (see Figure 33). A detailed comparison to Buta et al. (2014) classification will be presented in Paper 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size