2D stellar population and gas kinematics of the inner 1.5 kpc of the post-starburst quasar SDSS J0210-0903
Post-Starburst Quasars (PSQs) are hypothesized to represent a stage in the evolution of massive galaxies in which the star formation has been recently quenched due to the feedback of the nuclear activity. In this paper our goal is to test this scenario with a resolved stellar population study of the PSQ J0210-0903 ,as well as of its emitting gas kinematics and excitation. We have used optical Integral Field Spectroscopy obtained with the Gemini GMOS instrument at a velocity resolution of 120 km s and spatial resolution of 0.5 kpc. We find that old stars dominate the luminosity (at 4700 Å) in the inner 0.3 kpc (radius), while beyond this region (at 0.8 kpc) the stellar population is dominated by both intermediate age and young ionizing stars. The gas emission-line ratios are typical of Seyfert nuclei in the inner 0.3 kpc, where an outflow is observed. Beyond this region the line ratios are typical of LINERs and may result from the combination of diluted radiation from the nucleus and ionization from young stars. The gas kinematics show a combination of rotation in the plane of the galaxy and outflows, observed with a maximum blueshift of km s. We have estimated a mass outflow rate in ionized gas in the range yr and a kinetic power for the outflow of . This outflow rate is two orders of magnitude higher than the nuclear accretion rate of yr, thus being the result of mass loading of the nuclear outflow by circumnuclear galactic gas. Our observations support an evolutionary scenario in which the feeding of gas to the nuclear region has triggered a circumnuclear starburst 100’s Myr ago, followed by the triggering of the nuclear activity, producing the observed gas outflow which may have quenched further star formation in the inner 0.3 kpc.
keywords:galaxies: active – galaxies: kinematics and dynamics – galaxies: stellar content – galaxies: starburst
Our understanding of galaxy evolution has undergone a revolution in the past decade. Of particular interest is the correlation that links the nuclear supermassive black holes (SMBH) mass to the bulge stellar velocity dispersion, the relation (Gebhardt et al., 2000a, b; Merritt & Ferrarese, 2001; Tremaine et al., 2002; Ferrarese & Ford, 2005; Gültekin et al., 2009; Graham et al., 2011). In this scenario, Post-Starburst Quasars (PSQs) may acquire a great importance, since they seem to represent a critical phase in the secular evolution of galaxies that links the growth of the stellar bulge and that of the SMBH.
PSQs are broad-lined AGN that show Balmer jumps and high-order Balmer absorption lines from A stars characteristic of massive post-starburst populations with age of a few hundred Myrs. PSQs are hypothesized to represent a stage in the evolution of massive galaxies in which both star formation and nuclear activity have been triggered and are visible simultaneously before one or the other fades. Brotherton et al. (1999) reported the discovery of a spectacular PSQ that they propose to represent a stage in hierarchical galaxy evolution following a merger and enshrouded Ultra Luminous Infrared Galaxy (ULIRG) phase. HST images have shown that the host galaxy is in fact a post-merger remnant, and its spectrum shows evidence for multiple stages of star-formation as might be expected through such a process (Canalizo et al., 2000; Brotherton et al., 2002; Cales et al., 2011).
There are at least two possibilities to connect the presence of a post-starburst population and nuclear activity: (1) the flow of gas towards the nucleus first triggers star formation in the circumnuclear region and the nuclear activity is triggered after hundreds of Myrs; in the meantime, the star formation may cease due to exhaustion of the gas (Storchi-Bergmann et al., 2001; Davies et al., 2007); (2) the flow of gas towards the nucleus triggers star formation in the circumnuclear region and the nuclear activity, when triggered, quenches the star formation (Granato et al., 2004; Di Matteo et al., 2005; Hopkins et al., 2006; Cano-Díaz et al., 2012). In the first scenario, the fueling of the black hole usually occurs with a delay of a few 100 Myrs after the beginning of the star formation – in order to explain the systematic excess of intermediate age stellar population contribution in AGN hosts, although there are a few cases in which simultaneous on-going star formation and nuclear activity are observed (Storchi-Bergmann et al., 2001). In the scenario proposed by Davies et al. (2007) the gas fueling the nuclear activity originates from winds of AGB stars, which can be accreted more efficiently than the ejected mass by O and B stars and/or supernova. In these two versions of the first scenario, the star formation ceases due to the exhaustion of gas, and the post-starburst population is distributed in the circumnuclear region extending even to the nucleus. In the quenching scenario, the star formation occurring in the vicinity of the nucleus is abruptly interrupted due to feedback from the AGN (e.g., nuclear outflows, radiation and jets). Our proposition is that, in this case, the post-starburst population should be located in regions suffering the feedback effects – e.g., over regions where a gaseous outflow is observed, which can be mapped via the gaseous kinematics. With the goal of investigating the nature of the connection between the post-starburst stellar population and nuclear activity in PSQs, we began a program to map the stellar population and the manifestations of nuclear activity in the inner few kpc of PSQs using integral field spectroscopy.
In this paper we present the results obtained from observations of the PSQ SDSS J021011.5-090335.5 (hereafter PSQ J0210-0903), selected for its proximity, at a distance of
only 170 Mpc (from NED
This paper is organized as follows. In Section 2 we describe the observations and reduction processes. In Section 3 we present the methodology and results of the stellar population analysis. In Section 4 we present the methodology of analysis of the gas emission and the gas emission-line flux distributions and ratios. In Section 5 we report our results for the gas kinematics. In Section 6 we discuss and interpret our results as well as discuss the technique of Principal Component Analysis to study the gas kinematics. In sections 6.4 and 6.5 we present an estimate to the mass of the emitting gas and the mass outflow rate, respectively. In Section 7 we present a summary of our results as well as our conclusions.
2 Observations and data reduction
Two-dimensional spectroscopic data of the PSQ J0210-0903 were obtained in 2008 December 7 using the Gemini Multi-Object Spectrograph Integral Field Unit (GMOS-IFU) (Allington-Smith et al., 2002), in the program GN-2008B-Q-45. The observations were obtained in two-slit mode, covering a field of view (FOV) of 7.0 5.2 arcsec, using the B600G5303 grating and gG0301 filter, resulting in a spectral range of Å and a wavelength sampling of 0.913 Å pixel at a velocity resolution R 3000 (FWHM 120 km s). The total exposure time was 81700 s. The seeing during the observation was approximately 06, corresponding to a spatial resolution at the distance of the galaxy of 0.5 kpc.
Data reduction was accomplished using generic IRAF
In the top-left panel of Figure 1 we present the Gemini-GMOS -band acquisition image of the PSQ, which we have rotated to the same orientation of the IFU observations. In the top-right panel we present the H flux distribution within the IFU FOV, including the outline of a set of circular regions with 08 diameter from which we have integrated the spectra used in the stellar population synthesis (see next section). In the middle panels of the figure we present rest-frame spectra from the nucleus (the AGN, at region 1) together with those from two circumnuclear regions (regions 4 and 5) and from an HII region, at the more distant region 7. In the bottom panel the SDSS rest-frame spectrum shows the characteristic emission lines of the optical range, as well as the absorption lines from intermediate age stars, namely, the high-order Balmer lines between Å.
3 Stellar population
We have analyzed the stellar population using both the SDSS spectrum — which corresponds to an aperture of 3 ( 2.5 kpc at the distance of the galaxy) – covering the wavelength range Å and our integral field spectroscopy through the use of the technique of spectral synthesis.
The Integral Field Spectroscopy could in principle be used to obtain information on the spatial distribution of the stellar population applying the synthesis to the spectra of each pixel. However, the signal-to-noise ratio (SNR) of individual spectra was not high enough, thus, in order to increase it, we integrated the spectra within the 08 apertures described above and indicated in Figure 1. These regions were selected as sampling extranuclear regions with high enough S/N ratio in the absorption spectra to allow performing spectral synthesis and being separated enough from the nucleus in order to avoid possible contamination from the nuclear light due to the PSF wings.
In order to perform the spectral synthesis we used the STARLIGHT code (Cid Fernandes et al., 2004, 2005, 2009), which searches for the linear combination of simple stellar population (SSP), from a user-defined base, that best matches the observed spectrum. Basically, the code fits an observed spectrum solving the following equation for a model spectrum (Cid Fernandes et al., 2004):
where is the synthetic flux at the normalization wavelength, x is the population vector, whose components represent the fractional contribution of the each SSP to the total synthetic flux at , is the spectrum of the th SSP, with age and metallicity , normalized at , is the reddening term, and is the Gaussian distribution, centred at velocity with dispersion , used to model the line-of-sight stellar motions. The reddening term is modeled by STARLIGHT due to foreground dust and parametrized by the V-band extinction so that all components are equally reddened and to which we have adopted the Galactic extinction law of Cardelli, Clayton & Mathis (1989) with .
The fit of the model to the observed spectrum is carried out minimizing the following equation (Cid Fernandes et al., 2004):
where is the weight spectrum, defined as the inverse of the noise in . Emission lines and spurious features are masked out by fixing at the corresponding . The minimum of equation 2 corresponds to the best parameters of the model and the search for them is carried out with a simulated annealing plus Metropolis scheme. A detailed discussion of the Metropolis scheme applied to the stellar population synthesis can be found in Cid Fernandes et al. (2001).
We constructed the spectral base with the high spectral resolution evolutionary synthesis models of Bruzual & Charlot (2002) (BC03), where the SSPs cover 11 ages, , , , , , , , , , and yrs, assuming solar metallicity (). We have used the SSP spectra constructed from the STELIB library (Le Borgne et al., 2003), Padova-1994 evolutionary tracks and Chabrier (2003) Initial Mass Function (IMF). In order to account for the AGN featureless continuum (FC) a non-stellar component was also included, represented by a power law function (). In accordance to Cid Fernandes et al. (2004) we have binned the contribution of the SSPs into three age ranges: young (, where (t5 Myrs) represents the ionizing population and ( Myrs) represents the non-ionizing one), intermediate ( to ) and old ( to Gyr). Cid Fernandes et al. (2005) have shown that the star formation history of a galaxy may be very well recovered by this condensed population vector, as the individual contributions of each SSP are very uncertain. We have used this conservative spectral base, because we do not have sufficient constraints in our spectra to allow the use a broader spectral base. This procedure has been extensively used in previous studies (e.g., Cid Fernandes et al., 2004, 2005) and the investigation of the impact of the corresponding assumptions in the fit of the observed spectra is beyond the scope of the present paper. Our goal here is to investigate the relative contribution of broad age components as a function of location in the observed field.
The results of the synthesis are shown in Table 1 and Figure 2; they are presented in two sets, one for the SDSS spectrum and another for the IFU integrated spectra, which are discussed in the next two sections.
3.1 Synthesis of the SDSS spectrum
The SDSS spectrum shows clear signatures of intermediate age stars which contribute approximately 13% of the total luminosity at 4700 Å according to the synthesis (first line of Table 1). The old population contributes with 36% while the young and FC contribute with 18% and 32% respectively.
Our IFU data cover a limited spectral range, from 4260 to 5230 Å, thus we investigate the effect on the synthesis of this restricted range using the SDSS spectrum. We have labeled as the spectrum covering this range in Table 1. For this spectral range the synthesis results show an increase in the contribution of the old population, which is 41%, as well as of the young and FC components, which contribute with 33% and 27%, respectively, while the contribution of intermediate age stars disappears. Outside the nucleus, the blue end of the spectra are quite noisy, and we had to restrict the spectral range to the interval 4630–5230 Å, which we have labeled in the Table 1. In order to verify the effect of restricting the spectral range even further, we performed the synthesis of the SDSS spectrum also for this range. The result is a further increase of the old stellar population contribution to 59%, and a decrease of the young population and the FC component, while the intermediate age population continues to be absent.
In summary, using the wider spectral range we find that the intermediate age population contributes with approximately 13% of the flux at 4700 Å within the inner 3 arcsec of our PSQ. However, when we restrict the spectral range to the intervals labeled as and we lose the sensitivity to the intermediate age stellar population, most probably due to the fact that the main signatures of intermediate age stars (high order Balmer absorption lines) have been removed.
3.2 Synthesis of the IFU spectra
We have first constructed an “integrated” spectrum by adding the contribution of all spectra within the SDSS aperture (3). We have then performed the spectral synthesis in order to check the consistency of our results. These results are shown in the fourth line of Table 1 for the m wavelength interval, and they approximately agree with those for the SDSS, but presenting a difference of 10% in the population vectors and . According to the adev parameter, the percent difference between the model and the observed spectrum is 2% for both the SDSS and the IFU spectra, which corresponds to a good fit. However, because the stellar population synthesis method does not provide error estimates for the population vector, we can use the difference of 10% for the contribution of each component to the total light at 4700 Å as an error estimate for the population vector. In the 5 to the 12 line of Table 1 we show the results of the synthesis for the set of circular regions described in section 2 corresponding to apertures of 08 diameter. As discussed above, our IFU data spectra cover a limited spectral range, namely the ranges labeled as and . Only at the nucleus (region 1) we had enough SNR to cover the range, while in all other regions (2 to 7) we use the range .
In the nucleus (region 1), we have a predominance of the old population, contributing with 59%, while the intermediate age population contributes with 27% and the young + FC contribute with 14%. For the circumnuclear regions there is an increase in the contribution of the young populations, as we can see in the results for regions 2, 3, 4 and 5. The contribution of the old population decreases from 59% — at the nucleus — to 54%, 48%, 29% and 37% respectively, showing that the old population has the largest contribution at the nucleus. At region 6, only the intermediate age component is present. The young and intermediate age stellar populations are thus distributed preferentially in the circumnuclear region, at about 0.8 kpc from the central AGN.
Figure 2 shows the results of the synthesis for the SDSS spectrum and for five of the IFU integrated spectra of the circumnuclear region of the galaxy. The black dashed-dotted line, at the top of each panel, shows the observed spectrum. The red solid line shows the fitted model while the dotted line shows the masked regions of the spectrum. We have masked out emission lines and spurious features, such as the CCD gap. In each panel, at the bottom, we show the spectra corresponding to each age bin (identified as young, intermediate age and old) scaled to its percent contribution to the total continuum light at 4700 Å.
4 Emitting Gas
4.1 Measurement of the emission lines
The kinematics and the flux distribution of the emitting gas in the central region of active galaxies have typically been measured using a single Gaussian
function for each emission line profile. Nevertheless, a close inspection of the emission lines from the Narrow Line Region (NLR) of AGNs (e.g., Riffel et al., 2009a, 2010a)
shows that the profiles present deviations from a purely Gaussian profile, having “blue” or “red” wings, that are usually modeled by fitting a set of
two or more Gaussians. In order to take into account the presence of these asymmetries, we used instead of Gaussians the Gauss-Hermite series (hereafter, GH),
with which we have fitted the emission-line profiles of H and [O III ]
The GH series allows us to fit the wings of the emission lines via the and moments, in addition to the centroid velocities and the velocity dispersion. Following van der Marel & Franx (1993) and Cappellari & Emsellem (2004), the GH series can be written as
where is the amplitude of the GH series, is the central wavelength, is the velocity dispersion, are the GH moments and are the Hermite polynomials and the function is the standard Gaussian, given by
For any choice of free parameters of the Gaussian (, and ) there is a GH series described by the equation (3) that fits the emission line profile. A series that converge rapidly to the desired solution and that is very similar to the standard Gaussian function can be obtained by summing up to N=4 and taking the assumption that and (van der Marel & Franx, 1993). In this case, the GH series given by equation (3) can be approximated by
where and are the 3rd and 4th order Hermite polynomials.
The moment measures asymmetric deviations of a Gaussian profile, such as “blue” () or “red” () wings. The measures the symmetric deviations, where a wider and flatter profile than a Gaussian has and narrower and more peaked profile than a Gaussian has . In general, the use of the GH series allows a better fit to any emission-line profile.
The H and [O III ] emission-line profiles were measured using an IDL routine that solves for the best solution of parameters of the GH series using the nonlinear
least-squares Levenberg–Marquardt method. The IDL routine used is similar to that detailed in Riffel (2010b), which makes use of the routine MPFIT
4.2 Error estimates
Error estimates were obtained using Monte Carlo simulations adding to the original spectrum an artificial Gaussian noise normally-distributed, whose mean is equal to zero and the standard deviation equal to one. We have performed the Monte Carlo simulations with 200 realizations computing the best-fitting parameters of each realization so that we have the standard deviation of the best model in the end of the simulation. In principle, the larger the number of realizations, more reliable are the results, but in our case we have verified that with approximately 200 realizations the standard deviations are already of the same order as those obtained with 500 to 1000 realizations. In this way, we have obtained the error estimates for the free parameters of the GH function: the centroid velocity, velocity dispersion, the moments and and the flux of the emission line. The error values in each measured parameter are given together with the description of the corresponding maps in the following sections.
4.3 Emission line flux distributions
In Figure 3 we present the flux distributions in the H and [O III ]5007 emission lines, both given in logarithmic scale and erg cm s units. The location of the maximum brightness of the [O III ] flux distribution coincides — considering the spatial resolution of our data — with the position of the peak of the continuum. This position was determined as the centroid of the flux distribution in the continuum, obtained by collapsing the data cube between two continuum wavelengths. We have adopted this position as the galaxy nucleus and its location in the figures is identified either by two perpendicular dotted lines or a cross.
The H map shows a flux distribution elongated towards the South of the nucleus with a detached patch to the North. In the lower right corner there is another patch of emission which we attribute to an HII region located within what appears to be a spiral arm of the galaxy. The [O III ] flux distribution is also elongated along the North-South, but along the East-West is broader and more symmetrically distributed around the nucleus than H. Graham & Li (2009), using SDSS images, have reported the presence of a bar along the North-South direction, where we have reported the elongations in [O III ] and H flux distributions. The errors in the [O III ] emission line fluxes are 4% within the inner arcsecond, increasing to 13% at the borders of the mapped field. In the case of H the errors within the inner one arcsecond are 11% while at the borders of the mapped field they reach 25%.
In order to map the gas excitation, we show in Figure 4 the ratio between the [O III ] and H flux maps, whose distribution presents the highest ratios — and thus the highest excitation — at the nucleus and up to 1 to the North, roughly within the circle (with diameter of 08) labeled by the number 1. In this region the [O III ]/H ratio reaches a value of 20, which decreases steeply outwards, down to values as low as 3, which is typical of LINERs and may include also contribution of ionization by young stars. The lowest values of [O III ]/H are found at the regions labeled 4, 5 and 6, in regions where the stellar population synthesis results showed the largest contribution of young stellar populations.
5 Gas Kinematics
As described in the previous section, we have obtained the gas kinematics — centroid velocities, velocity dispersions and GH moments and — using the GH series to fit the emission lines H and [O III ]. In Figure 5 we present, to the left, the centroid velocity maps and, to the right, the corresponding velocity dispersion maps. The velocity maps are masked according to the error estimates, which are on average in the range km s. Pixels for which the error values are larger than these were masked. Centroid velocities are shown relative to the systemic velocity of the galaxy, 12480 km s, whose value was measured from the velocity channel maps in the H emission-line, as explained in section 5.1.
The [O III ] velocity field (top left panel of Fig. 5) shows positive values towards the South and negative values towards the North, suggesting a rotation pattern with the line of nodes approximately along the North–South direction and a velocity amplitude of 150 km s. Alternatively, this pattern could also be due to a bi-polar outflow, frequently observed in the Narrow-Line Region of AGN. Considering that the velocity pattern is distorted relative to that corresponding to simple rotation and it’s not centered at the nucleus, it is more likely that both rotation and outflows are present, as evidenced also by the blueshifts around the nucleus. The average errors in the velocity values increase from 6 km s within the inner arcsecond to 15 km s at the borders of the mapped field.
The [O III ] velocity dispersion map (top right panel of Fig. 5) presents a range of values from 130 km s to 300 km s. Within a radius of 05 from the nucleus we observe the lowest ’s, of 130 km s. This region is surrounded to the West by a semi-circle where the values are slightly higher, 175 km s, while to the SE, East and NE there is another semi-circle where the values reach even higher values of 250-300 km s. The errors and their distribution are similar to the ones for the [O III ] velocity field.
The H velocity field (bottom left panel of Fig. 5) also presents a pattern that suggests rotation plus outflows with similar orientation to that observed in [O III ], but reaching a higher velocity amplitude of 200 km s. As in the case of [O III ], we also observe a rotational pattern with kinematic center displaced from the continuum peak, showing blueshift at the nucleus. The average errors in the H velocity field range from 10 km s for the inner arcsecond to 25 km s at the borders of the mapped field.
The H velocity dispersion map (bottom right panel of Fig. 5) shows a smaller region of 02 near the the nucleus where the has values of 130 km s surrounded by a semi-circle to the West where the is larger than in the case of [O III ] (reaching 200 km s). Outwards, the values decrease to 100 km s. Although the H emission could not be measured as far as that of [O III ] towards the East, the small region where the could be measured also shows an increase there, similar to what we have observed for .
In the top panels of Fig. 6 we show the (left) and (right) GH moments of the [O III ] emission-line profiles. The moment, over most of the IFU field, presents values between 0 and , indicating that there isn’t significant asymmetric deviations with relation to Gaussian profiles. However, it reaches values of up to to the NE, indicating the presence of blue wings in this region, giving additional support to the presence of outflows. The moment shows only positive values, meaning that the profiles are narrower and more “pointy” than Gaussians. In the central region the moment has values around 0.2, increasing to 0.3 to the NE and decreasing to about 0.1 to the SW. The average errors obtained from the Monte Carlo simulations for the and are 0.02 within the inner arcsecond, increasing to 0.08 at the borders of the mapped field.
In the bottom panels of Fig. 6 we show the (left) and (right) GH moments of the H emission-line profiles. The moments present values of about around the nucleus, suggesting the presence of outflows towards us in this region, whereas in other regions the values are close to zero. The moment presents an irregular distribution, with values not deviating much from zero, indicating there are no significant symmetric deviations from a Gaussian profile for the H emission line. The average errors in and are larger than for [O III ], increasing from 0.04 within the inner arcsec to 0.2 towards the borders of the mapped field.
5.1 Channel maps
We have obtained channel maps by integrating the flux distribution in different slices of velocity along the [O III ] and H emission-line profiles. We have integrated the flux in each velocity channel after subtraction of the continuum contribution from both sides of the profile. This technique allows us to map the velocity field of the emitting gas throughout the emission line and not only at the central wavelength (like the velocity maps shown in Fig. 5).
In Figure 7, we present the channel maps of the H emission line for 16 velocities. The systemic velocity has been subtracted from the maps and its value was obtained as described below. The highest blueshifts, of km s, are observed at the nucleus. As the velocities become less negative (channel km s), it is possible to observe emission outside the nucleus towards the North. For velocities around zero, we can observe emission distributed over most parts of the field, including in the spiral arm to the West of the nucleus and in a similar structure to the East. As the velocities become more positive, the emission distribution moves to the South side. The highest positive velocities reach 290 km s to the South of the nucleus. The kinematics observed in the velocity channels are consistent with a rotation pattern with line of nodes approximately along North-South and an outflow observed in negative velocities within 1 arcsec from the nucleus. As the nuclear outflow disturbs the rotational pattern in the nuclear region in the centroid velocity maps (Fig. 5), we have used the H velocity channels to estimate the systemic velocity, observing that the flux distribution pattern becomes symmetric relative to the nucleus along the galaxy major axis for a systemic velocity of km s. A comparison between the panels at and km s, for example, (approximately symmetric relative to the zero velocity) indeed shows that the extranuclear emission is approximately symmetric relative to the nucleus (within 20 km s). In addition, most of the extranuclear emission along the minor axis is found at zero velocity, supporting our adopted systemic velocity.
In Figure 8, we present the channel maps of the [O III ] emission line in 25 different velocity channels, with an increase in velocity of
55 km s from one panel to another. The highest negative velocities, that reach km s — therefore higher than the ones observed
for H — are also observed at the nucleus, but extending circa 10 towards Northeast (NE) of the nucleus. As the velocities become less negative,
the emission becomes more extended toward the North. For velocities around zero we can observe the emission extending from North to South. As the velocities
become more positive, the emission moves towards South, in the same way as observed for H (up to 275 km s). For [O III ] even higher positive
velocities are observed, however the peak values remain centered on the nucleus.
6.1 Gas Excitation
We have investigated the gaseous excitation by mapping the ratio of [O III ]/H and by measuring its value in the integrated spectra of the circumnuclear region, as shown in column 6 of Table 2. Additionally, we have separated the synthesis results for the young stellar population in terms of ionizing () and not ionizing population (), given in percent contribution in columns 5-6 of Table 1.
Using as reference the BPT diagnostic diagrams (Baldwin, Phillips & Terlevich, 1981; Kauffmann et al., 2003; Kewley et al., 2006; Stasińska et al., 2006) we conclude that the emission-line ratios within the inner 0.3 kpc (radius) are typical of Seyfert nuclei with a mean value of [O III ]/H17.9, but reaching up to 25 at locations closest to the nucleus. In the circumnuclear region smaller values are observed, decreasing from 8.1 (at region 2) to 0.2 in the detached region 7. At these locations (regions 2-7), the synthesis shows that the young population is essentially comprised by ages younger than 5 Myrs. The BPT diagrams reveals additionally LINER-like emission-line ratios, which we have interpreted as due to excitation by the combined contribution of radiation from the AGN and hot stars.
6.2 Gas Kinematics
The gaseous velocity fields of [O III ] and H presented in Fig. 5 both show similar rotation patterns, but with the rotation displaced from the continuum peak position by circa 10 to the South. Such velocity field suggests that a possible rotation pattern around the nucleus is disturbed by a kinematic component not due to the gravitational potential of the bulge. In this region the emitting gas is blueshifted, with the velocity centroid reaching values of the order of km s. The velocities observed in the wings of the emission-line profiles are, nevertheless, much higher, reaching km s for [O III ] and km s for H, as observed in the channel maps of Figs. 8 and 7. These blueshifts support the presence of a nuclear outflow, in agreement also with the high velocity dispersions observed in this region (right panels of Fig. 5). Considering both the centroid velocities and the channel maps, the kinematic of the emitting gas can be understood as the combination of gas in rotation in the plane of the galaxy combined with an outflow component.
The orientation of the rotation axis is approximately along the North-South direction, although the true orientation is not clear due the presence of the outflow in the central region. The rotation in the H emission is seen in the channel maps between and km s, with blueshifts observed along the major axis to the North of the nucleus. These blueshifts decrease gradually as the emission shifts towards the South side of the major axis and the velocity increases to redshifts.
Additional support for the rotation in the [O III ] centroid velocity field is given by the Principal Component Analysis of the data cube, which is described in the next section. The Principal Component 4 (PC4) together with its tomogram seems to be consistent with the presence of a rotational component, as shown in the channel maps of the [O III ] emission-line.
6.3 PCA tomography
Astronomical observations with IFU’s provide a large amount of information, with two spatial dimensions and one spectral dimension, so that there are tens of millions of pixels to be analyzed. An efficient method to analyze all this information in an optimized way is the Principal Component Analysis (PCA), which has been adapted by Steiner et al. (2009) for application in IFU data. The PCA technique has been previously employed in astronomy for morphological classification of galaxies (Lahav et al., 1996) and for distinguishing radio-loud from radio-quiet QSOs (Boroson, 2002).
In order to extract additional information from our data and its complex gas kinematics, we have employed the PCA technique as described by Steiner et al. (2009), to the calibrated data cube. The PCA allows us to separate the information initially contained in a system of correlated coordinates into a new system of uncorrelated coordinates, whose components are ordered in accordance with decreasing values of the variance of the principal components. The new coordinates are called “eigenspectra” — i.e., eigenvectors given as a function of wavelength — so that their projection on the data produces images which are called “tomograms”. While the eigenspectra show the correlations and anti-correlations in the data, the tomograms will map them spatially. The tomograms are like slices of data in the space of eigenspectra, whose interpretation can reveal correlations or anti-correlations between spectral and spatial characteristics of the data.
In the panels of Figure 9 we present the results of the PCA applied to our data cube, where we show the principal components PC1, PC2 and PC4 (which show the largest variance of the data), represented through the corresponding tomograms and eigenspectra. We are not presenting the PC3, which has the third largest variance, because it contains mainly noise, mainly due to the difference in response between the blue and red slit detectors of the IFU. The other components do not present a significant variance. Correlations in the emission lines appear as all of them having positive values in the eigenspectrum (e.g., bottom left panel of Fig. 9). Anti-correlations appear as one line having positive value and the anti-correlated one having negative value (e.g., H negative and [O III ] positive in the bottom middle panel of Fig. 9).
In the left panels of Figure 9 we show tomogram 1 and its eigenspectrum (PC1), which represents the principal component of the PCA with circa 93.4% of the total variance. It shows a positive correlation between the emission lines, indicating that they are all generated in the same spatial region — the brighter region of the tomogram. Since the emission lines that are present in the eigenspectrum have [O III ]/H ratio typical of ionization by AGN, as well as contribution of the stellar population, we conclude that PC1 reveals the dominance of the AGN and the stellar population over the entire central region.
The central panels of Figure 9 show eigenspectrum and tomogram 2, for which there exists an inverse correlation between [O III ] and both the H emission and the stellar component, whose features are inverted and below zero in the eigenspectrum. The tomogram indeed shows that the corresponding emission originates in different regions: [O III ] dominates in the inner 03 (positive values), while the H and stellar population dominate in the surrounding region (negative values). The interpretation that we propose for PC2 is that although PC1 shows that most of the H emission is due to the AGN, there is some H emission which has a different origin. This emission is thus not generated by the AGN, but possibly by young stars that ionize the gas, indicating the presence of recent star formation around the nucleus. This conclusion is in agreement with the results of the stellar population synthesis, which show, in fact, the presence of young stars in these regions.
In the right panel of Figure 9 we show PC4, which represents 0.09% of the variance. The eigenspectrum 4 shows that the [O III ] emission has part of the profile in redshift and part in blueshift, and that they are anti-correlated, since the blueshift profile has negative values and the redshift profile has positive values; i.e. they originate in different regions, as seen in tomogram 4. The part in redshift originates in a compact structure to the Southwest of the nucleus which curves and extends to the South, while the blueshifted part originates in a compact region to the Northeast of the nucleus. This structure may be identified with the rotation component seen in the [O III ] channel maps.
6.4 Mass of the emitting gas
We have estimated the total mass of ionized gas of the inner region (excluding the detached spiral arm) as:
where is the electron density, is the proton mass, is the filling factor and is the volume of the emitting region.
The filling factor can be calculated using the H luminosity, as follows. The emissivity of the H line is given by (Osterbrock & Ferland, 2006):
where is the proton density, is the H effective recombination coefficient and is the corresponding energy. The total luminosity in H is given by integrating the emissivity over the entire emitting volume and over all directions, assuming a pure hydrogen gas which is completely ionized ():
where we have adopted the for a temperature of K. From the above equation, we obtain the product :
where is the H luminosity, in units of ergs s and is the electron density in units of cm. Thus, we can estimate the mass of the emitting region by inserting the expression 9 for above in 6, obtaining:
given in units of solar masses (Peterson, 1997).
The H luminosity was calculated from the integrated H flux F(H) of the inner region corrected for the reddening obtained from the SDSS spectrum (), using the reddening law of Cardelli, Clayton & Mathis (1989) adopting the theoretical ratio of 3.0, corresponding to case B recombination (Osterbrock & Ferland, 2006). For the assumed distance of 170 Mpc, we obtain ergs s. Using the [S II] ratio, obtained from the SDSS spectrum, the electron density was obtained by solving numerically the equilibrium equation for a -level atom using the IRAF routine stsdas.analysis.nebular.temden (De Robertis et al., 1987; Shaw & Dufour, 1994). The assumed electron temperature was K (Peterson, 1997), resulting in . With these assumptions, we obtain a mass of ionized gas of .
6.5 Mass outflow rate
We have concluded from the gas kinematics that the blueshifted emission in the nuclear region is due to outflowing gas. Although the geometry of the outflow is not clear from our maps, previous kinematic studies (e.g., Crenshaw & Kraemer, 2000; Das et al., 2005, 2006; Storchi-Bergmann et al., 2010) of other active galaxies show that AGN outflows present a conical geometry. On the basis of these previous studies, we assume that the outflow in PSQ J0210-0903 also has a conical geometry, with the cone axis directed approximately towards us. We thus consider that the outflowing gas is crossing the base of a cone whose radius we estimate from the extent of the blueshifted region around the nucleus seen in the channel maps. We have measured the ionized gas mass outflow rate at negative velocities through a circular cross section with radius 6 around the nucleus, assuming that the height of the cone is equal to the diameter of its base. The mass outflow rate can be calculated as (Riffel & Storchi-Bergmann, 2011a, b):
where cm is the area of the circular cross section, g is the proton mass, is the electron density and is the filling factor. The filling factor can be obtained from equation 9 above taking into account only the luminosity of the gas in outflow. We have considered that the outflowing gas corresponds to that observed in the velocity channels between and km s. The total flux was obtained by adding the flux of the channels centered between and km s and by adding 50% of the flux of the channels between and km s, for which the remaining emission – the other 50% – is assumed, as a best guess, to originate in the plane of the galaxy. The total flux in these channels is erg cm s and erg s. From the adopted geometry, the volume of the cone whose height and radius are such that pc is cm. The velocity adopted for the outflow is the average of the velocities of each channel above weighted by the corresponding flux, which results km s. Adopting cm (obtained from the [SII] line ratio from the SDSS spectrum), , we obtain M yr. For a somewhat larger gas density, of cm (as the aperture of the SDSS is too large, and the inner gas density being probably larger), we obtain and M yr.
Riffel & Storchi-Bergmann (2011b) have obtained M yr for Mrk 1157. They have also revised previous results for another 6 active galaxies in Barbosa et al. (2009) and in Riffel et al. (2009b), obtaining values for the mass outflow rate in range of M yr. Similar values, of M yr, have been reported by Storchi-Bergmann et al. (2010) for NGC 4151 as well as by Veilleux et al. (2005), with M yr for a set of luminous active galaxies. The mass outflow rate that we have estimated for PSQ J0210-0903 is thus in good agreement with previous estimates for other AGNs.
Now we can compare the outflow mass rate with the accretion rate necessary to feed the AGN, which can be calculated as (Peterson, 1997):
where is the nuclear bolometric luminosity, is the light speed and is the efficiency of conversion of the rest mass energy of the accreted material into radiation power and is bolometric luminosity in units of 10 erg s. Following Riffel & Storchi-Bergmann (2011a, b), can be approximated as , where is the H nuclear luminosity. Using the interstellar extinction coefficient , the intrinsic ratio (Osterbrock & Ferland, 2006) and the reddening law of Cardelli, Clayton & Mathis (1989), we have obtained erg cm s within 04 of the nucleus. At the distance of 170 Mpc, erg s and the corresponding nuclear bolometric luminosity is erg s. If we assume an efficiency , for an optically thick and geometrically thin accretion disc (Frank, King & Raine, 2002) we derive an accretion rate of yr. Since is 2 orders of magnitude higher than , we conclude that the outflowing gas does not originate only from the AGN, but is composed mostly of interstellar gas from the surrounding region of the galaxy, which is swept away by the AGN outflow.
We can additionally estimate the kinetic power of the outflow, considering both the radial and turbulent component of the the gas motion. Using the H kinematics to estimate this power, we have
where is the velocity of the outflowing gas and is the velocity dispersion. From Fig. 5 we have km s and adopting km s and M yr, as discussed above, we obtain erg s, which is , so that between 0.03% and 0.1% of the mass accretion rate is transformed in kinetic power in the outflow. These values are one order of magnitude smaller than the AGN feedback derived by Di Matteo et al. (2005) in simulations to account for the co-evolution of black holes and galaxies and approximately at the same order of those estimated by Barbosa et al. (2009) and Storchi-Bergmann et al. (2010) using similar IFU data.
|Pos.||ap.||log||log||[O III ]/H|
6.6 Stellar Population
The stellar population synthesis results are affected by the restrict wavelength range of our data, which excludes the near-UV region (below 4200 Å) – the most sensitive to the intermediate age population. This limitation may have affected the absolute value of the contribution of this age component. In order to compensate for this, we have restricted the population vector to only 3 stellar components (young, intermediate age and old). The relative contribution of these different age components at the different locations seems to be robust, supported, for example, by the dominance of the young stellar population component in region 7, for which the emission-line spectrum supports ionization by stars.
From the synthesis, we conclude that within the inner 0.3 kpc the old population dominates the flux at 4700 Å ( 60%), but that there is also some contribution of the intermediate age one ( 30%). The gaseous kinematics has shown that the AGN feedback is concentrated approximately within the inner 0.3 kpc, coinciding with the region where the stellar population synthesis has revealed predominance of old population, but with some contribution of intermediate age population and FC. Beyond this region, at the 0.8 kpc ring, the young stellar population has a contribution ranging from 40% to 60%, reaching 100% at the detached region 7, located in an spiral arm. In regions 4, 5 and 6, there is also significant contribution of the intermediate age population.
We can represent the mixture of ages of the stellar population components in a more condensed form by its mean stellar age. Following Cid Fernandes et al. (2005) we can define the mean stellar age weighted by the light fraction as
and by the stellar mass fraction as
While the light-weighted mean stellar age is biased towards younger ages, as we can see in Table 2, the mass-weighted stellar age is more representative of the older population (which has a larger mass-to-light ratio). According to Cid Fernandes et al. (2005), although the latter is more physically representative, it has less relation with the observed spectrum. These authors argue that, in practice, is the most useful of the two indices, due to the large M/L ratio of old stars. In Table 2 we list these two mean ages for the integrated spectrum and for the regions . is indeed characteristic of a post-starburst population at most locations, where mean ages are smaller than 1 Gyr, except at the nucleus, where the mean age is larger than 3 Gyr.
7 Summary and Conclusions
We have used integral field optical spectroscopy in order to map the stellar population and the gas kinematics within the inner 1.5 kpc (radius) around the Post-Starburst Quasar J0212-0903 at spatial resolution of 0.5 kpc (0.6 arcsec) and velocity resolution of 120 km s. This is the first two-dimensional study of both the stellar population and gas kinematics of a Post-Starburst Quasar and our main conclusions are:
The stellar population is dominated by old stars within 0.3 kpc (radius) from the nucleus, while in the circumnuclear region both intermediate age (100 Myrs 2.5G Gyrs) and young stars ( Myrs) dominate the optical flux. Our results support a location of most of the post-starburst population in a ring with radius of 0.8 kpc, where active star formation is also occurring;
Extended emission up to 1.5 kpc from the nucleus is observed in both [O III ] and H. The [O III ] flux distribution is more centrally concentrated, while the H flux distribution is more extended approximately along the line of the nodes and is also observed in a detached patch which seems to belong to a spiral arm with active star formation;
The emission-line ratios are typical of Seyfert excitation within the inner 0.3 kpc, where an outflow is observed, while beyond this region the line ratios are typical of LINERs, which we attribute to a combination of diluted nuclear radiation and ionization by young stars. A Principal Component Analysis (PCA) supports the conclusion that recent star formation dominates the spectrum in the circumnuclear region;
The mass of ionized gas in the inner 1.25 kpc (radius) is ;
The gas kinematics can be reproduced by a combination of (1) rotation in the plane of the galaxy whose line of the nodes runs approximately North-South with amplitude km s (supported also by the PCA) and (2) an outflow observed within the inner 0.3 kpc, with blueshifts reaching up to km s;
The mass outflow rate is in the range of 0.3-1.1 M yr, which is times the AGN mass accretion rate of yr, implying that most of the outflow originates via mass-loading in the surrounding interstellar medium of the galaxy, swept away by the AGN outflow.
We conclude that our observations support both the evolutionary and quenching scenario for this galaxy, as follows. The feeding of gas to the nuclear region has triggered a circumnuclear starburst a few 100’s Myr ago, extending all the way to the nucleus. The remaining gas from this inflow, combined with the mass loss from the newly formed stars may have then triggered the nuclear activity, producing the observed gas outflow. This outflow, observed within the inner 0.3 kpc, has then quenched further star formation at this location, in agreement with the observed absence of young stars and the contribution of intermediate age stars in the inner 0.3 kpc. Beyond this region affected by the outflow, star-formation seems to be still active in a 0.8 kpc ring. The presence of a delay between the triggering of star formation and the nuclear activity is supported by the presence of the intermediate age stellar population at the nucleus.
We thank an anonymous referee for the valuable suggestions which helped to improve the paper. We thank Sabrina Lyn Cales for valuable suggestions which helped to improve the present paper and Rajib Gangulythe who has helped with the Gemini Proposal. Based on observations obtained at the Gemini Observatory, which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the NSF on behalf of the Gemini partnership: the National Science Foundation (United States), the Science and Technology Facilities Council (United Kingdom), the National Research Council (Canada), CONICYT (Chile), the Australian Research Council (Australia), Ministério da Cin̂cia, Tecnologia e Inovação (Brazil) and Ministerio de Ciencia, Tecnología e Innovación Productiva (Argentina). This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. The STARLIGHT project is supported by the Brazilian agencies CNPq, CAPES and FAPESP, and by the France-Brazil CAPES/Cofecub program. MSB wishes to thank Brazilian agencie CAPES for the visiting professor Fellowship. This work has been partially supported by the Brazilian institution CNPq.
- pagerange: 2D stellar population and gas kinematics of the inner 1.5 kpc of the post-starburst quasar SDSS J0210-0903–2D stellar population and gas kinematics of the inner 1.5 kpc of the post-starburst quasar SDSS J0210-0903
- pubyear: 2012
- The NASA/IPAC Extragalactic Database (NED) is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.
- IRAF is distributed by the National Optical Astronomy Observatories, which is operated by the Association of Universities for Research in Astronomy, Inc. (AURA) under cooperative agreement with the National Science Foundation.
- We denote the [O III ] emission line simply by [O III ]
- MPFIT provides a robust and relatively rapid way to run the least squares method and is available in http://purl.com/net/mpfit as part of Markwardt IDL Library
- Allington-Smith J., Murray G., Content R., Dodsworth G., Davies R., Miller B.W., Jorgensen I., Hook I. et al., 2002, PASP, 114, 892
- Baldwin J. A., Phillips M. M., Terlevich R., 1981, PASP, 93, 5
- Barbosa F. K. B., Storchi-Bergmann T., Cid Fernandes R., Winge C., Schmitt H., 2009, MNRAS, 396, 2
- Boroson, T.A., 2002, ApJ, 565, 78
- Bruzual G., Charlot S., 2003, MNRAS, 344, 1000
- Brotherton M.S., van Breugel W., Stanford S.A., Smith R.J., Boyle B.J., Miller L., Shanks T., Croom S.M., Filippenko A.V., 1999, ApJL, 520, 87
- Brotherton M. S., Grabelsky M., Canalizo G., van Breugel W., Filippenko A. V., Croom S., Boyle B., Shanks T., 2002, PASP, 114, 593
- Cales S.L., Brotherton M. S., Shang Z., Bennert V. N., Canalizo G., Stoll R., Ganguly R., Vanden Berk D., Paul C., Diamond-Stanic A., 2011, ApJ, 741, 106
- Canalizo G., Stockton A., Brotherton M.S., van Breugel W., 2000, AJ, 119, 59
- Cano-Díaz M., Maiolino R., Marconi A., Netzer H., Shemmer O., Cresci G., 2012, A&A, 537, 8
- Cappellari M. and Emsellem E., 2004, PASP, 116, 138
- Cardelli J. A., Clayton G.C., Mathis J.S., 1989, ApJ, 345, 245
- Chabrier G., 2003, PASP, 115, 763
- Cid Fernandes R., Sodré L., Schmitt H. R., Leão J. R. S., 2001, MNRAS, 325, 60
- Cid Fernandes R., Gu, Q. Melnick J., Terlevich E., Terlevich R., Kunth D., Rodrigues L.R., Joguet B., 2004, MNRAS, 355, 273
- Cid Fernandes R., Mateus A., Sodré L., Stasińska G., Gomes J.M., 2005, MNRAS, 358, 363
- Cid Fernandes R., Schoenell W., Gomes J. M., Asari N. V., Schlickmann M., Mateus A., Stasińska G., Sodré Jr. L., Torres-Papaqui J. P. and Seagal Collaboration, 2009, RMXAC, 35, 127
- Crenshaw, D.M. and Kraemer S.B., 2000, ApJ, 532, 101
- Das V., Crenshaw D.M., Hutchings J.B., Deo R.P., Kraemer S.B., Gull T.R., Kaiser M.E., Nelson C.H., Weistrop D., 2005, AJ, 130, 945
- Das V., Crenshaw D.M., Kraemer S.B., Deo R.P., 2006, AJ, 132, 620
- Davies R. I., Müller Sánchez F., Genzel R., Tacconi L. J., Hicks E. K. S., Friedrich S., Sternberg A., 2007, ApJ, 671, 1388
- Di Matteo T., Springel V., Hernquist L., 2005, Nat, 433, 604
- De Robertis M.M., Dufour R.J., Hunt R.W., 1987, JRASC, 81, 195
- Ferrarese L., Ford H., 2005, Space Sci. Rev., 116, 523
- Frank J., King A., Raine, D.J., 2002, Accretion Power in Astrophysics, 3rd ed., Cambridge Univ. Press, Cambridge
- Gebhardt K., Bender R., Bower G., Dressler A., Faber S. M., Filippenko A. V., Green R., Grillmair C., Ho L. C., Kormendy J., Lauer T. R., Magorrian J., Pinkney J., Richstone D., Tremaine S., 200, ApJ, 539, 13
- Gebhardt K., Kormendy J., Ho L. C., Bender R., Bower G., Dressler A., Faber S. M., Filippenko A. V., Green R., Grillmair C., Lauer T. R., Magorrian J., Pinkney J., Richstone D. Tremaine S., 2000, ApJ, 543, 5
- Graham A.W., Li I-hui., 2009, ApJ, 698, 812
- Graham A. W., Onken C. A., Athanassoula E., Combes F., 2011, MNRAS, 412, 2211
- Granato G. L., De Zotti G., Silva L., Bressan A., Danese L., 2004, ApJ, 600, 580
- Gültekin K., Richstone D. O., Gebhardt K., Lauer T. R., Tremaine S., Aller M. C., Bender R., Dressler A., Faber S. M., Filippenko A. V., Green R., Ho L. C., Kormendy J., Magorrian J., Pinkney J., Siopis C., ApJ, 698, 198
- Heckman T.M., Kauffmann G., Brinchmann J., Charlot S., Tremonti C., White S. D. M., 2004, 613, 109
- Hopkins P. F., Hernquist L., Cox T. J., Robertson B., Springel V., 2006, ApJS, 163, 50
- Kauffmann G., Heckman T. M., Tremonti C., Brinchmann J., Charlot S., White S. D. M., Ridgway S. E., Brinkmann J., Fukugita M., Hall P. B., Ivezić Ž., Richards G. T., Schneider D. P., 2003, MNRAS, 346, 1055
- Kewley L. J., Groves B., Kauffmann G., Heckman T., 2006, MNRAS, 372, 961
- Lahav O., Naim A., Sodré Jr.L., Storrie-Lombardi M.C., 1996, MNRAS, 283, 207
- Le Borgne J.F., Bruzual G., Pelló R., Lançon A., Rocca-Volmerange B., Sanahuja B., Schaerer D., Soubiran C., Vílchez-Gómez R., 2003, A&A, 402, 433
- Markwardt, C.B., 2009, ASPCS, 411, 251
- Merritt D. and Ferrarese L., 2001, MNRAS, 320, 30
- Osterbrock D.E. and Ferland G.J., 2006, Astrophysics of Gaseous Nebulae and Active Galactic Nuclei, 2nd. ed., University Science Books, California
- Osterbrock Peterson B.M., 1997, An Introduction to Active Galactic Nuclei, Cambridge Univ. Press, Cambridge
- Riffel R.A., Storchi-Bergmann T., McGregor, P.J., 2009, ApJ, 698, 1767
- Riffel R. A., Storchi-Bergmann T., Dors O. L., Winge C., 2009, MNRAS, 393, 783
- Riffel R.A., Storchi-Bergmann T., Nagar N.M., 2010, MNRAS, 404, 166
- Riffel, R.A., 2010, APSS, 327, 239
- Riffel R.A., Storchi-Bergmann T., 2011, MNRAS, 411, 469
- Riffel R.A., Storchi-Bergmann T., 2011, MNRAS, 417, 275
- Schlegel D. J., and Finkbeiner D. P., Davis, M., 1998, ApJ, 500, 525
- Shen Y. et al., 2011, ApJS 194, 45
- Stasińska G., Cid Fernandes R., Mateus A., Sodré L., Asari N. V., MNRAS, 2006, 371, 972
- Steiner J.E., Menezes R.B., Ricci, T.V., Oliveira, A.S., 2009, MNRAS, 395, 64
- Shaw R.A. and Dufour R.J., 1994, ASPC Ser. 61: Astronomical Data Analysis Software and Systems III, 61, 327
- Storchi-Bergmann T., González Delgado R.M., Schmitt, H.R., Cid Fernandes R., Heckman T., 2001, ApJ, 559, 147
- Storchi-Bergmann T., Lopes R.D.S., McGregor P.J., Riffel R.A., Beck T., Martini P., 2010, MNRAS, 402, 819
- Tremaine S., Gebhardt K., Bender R., Bower G., Dressler A., Faber S. M., Filippenko A. V., Green R., Grillmair C., Ho L. C., Kormendy J., Lauer T. R., Magorrian J., Pinkney J., Richstone, D., 2002, ApJ, 574, 740
- van der Marel R.P. and Franx M., 1993, ApJ, 407, 525
- van Dokkum, P. G., 2001, PASP, 113, 1420
- Veilleux S., Cecil G., Bland-Hawthorn J., 2005, AR&A, 43, 769