Äîêóìåíò âçÿò èç êýøà ïîèñêîâîé ìàøèíû. Àäðåñ îðèãèíàëüíîãî äîêóìåíòà : http://www.stsci.edu/~marel/abstracts/psdir/2003ApJ...586..868M.ps Äàòà èçìåíåíèÿ: Sun Mar 30 23:31:02 2003 Äàòà èíäåêñèðîâàíèÿ: Sat Dec 22 18:35:02 2007 Êîäèðîâêà: Ïîèñêîâûå ñëîâà: ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï ï

IS THERE REALLY A BLACK HOLE AT THE CENTER OF NGC 4041?
CONSTRAINTS FROM GAS KINEMATICS 1
A. Marconi, 2 D. J. Axon, 3 A. Capetti, 4 W. Maciejewski, 2,5 J. Atkinson, 3 D. Batcheldor, 3
J. Binney, 6 M. Carollo, 7 L. Dressel, 8 H. Ford, 9 J. Gerssen, 8 M. A. Hughes, 3
D. Macchetto, 8,10 M. R. Merrifield, 11 C. Scarlata, 8 W. Sparks, 8
M. Stiavelli, 8 Z. Tsvetanov, 9 and R. P. van der Marel 8
Received 2002 July 19; accepted 2002 November 27
ABSTRACT
We present Space Telescope Imaging Spectrograph spectra of the Sbc spiral galaxy NGC 4041, which were
used to map the velocity field of the gas in its nuclear region. We detect the presence of a compact
(r ' 0>4 ' 40 pc), high surface brightness, rotating nuclear disk cospatial with a nuclear star cluster. The disk
is characterized by a rotation curve with a peak­to­peak amplitude of #40 km s #1 and is systematically
blueshifted by #10--20 km s #1 with respect to the galaxy systemic velocity. With the standard assumption of
constant mass­to­light ratio and with the nuclear disk inclination taken from the outer disk, we find that a
dark point mass of Ï1 ×0:6
#0:7 ÷ # 10 7 M # is needed to reproduce the observed rotation curve. However, the
observed blueshift suggests the possibility that the nuclear disk could be dynamically decoupled. Following
this line of reasoning, we relax the standard assumptions and find that the kinematical data can be accounted
for by the stellar mass provided that either the central mass­to­light ratio is increased by a factor of #2 or the
inclination is allowed to vary. This model results in a 3 # upper limit of 6 # 10 6 M # on the mass of any nuclear
black hole (BH). Overall, our analysis only allows us to set an upper limit of 2 # 10 7 M # on the mass of the
nuclear BH. If this upper limit is taken in conjunction with an estimated bulge Bmagnitude of #17.7 and with
a central stellar velocity dispersion of '95 km s #1 , then these results are not inconsistent with both the M BH ­
L sph and the M BH ­# * correlations. Constraints on BH masses in spiral galaxies of types as late as Sbc are still
very scarce; therefore, the present result adds an important new data point to our understanding of BH
demography.
Subject headings: black hole physics --- galaxies: individual (NGC 4041) ---
galaxies: kinematics and dynamics --- galaxies: nuclei --- galaxies: spiral
1. INTRODUCTION
It has long been suspected that the most luminous active
galactic nuclei (AGNs) are powered by accretion of matter
onto massive black holes (BHs; e.g., Lynden­Bell 1969). This
belief, combined with the observed evolution of the space den­
sity of AGNs (Sootan 1982; Chokshi & Turner 1992; Marconi
& Salvati 2002) and the high incidence of low­luminosity
nuclear activity in nearby galaxies (Ho, Filippenko, & Sargent
1997), implies that a significant fraction of luminous galaxies
must host BHs of mass 10 6 --10 10 M # .
It is now clear that a large fraction of hot spheroids (E--
S0) contain a BH (Harms et al. 1994; Kormendy & Rich­
stone 1995; Kormendy et al. 1996; Macchetto et al. 1997;
van der Marel & van den Bosch 1998; Bower et al. 1998;
Marconi et al. 2001) with mass proportional to the mass (or
luminosity) of the host spheroid (M BH =M sph # 0:001; e.g.,
Merritt & Ferrarese 2001a). Recently, Ferrarese & Merritt
(2000) and Gebhardt et al. (2000) have shown that a tighter
correlation holds between the BH mass and the velocity dis­
persion of the bulge. Clearly, any correlation of BH and
spheroid properties would have important implications for
theories of galaxy formation in general and bulge formation
in particular. However, to date, there are very few secure
BH measurements or upper limits in spiral galaxies even
though we know that AGNs are common in such systems
(Maiolino & Rieke 1995). In total, there are 37 secure BH
detections according to Kormendy & Gebhardt (2001), or
just 22 according to Merritt & Ferrarese (2001b), depending
on the definition of `` secure.'' Only #20% of these BH
detections (7/37 or 4/22, respectively) are in galaxy types
later than S0, and only three in Sbc types and later (the
Milky Way, Genzel et al. 2000; NGC 4258, Miyoshi et al.
1995; NGC 4945, Greenhill, Moran, & Herrnstein 1997). It
is therefore important to directly establish how common
BHs are in spiral galaxies and if they follow the same M BH ­
M sph and M BH ­# * correlations as elliptical galaxies.
1 Based on observations made with the NASA/ESA Hubble Space Tele­
scope, obtained at the Space Telescope Science Institute, which is operated
by the Association of Universities for Research in Astronomy, Inc., under
NASA contract NAS 5­26555. These observations are associated with
proposal 8228.
2 INAF, Osservatorio Astrofisico di Arcetri, Largo Enrico Fermi 5,
I­50125 Florence, Italy.
3 Department of Physical Sciences, University of Hertfordshire, Hatfield
AL10 9AB, UK.
4 INAF, Osservatorio Astronomico di Torino, Strada Osservatorio 20,
10025 Pino Torinese, Torino, Italy.
5 Also at Obserwatorium Astronomiczne Uniwersytetu Jagiellon ’ skiego,
Poland.
6 Theoretical Physics, University of Oxford, 1 Keble Road, Oxford
OX1 3NP, UK.
7 Eidgenoessische Technische Hochschule Zuerich, Hoenggerberg HPF
G4.3, CH­8092 Zurich, Switzerland.
8 Space Telescope Science Institute, 3700 San Martin Drive, Baltimore,
MD 21218.
9 Department of Physics and Astronomy, Johns Hopkins University,
3400 North Charles Street, Baltimore, MD 21218.
10 ESA Space Telescopes Division.
11 School of Physics and Astronomy, University of Nottingham,
Nottingham NG7 2RD, UK.
The Astrophysical Journal, 586:868--890, 2003 April 1
# 2003. The American Astronomical Society. All rights reserved. Printed in U.S.A.
868

This can be achieved only with a comprehensive survey
for BHs that covers quiescent and active spiral galaxies of
all Hubble types. Such a survey would pin down the mass
function and space density of BHs and their connection with
host galaxy properties (bulge mass, disk mass, etc.).
To detect BHs, one requires spectral information at the
highest possible angular resolution: the `` sphere of influ­
ence '' (Bahcall &Wolf 1976) of BHs is typically #1 00 even in
the most nearby galaxies. Nuclear absorption line spectra
can be used to demonstrate the presence of a BH (Kor­
mendy & Richstone 1995; Richstone 1998; van der Marel et
al. 1998), but the interpretation of the data is complex
because it involves stellar dynamical models that have many
degrees of freedom, which can be pinned down only when
data of very high signal­to­noise ratio (S/N) are available
(Binney & Mamon 1982; Statler 1987; Merritt 1997; Binney
&Merrifield 1998). Radio­frequency measurements of mas­
ers in disks around BHs provide some of the most spectacu­
lar evidence for BHs but have the disadvantage that only a
small fraction of the disks will be inclined such that their
maser emission is directed toward us (Braatz, Wilson, &
Henkel 1997). Studies of ordinary optical emission lines
from gas disks provide an alternative and relatively simple
method to detect BHs.
Hubble Space Telescope (HST) studies have discovered
many cases of such gas disks in early­type galaxies (M87,
Ford et al. 1994; NGC 4261, Ja#e et al. 1996; NGC 5322,
Carollo et al. 1997; Cen A, Schreier et al. 1998) and have
demonstrated that both their rotation curves and line pro­
files are consistent with thin disks in Keplerian motion
(Ferrarese, Ford, & Ja#e 1996; Macchetto et al. 1997; Ford
et al. 1998; van der Marel & van den Bosch 1998; Bower et
al. 1998; Marconi et al. 2001).
In early­type galaxies there are still worrisome issues
about the dynamical configuration of nuclear gas (misalign­
ment with the major axis, irregular structure, etc.). By con­
trast, nuclear gas in relatively quiescent spirals is believed to
be organized into well­defined rotating disks seen in optical
line images (e.g., M81; Devereux, Ford, & Jacoby 1997). Ho
et al. (2002) recently found that the majority of spiral gal­
axies in their survey have irregular velocity fields in the
nuclear gas, not well suited for kinematical analysis. Still,
25% of the galaxies where H# emission was detected all the
way to the center have velocity curves consistent with circu­
lar rotation, and the galaxies with more complicated veloc­
ity curves can also be useful for BHmass measurement after
detailed analysis of the spectra. Indeed, even in the most
powerful Seyfert nuclei such as NGC 4151, where the gas is
known to be interacting with radio ejecta, it may be possible
to get the mass of the BH from spatially resolved HST spec­
troscopy by careful analysis of the velocity field to separate
the underlying quiescently rotating disk gas from that dis­
turbed by the jets (Winge et al. 1999).
Prompted by these considerations, we have undertaken a
spectroscopic survey of 54 spirals using the Space Telescope
Imaging Spectrograph (STIS) on HST. Our sample was
extracted from a comprehensive ground­based study by
D. J. Axon et al. (2003, in preparation), who obtained H#
and [N ii] rotation curves, at a seeing­limited resolution of
1 00 , of 128 Sb, SBb, Sc, and SBc spiral galaxies from RC3.
By restricting ourselves to galaxies with recession velocities
V < 2000 km s #1 , we obtained a volume­limited sample of
54 spirals that are known to have nuclear gas disks and span
wide ranges in bulge mass and concentration. The systemic
velocity cuto# was chosen so that we can probe close to the
nuclei of these galaxies and detect even lower mass BHs.
The frequency of AGNs in our sample is typical of that
found in other surveys of nearby spirals, with comparable
numbers of weak nuclear radio sources and LINERs. The
sample is described in detail by D. J. Axon et al. (2003, in
preparation).
This paper presents the observations of NGC 4041, the
first object observed, and a detailed description of the analy­
sis and modeling techniques that will be applied to the other
galaxies in the sample. From the Lyon/Meudon Extra­
galactic Database 12 (LEDA), NGC 4041 is classified as an
Sbc spiral galaxy with no detected AGN activity. Its average
heliocentric radial velocity from radio measurements is
1227 # 9 km s #1 , becoming '1480 km s #1 after correction
for Local Group infall onto Virgo. With H 0 ¼ 75 km s #1
Mpc #1 this corresponds to a distance of '19.5 Mpc and a
scale of 95 pc arcsec #1 .
The outline of the paper is as follows. In x 2 we present
the adopted observational strategy and data reduction tech­
niques. In x 3 we present the rotation curves of the ionized
gas and the broadband images of the nuclear region of the
galaxy. In x 4.1 we derive the stellar luminosity density from
the observed surface brightness distribution, which is then
used in the model fitting of the kinematical data. All the
details of the inversion procedure are described in Appendix
A. The model fitting of the kinematical data is described in
x 4.2, and all the details of the model computation are
described in Appendix B. In particular, x 4.2.2 describes the
standard approach, while in x 4.2.3 an alternative approach
is considered. In x 5 we discuss the e#ects of our assumptions
on the derived value of the BH mass for which, in the
present case, only an upper limit can be set. We then com­
pare this upper limit with the M BH ­L sph and M BH ­# * corre­
lations. Finally, our conclusions are presented in x 6.
2. OBSERVATIONS AND DATA REDUCTION
2.1. STIS Observations
NGC 4041 was observed with STIS on HST in 1999 July
7. An #5 00 # 5 00 acquisition image was obtained with the
F28X50LP filter, and the galaxy nucleus, present within the
field of view, was subsequently centered and reimaged fol­
lowing the ACQ procedure. The exposure time of the
acquisition images was 120 s.
The observational strategy consisted in obtaining spectra
at three parallel positions with the central slit centered on
the nucleus and the flanking ones at a distance of 0>2. The
slit positions are overlaid on the acquisition image in
Figure 1, and their position angle (P.A.) is 43 # . At each slit
position we obtained two spectra with the G750M grating
centered at H#, with the second spectrum shifted along the
slit by an integer number of detector pixels in order to
remove cosmic­ray hits and hot pixels. The nuclear spec­
trum (NUC) was obtained with the 0>1 slit and no binning
of the detector pixels, yielding a spatial scale of 0>0507
pixel #1 along the slit, a dispersion of D# ¼ 0:554 A š pixel #1 ,
and a spectral resolution of R ¼ #=Ï2D#÷ ' 6000. The o#­
nuclear spectra (POS1 and POS2) were obtained with the
0>2 slit and 2 # 2 on­chip binning of the detector pixels,
yielding 0>101 pixel #1 along the slit, 1.108 A š pixel #1 along
12 Available at http://leda.univ­lyon1.fr.
BLACK HOLE AT CENTER OF NGC 4041 869

the dispersion direction, and R ' 3000. Total exposure
times were 950 s for the NUC position and 420 and 500 s for
POS1 and POS2, respectively.
The acquisition images were flat­fielded, realigned, and
co­added in order to improve the S/N. The pixel scale
is 0>0507. The flux calibration was first obtained using
the PHOTFLAM header in the image. We subsequently
applied a color correction to convert to Johnson R magni­
tudes. In order to do this, we used the average spectra of Sb
and Sc spiral galaxies from Kinney et al. (1996) as spectral
templates.
The raw spectra were reprocessed through the calstis
pipeline using the darks obtained daily for STIS. Standard
pipeline tasks were used to obtain flat­field--corrected
images. The two exposures taken at a given slit position
were then realigned with a shift along the slit direction (by
an integer number of pixels), and the pipeline task ocrreject
was used to reject cosmic rays and hot pixels. Subsequent
calibration procedures followed the standard pipeline
reduction described in the STIS Instrument Handbook
(Leitherer et al. 2001); i.e., the spectra were wavelength cali­
brated and corrected for two­dimensional distortions. The
expected accuracy of the wavelength calibration is 0.1--0.3
pixels within a single exposure and 0.2--0.5 pixels among
di#erent exposures (Leitherer et al. 2001), which converts
into #3--8 km s #1 (relative) and #5--13 km s #1 (absolute).
The relative error on the wavelength calibration is negligible
for the data presented here because our analysis is restricted
to the small detector region including H# and [N ii]
(D# < 100 A š ).
The nominal slit positions obtained as a result of the STIS
ACQ procedure were checked by matching the light profiles
measured along the slit with the synthetic ones derived from
the acquisition image: we collapsed the spectra along the
dispersion direction and compared the resulting light pro­
files with the ones extracted from the acquisition image for a
given slit position. The agreement is good for all slit posi­
tions and the center of the NUC slit if o#set by only #0>03
with respect to the position of the target determined by the
STIS ACQ procedure (0, 0 position in Fig. 1).
We selected the spectral regions containing the lines of
interest and subtracted the continuum by fitting a linear pol­
ynomial row by row along the dispersion direction. The
continuum­subtracted lines were then fitted row by row
with Gaussian functions using the task LONGSLIT in the
TWODSPEC FIGARO package (Wilkins & Axon 1992)
and the task specfit in the IRAF 13 STSDAS package. When
the S/N was insu#cient (thin line), the fitting was improved
by co­adding two or more pixels along the slit direction.
2.2. WFPC2 Images
WFPC2 images in the F450W (#B), F606W (#R), and
F814W (#I ) filters were retrieved from the archive. These
images encompass the entire galaxy with the nucleus located
in the WF3 chip. The data were automatically reprocessed
with the best calibration files available before retrieval. Two
exposures were performed in each of the three filters, in
order to remove cosmic rays. Warm pixels and cosmic rays
were removed (STSDAS tasks warmpix and crrej) and
mosaicked images with spatial sampling of 0>1 pixel #1 were
obtained using the wmosaic task, which also corrects for the
optical distortions in the four WFPC2 chips. The back­
ground was estimated from areas external to the galaxy
where no emission is detected. Flux calibration to Vega
magnitudes was performed using the zero points by
Whitmore (1995). To convert to standard filters in the John­
son­Cousins system, we estimated the color correction using
various spectral templates, A0 V and K0 V stars, and the Sb
and Sc spiral spectra from Kinney et al. (1996). The color
corrections are the following: I#F814W ' #0:1 (Sb, K0 V)
and #0.005 (A0 V); R#F606W ' #0:4 (Sb, K0 V) and
#0.05 (Sc, A0 V); and B#F450W ' 0:1 (Sb, K0 V) and 0.0
(Sc, A0 V). With these corrections, the di#erences between
colors from Johnson­Cousins and HST instrumental mag­
nitudes are DÏR#I÷ ' #0:3, #0.04; DÏB#I÷ ' 0:2, 0.0; and
DÏB#R÷ ' 0:5, 0.0 according to the templates used. Unless
stated otherwise, we have applied the color corrections
derived from the K0 V star and Sb spiral templates, which
are most suitable to the present data.
2.3. Ground­based Observations
NGC 4041 was also observed with the Near Infrared
Camera Spectrometer (NICS; Ba#a et al. 2001) at the Tele­
scopio Nazionale Galileo (TNG) in 2001 February 12 using
the K 0 filter and the small­field camera that yields a 0>13
pixel size. The night was not photometric and seeing during
the observations was #0>8. Observations consisted of sev­
eral exposures, each with the object at a di#erent position
on the array. Data reduction consisted of flat­fielding and
sky subtraction. The frames were then combined into a
Fig. 1.---Slit positions overlaid on the acquisition image. The 0, 0 posi­
tion is the position of the target derived from the STIS ACQ procedure.
The white cross is the kinematic center derived from the fitting of the rota­
tion curves (see x 4.2).
13 IRAF is distributed by the National Optical Astronomy Observ­
atories, which are operated by the Association of Universities for Research
in Astronomy, Inc., under cooperative agreement with the National Science
Foundation.
870 MARCONI ET AL. Vol. 586

mosaic. In order to flux­calibrate the NICS image, we have
used the Two Micron All Sky Survey (2MASS) image avail­
able on the World Wide Web. The 2MASS K­band image is
flux­calibrated to an accuracy of less than 0.1 mag and has a
spatial resolution of #3>5 with a pixel size of 1 00 . We have
therefore degraded the spatial resolution of the NICS image
to the 2MASS image and rebinned to 1 00 pixels. We have
performed ellipse fitting to both images and rescaled the
NICS image to the flux­calibrated 2MASS data by compar­
ing the light profiles.
The NICS image and WFPC2 images were then realigned
by cross­correlating common features in the nuclear region
because no point sources are present in the near­infrared
image. The main feature that drives the correlation is the
strong nuclear peak; however, we have checked that, even
excluding the central peak, the shift between the images can
be determined with an accuracy of #0>2 (i.e., #2 pixels of
WFPC2) in both directions. With this accuracy the nuclear
continuum peaks present in both optical and infrared
images are consistent with being at the same position.
3. RESULTS
3.1. Kinematics
3.1.1. Line­fitting Procedures
Line­of­sight velocities, FWHMs, and surface bright­
nesses along each slit were obtained by fitting single Gaus­
sians to H# and [N ii] emission lines in each row of the
continuum­subtracted two­dimensional spectra. In the left­
hand panels of Figure 2 we plot the measured velocities.
Beyond 0>5 of the nucleus velocities show considerable
small­scale variations likely due to local gas motions that do
not reflect the mass distribution. Therefore, we averaged
velocities by binning the spectra in steps of 1 00 (10 rows
in NUC and five in POS1 and POS2), while within 0>5 of
the nucleus velocities are measured along each row to take
advantage of all the spatial information. Similarly, the
middle and right­hand panels of Figure 2 display FWHMs
and H# surface brightnesses as measured along the slit.
In Figure 2 the dotted line superimposed on the NUC
rotation curve represents the H# surface brightness shown
in the right middle panel. This helps to distinguish the pres­
ence of two components: an extended one, characterized by
a low surface brightness (<10 #13 ergs cm #2 s #1 arcsec #2 ),
roughly constant along the slit, and a nuclear one, compact
(within #0>4 and 0>4), bright, and cospatial with the posi­
tion of the nuclear continuum peak. This might be inter­
preted as a nuclear disk of the same extent as the nuclear
stellar cluster described in x 3.2. As shown in more detail in
the right­hand panel of Figure 3, the emission­line surface
brightness of the nuclear disk is double peaked while the
nuclear continuum source, roughly coincident with the cen­
ter of rotation, is located in between the two peaks.
Fitting single unconstrained Gaussians is acceptable for
the extended component, but it does not produce good
results for the points in the nuclear region. In particular, a
single Gaussian fit produces velocities of H# and [N ii] that
di#er by as much as #30--40 km s #1 (see Fig. 3). This is not a
worrisome issue if the amplitude of the rotation curve is a
few hundred kilometers per second but makes interpreta­
tion of the data uncertain in the present case in which the
amplitude of the nuclear rotation curve is only #40 km s #1 .
A careful analysis of the line profiles in the nuclear region
shows that they are persistently asymmetric with the pres­
ence of a blue wing (see Fig. 4).
A fit row by row with two Gaussian components, the
main component and the blue one, with the constraint that
they have the same velocities and widths for H# and [N ii],
shows that, within the large uncertainties, the `` blue wing ''
has always the same velocity and width. We have then
deblended the `` blue '' component in the spectrum obtained
by co­adding the central five rows. The velocity and width
of the blue component were then used in the row­by­row fit.
The constrained fit is good and H# and [N ii] now have
the same velocity in the main component. The measured
velocities, FWHMs, and H# surface brightnesses in the
Fig. 2.---Velocity (left), FWHM (middle), and surface brightness (right) measured along the slit at POS1 (top), NUC (middle), and POS2 (bottom). Vertical
bars are 1 # errors, while horizontal bars indicate the size of the aperture over which the quantity was measured. In the left­hand panel, the dashed line is the
velocity gradient measured in the ground­based data. The dotted line in the NUC rotation curve panel is the line surface brightness along the slit, drawn in
order to pinpoint the high surface brightness nuclear disk. The dotted lines in the FWHM panels correspond to instrumental widths. The 0 position along the
slit corresponds to the position of the nuclear continuum peak.
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 871

Fig. 3.---Same as previous figure, but for the points related to the nuclear disk. The points with the error bars are the quantities derived in the fit that takes
into account the presence of the blue wing and in which H#­ and [N ii]--emitting media are constrained to have the same velocity and FWHM.Conversely, the
points without the error bars are derived with unconstrained, single Gaussian fits of H# (open squares) and [N ii] (crosses). For more details see x 3.1.1.
Fig. 4.---Fits of the line profiles at the three slit positions. The dotted lines identify the two components of the fit. Each component is characterized by the
same velocity and FWHM for H# and [N ii]. The ratio between the two [N ii] lines is that fixed by atomic physics. At every slit position the blue component
has the same velocity and width at each row and their values are determined in the fit of the overall nuclear spectrum shown in the top panels. The numbers in
the upper right corners of each panel are the row or the range of rows where the fit was performed. The strong narrow lines in the right­hand panels are residual
cosmic rays that have been excluded from the fit.

nuclear region are shown in Figure 3, where we also plot, as
a comparison, the values obtained from unconstrained sin­
gle Gaussian fits of H# and [N ii]. Given the S/N of the
present data, the nature of the blue wing is, as yet, unclear.
3.1.2. Velocity Curves
The velocity field of the extended component shows a
quasi­linear gradient of #8 km s #1 arcsec #1 , which agrees
with the velocity gradient (Fig. 5) measured from ground­
based observations at the same P.A. as our STIS spectra
(D. J. Axon et al. 2003, in preparation). The velocity gra­
dient measured from the ground­based rotation curve is
shown as a straight line in the figures. The large­scale trend
observed from HSTmatches very well the expected gradient
from the solid­body part of the ground­based rotation curve
although it presents structures at small scales.
Within 0>5 of the nucleus, the velocity field shows a
smooth S­shaped curve with a peak­to­peak amplitude of
only 40 km s #1 (left­hand panel of Fig. 3). Overall there is
no hint of a steep Keplerian rise around a pointlike mass.
The central velocities of the nuclear curves are systemati­
cally o#set from the large­scale velocity field by #10--20 km
s #1 (compare with the dashed line that is the solid­body part
of the rotation curve; left­hand panel of Fig. 3). The o#­
nuclear slits show essentially the same velocity field as the
on­nucleus slit. Thus, this blueshift must be real and not an
instrumental artifact generated by light entering the slit o#­
center (for a detailed discussion on the e#ects of light enter­
ing the slit o#­center see Maciejewski & Binney 2001). For
example, consider the case in which the peak of line emis­
sion is on the blue side (left on Fig. 1) of the NUC slit. Then
the measured velocity will be blueshifted with respect to the
true value. The same would happen for POS1. However, in
POS2 the peak of the line emission will be on the red side of
the slit and the measured velocities will be shifted to the red
with respect to the true value. This is not the case for our
data. The possible origin and implications of this blueshift
are discussed in x 5.
3.2. Morphology
The acquisition image, shown in Figure 1, has a field of
view of 5 00 # 5 00 . In Figure 6 we show the inner 30 00 # 30 00
and 5 00 # 5 00 of the WFPC2 and NICS images. Finally, in
Figure 7 we plot the color maps.
The acquisition image (Fig. 1) shows a compact but re­
solved (FWHM ' 0>2) bright central feature superimposed
to the central region of the bulge. A second bright feature is
present #0>4 southwest of the nucleus. Figure 8 represents
the radial light profile from the STIS acquisition image. As
already clear from the acquisition image, the light profile
indicates the presence of a bright central feature. This fea­
ture is spatially extended as can be seen from a comparison
with the light profile of an unresolved source (NGC 4051)
obtained with the same instrumental setup. This nuclear
feature is likely to be a star cluster; a photometrically dis­
tinct star cluster is often present in the dynamical center of
spiral galaxies of all Hubble types (e.g., Carollo, Stiavelli,
& Mack 1998; Carollo et al. 2002; Bo ˜ ker et al. 2002 and
references therein).
The same central source is visible in all the WFPC2
images, and an analysis of its colors (Fig. 7) indicates that it
is bluer than the surrounding regions but still redder than
the galactic disk and bulge. Its Vega magnitudes in the three
WFPC2 filters are 19.10 (F450W), 18.05 (F606W), and
17.16 (F814W) with formal uncertainties of #0.02 mag.
These values, converted to Johnson magnitudes, become
B ¼ 19:2, R ¼ 17:7, and I ¼ 17:1 for the K0 V and Sb spiral
spectra and B ¼ 19:1, R ¼ 18:0, and I ¼ 17:15 for the A0 V
and Sc spiral spectra. The star cluster emission is very weak
in the K band and not readily visible, but its location is coin­
cident with the location of the K­band peak within the align­
ment uncertainties (#0>2; see Fig. 6).
The color images show that the inner few arcseconds are
redder than the galactic disk and bulge and this could be
due to either the presence of obscuration by dust or a change
in stellar population.
From several galaxy catalogs it can be inferred that the
inclination of the large­scale disk is #20 # . This is supported
by the nearly circular isophotes observed in the NICS and
WFPC2 images at large scales greater than 20 00 --30 00 . How­
ever, in the inner regions, a barlike structure is visible in the
K­band image oriented approximately east­west. It is
extended for #15 00 , i.e., #1.4 kpc with isophotes symmetri­
cally twisted on both ends. The twisting is most likely due to
spiral arms just outside (and apparently connecting to) the
bar. These give the bar the appearance of being larger than
it is. It is well known (see, e.g., Shaw et al. 1995; Knapen et
al. 1995) that strong concentrations of young stars can dom­
inate the K band. This likely happens in the region where
spiral arms emanate from the bar in NGC 4041: blue­star
complexes are seen in color images there, particularly at the
west end (Fig. 6). The bar itself is rather short, with semi­
major axis length at most #3 00 and a P.A. of about 60 # . It is
a weak bar, given the fatness of the observed structure and
that it is lacking the pair of dust lanes characteristic for
strong bars (Athanassoula 1992). Instead, a minispiral
structure is seen in the STIS ACQ image (Fig. 1), extending
at least out to 3 00 from the nucleus, i.e., throughout the entire
IR bar. Such a morphology indicates that we should not
Fig. 5.---Ground­based rotation curve obtained at P:A: ¼ 43 # . The
dashed line is the constant velocity gradient from the solid­body part of the
rotation curve also shown in Figs. 2 and 3.
BLACK HOLE AT CENTER OF NGC 4041 873

Fig. 7.---20 00 # 20 00 view of the color images of the central region of NGC 4041. The center of the images coincides with the location of the nucleus. North is
up and east is to the left. From left to right: B#I, B#I degraded to NICS resolution, B#K with B (from WFPC2/F450W) degraded to NICS resolution. The
dark regions have redder colors.
Fig. 6.---Overlay of the K­band isophotes on the gray scales of the F450W (left) and F814W (right) WFPC2 images. North is up and east is to the left. The
bottom panels show an expanded view of the nuclear region.

expect strong departures from the circular motion in the gas
flow. In fact, the CO velocity field observed by Sakamoto et
al. (1999; see their Fig. 1f ) shows no such deviations down
to the 4 00 beam size of the observations. The P.A. of the kine­
matical line of nodes is '85 # , close to the P.A. of the bar.
3.2.1. Light Profiles
The canonical way to derive light profiles from images is
to fit ellipses to the isophotes and to measure encircled
fluxes. However, inspection of Figures 1 and 6 indicates the
presence of small­scale structures invalidating the ellipse fits
of the WFPC2 and STIS images. The NICS image is much
smoother, but the ellipticity of its inner isophotes is caused
by the presence of a bar and spiral arms and not by geomet­
rical projection e#ects. Since the inclination of the large­
scale disk is only '20 # , one expects almost circular iso­
photes, and it is reasonable to obtain the light profiles from
azimuthal averages instead of canonical ellipse fitting.
In Figure 9 we plot the azimuthally averaged light profiles
extracted from the WFPC2 and NICS images. The light
profiles have all been rescaled to match the external one of
the F814W image. The profiles all match well beyond
r > 10 00 , and the colors of that region can be estimated as
B#I # 1:7, R#I # 0:4, and I#K # 2:1. These colors are
consistent with other spiral galaxies: for instance, de Jong
(1996) has shown that the spirals of this Hubble type on
average have integrated colors B#I ¼ 1:8 # 0:2 (observed
1.7) and B#K ¼ 3:5 # 0:3 (3.8).
The red feature observed in the color maps around the
nucleus results in a flattening of the light profiles, which
increases at decreasing wavelength. From the figure one
can immediately estimate the color excess or reddening for
each pair of bands with respect to the external points.
Roughly EÏB#I÷ # 0:8, EÏR#I÷ # 0:4, and EÏI#K÷ # 0:5.
Assuming the Galactic extinction law by Cardelli, Clayton,
& Mathis (1989), one finds that with A V # 1:2 mag one can
approximately match the observed color excesses [AV # 1:2
implies EÏB#I÷ # 1:03, EÏR#I÷ # 0:32, and EÏI#K÷ #
0:44]. Thus, the red feature observed in the color maps is
likely to be due to extinction by dust with an average value
of AV # 1:2 mag.
4. MODEL FITTING
4.1. The Stellar Luminosity Density
In order to estimate the stellar contribution to the gravita­
tional potential in the nuclear region, we have to derive the
stellar luminosity density from the observed surface bright­
ness distribution. This inversion is not unique if the gravita­
tional potential does not have a spherical symmetry. Here we
assume that the gravitational potential is an oblate spheroid
(e.g., van der Marel & van den Bosch 1998). We further
assume that the principal plane of the potential has the same
inclination as the large­scale disk (i ' 20; LEDA; de
Vaucouleurs et al. 1991). The knowledge of i, combined with
the observed axial ratio of the isophotes, should directly pro­
vide the axial ratio of an axisymmetric light distribution, but
the presence of a barlike structure invalidates this approach
since the flattening of the isophotes cannot be ascribed to
simple projection e#ects. Therefore, we consider two extreme
cases: a spherical and a disk light distribution. The model
light profiles are computed taking into account the finite
spatial resolution and pixel size of the detector. A detailed
description of the relevant formulae for the inversion and fit­
ting procedure is presented in Appendix A.
The light profiles derived from HST images have
better spatial resolution than the one derived from the
Fig. 9.---Azimuthally averaged light profiles for the WFPC2 and NICS
images. All images have been rescaled to match the external light profile of
the F814W image. The scaling factors applied are indicated in the figure
and directly show the color in the outer regions where the profiles match.
`` NICS res.'' indicates that the F814W image has been degraded to the
NICS resolution.
Fig. 8.---Radial light profile from the STIS acquisition image. The open
diamonds represent the radial light profile from a galaxy, NGC 4051,
having a strong central unresolved source; the profile was rescaled to match
the NGC 4041 one in the central pixel.
BLACK HOLE AT CENTER OF NGC 4041 875

ground­based NICS image. However, as shown above, the
optical HST images su#er from AV # 1:2 mag extinction in
the nuclear region. It is also well known that K­band light is
a better tracer of mass. Therefore, it is useful to use both the
HST and ground­based light profiles in order to infer the
mass profiles. The star cluster is probably washed out by
the lower spatial resolution in the NICS image; thus, in
order to account for its presence, one can determine its
geometrical parameters from the fit of the HST light profiles
and use them in fitting the NICS image.
In principle, the WFPC2 images could be reddening cor­
rected by using the color maps in Figure 7 to derive the
E(B#V ) map. In order to do this, one might assume that
the colors in the nuclear region are constant and equal to
those beyond 10 00 (see Fig. 9). For a more detailed descrip­
tion of this procedure see, for instance, Marconi et al. 2000.
However, the cluster would be forced to have the same
intrinsic color as the galactic disk, which is probably not
realistic. Moreover, its shape would be a#ected by point­
spread function (PSF) di#erences among the images.
Finally, it is not clear if the color di#erences are due entirely
to reddening. We therefore decided not to apply any redden­
ing correction, and the consequences of this will be dis­
cussed in x 5.
We first consider the light profile derived from WFPC2/
F814W because it is the least sensitive to reddening among
the HST images. We then fit a model light distribution com­
posed of the central star cluster and the more extended com­
ponent. We assume that the star cluster luminosity density
is spherically symmetric with a density law:
# # Ïr÷ ¼ # # b 1 ×
r
r # b
# # 2
" # ## #
: Ï1÷
The extended component has a functional form of the type
#Ïm÷ ¼ # b
m
m b
# # ##
1 ×
m
m b
# # 2
" # ##
; Ï2÷
where m 2
¼ x 2
× y 2
× z 2 =q 2 corresponds to the radius in the
spherical case and q is the axial ratio of the mass distribu­
tion. This functional form of the extended component is a
reasonable description of the central part of the galaxy
where the bulge dominates. It is su#cient for our purposes
because we are interested in modeling only the inner r < 5 00 .
In equation (2), q ¼ 1 for a spherical light distribution,
and q ¼ 0:1 for a disklike one. These are the two extreme
cases. The best­fit parameters for the HST light profile are
shown in the upper part of Table 1. The left­hand panel of
Figure 10 shows the fit in the spherical case, which is visually
indistinguishable from the disk case.
Using the geometrical parameters of the cluster deter­
mined from the fit to the HST data (r # b and #*), we fit the
NICS profile, both the spherical and disk case. The results
are shown in Table 1. The right­hand panel of Figure 10
shows the fit in the spherical case, which, as before, is visu­
ally indistinguishable from the disklike case.
The values of the reduced # 2 reported in Table 1 are much
larger than the expected value of 1, and the main reason can
be found in the systematic deviation of the residual around
3 00 (Fig. 10) where the data are systematically lower than the
model. This systematic deviation, whose maximum value is
#0.03 dex, i.e., 7%, is present in both the optical and near­
IR light profile and is caused by the presence of the barlike
feature.
The fits to the light profiles allow an estimate of the I and
Kmagnitudes of the star cluster. With the assumed distance,
the cluster luminosity is 10 7.3 L #wI , where wI refers to the
F814W band (in this band the Sun has M#wI ¼ 4:16). In the
K band, the cluster has 10 7.2 L #K for the spherical case and
10 7.1 L #K for the disk case (M#K ¼ 3:286). These results sug­
gest that the cluster could be blue, with estimated color
I#K ' 0:9 (or 0.5 in the disk case). The estimate of the K
magnitude of the cluster is, however, uncertain as a result of
the fact that the cluster is unresolved in the ground­based
observations. In order to set a limit to the amount of the
cluster contribution to the K­band light profile, in Figure 10
we also plot the results of the fit when the cluster normaliza­
tion has been fixed to a value 0.6 dex larger (i.e., the cluster
is 1.5 mag brighter) than the best­fit value (dashed line). The
fit is worse, as clearly indicated by the residuals, and we con­
sider this a firm upper limit for the cluster contribution to
the K­band light profile. This implies that the observed color
of the cluster is I#K < 2:4. If we assume that the cluster is
subjected to A V ¼ 1:2 mag of extinction (x 3.2.1), the upper
limit changes to I#K < 2:0.
Figure 11 shows the circular velocities in the principal
plane of the potential expected from the density distribu­
tions determined from the fit of the NICS and HST data. In
order to match the di#erent rotation curves, the circular
velocities are plotted for # K ¼ 0:5 or # I ¼ 2:3, the values
derived in x 4.2 from fitting the rotation curves. It is clear
TABLE 1
Best­Fit Parameters of the Stellar Light Density
q
# # b
a
(M # pc #3 )
r # b
(arcsec) #*
# b a
(M # pc #3 )
m b
(arcsec) # # # 2
red
I Band (WFPC2/F814W)
1.0.............. 1023 0.27 4.2 0.73 12.6 0.9 1.6 10.8
0.1.............. 951 0.31 5.0 5.8 12.6 0.9 1.6 10.5
KBand (NICS)
1.0.............. 910 0.27 b 4.2 b 0.44 28.0 1.4 3.0 9.8
0.1.............. 609 0.31 b 5.0 b 1.60 45.9 1.4 6.6 7.7
a Assuming # ¼ 1.
b Kept fixed at the value from the previous fit.
876 MARCONI ET AL. Vol. 586

that the rotation curves are very similar beyond r > 0>5.
However, they di#er at smaller radii because of the di#erent
contribution of the star cluster to the total mass budget.
This di#erence is due to the fact that the mass­to­light ratio
has been assumed constant in each band.
4.2. Kinematics
In order to model the gas kinematical data (velocities and
widths measured along the slit), we select the simplest possi­
ble approach and assume that the ionized gas is circularly
rotating in a thin disk located in the principal plane of the
galaxy potential. The latter assumption is not needed if the
galaxy potential has a spherical symmetry. We assume that
the disk is not pressure supported, and we neglect all hydro­
dynamical e#ects. Thus, the disk motion is completely deter­
mined by the gravitational potential, which is made of two
components: one is stellar and is completely determined by
the mass distribution, derived in the previous section, and
by its mass­to­light ratio #; the other one comes from a dark
mass concentration (the BH), spatially unresolved at
HST+STIS resolution and defined by its total mass M BH .
This is the same approach that has been followed in all pre­
vious gas kinematical studies of circumnuclear gas disks
(e.g., Macchetto et al. 1997; van der Marel & van den Bosch
1998; Barth et al. 2001; Marconi et al. 2001).
In principle, in order to constrain the BH mass, one
should compare the emission­line profile predicted by the
model with the observed one. Even in the simple potential
described above, the line profiles can be very complicated,
with multiple peaks (Maciejewski & Binney 2001). Never­
theless, as shown in Figure 4, the line profiles are well repre­
sented by single Gaussians (after subtracting the `` blue ''
component); therefore, we compare the moments of the line
profiles: the average velocity hvi and velocity dispersion #
defined as # 2
¼ hv 2
i # hvi 2 . In order to be compared with
the observations, the model hvi and # are computed taking
into account the spatial resolution of HST/STIS and the
size of the apertures over which the observed quantities are
measured. The formulae used to compute velocities, widths,
and line fluxes and the details of their derivation, num­
erical computation, and model fitting are described in
Fig. 10.---Left: Fit of the light profile obtained from the WFPC2/F814W data. Right: Fit of the NICS/K light profile with the central star cluster and the
extended component (the geometrical parameters of the star cluster were determined in the previous fit). In both cases the lower panel shows the fit residuals.
The dashed line in the right­hand panel shows the K­band fit obtained by fixing the cluster normalization to a value that is 0.6 dex larger than the best­fit value.
The open squares show the corresponding residuals.
Fig. 11.---Circular velocities in the principal plane of the potential for
the density distributions determined from the fit of the NICS and HST
data. The circular velocities are derived assuming a spherically symmetric
mass distribution and are plotted for #K ¼ 0:5 and # I ¼ 2:3. The dotted
lines represent the extreme case in which the cluster normalization in the K
band has been fixed to +0.6 dex of the best­fit value (see text).
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 877

Appendix B. To derive the observed hvi and #, our approach
is to parameterize the observed spectra by fitting them with
as many Gaussian components as required by the data.
Then, components that are not obviously coming from cir­
cularly rotating gas can be discarded (e.g., Winge et al.
1999) and the average velocities and line widths can be com­
puted from the remaining ones (this is trivial if, as in the
present case, only a single Gaussian component is left).
With this approach, the observed and model parameters are
compared. The advantages of our approach are twofold.
First, much computational time is saved in computing just
hvi and # instead of the whole line profile. Second, modeling
the observed emission­line profiles in detail requires very
high S/N data, and it is seldom possible to obtain such
high­quality data with HST except for a few galaxies (e.g.,
M87; Macchetto et al. 1997). Even with the high S/N avail­
able from 8 m class telescopes the matching of the line pro­
files remains a significant problem (e.g., Cen A; Marconi et
al. 2001).
Thus, the model hvi and # depend on #, the intrinsic sur­
face brightness distribution of the emission lines, and on the
following parameters:
1. s 0 , the position of the kinematical center along the slit
with respect to the position of the continuum peak.
2. b, the distance between the reference slit center of the
NUC slit and the kinematical center.
3. i, the inclination of the rotating disk.
4. h, the angle between the disk line of nodes and the slits.
5. V sys , the systemic velocity of the disk.
6. #, the mass­to­light ratio of the stellar population.
7. M BH , the BHmass.
Not all of these parameters can be independently deter­
mined by fitting the observations.
It can be inferred from the equations in Appendix B that
M BH , #, and i are directly coupled. This fact has already
been discussed in detail by Macchetto et al. (1997), but here
we briefly indicate how this coupling can be used to pose a
lower limit on the disk inclination i. Considering the case in
which no BH is present, the amplitude of the rotation curves
depends linearly on # 0:5 sin i, and its value can be deter­
mined by fitting the observations. Therefore, # will increase
with decreasing i. An upper limit on # can be certainly set
by the largest reasonable value for a very old stellar popula­
tion (e.g., Bruzual & Charlot 1993; Maraston 1998), and
this will immediately fix the lower limit for i. If a BH is
present, # 0:5 sin i can still be determined from the kinemati­
cal data obtained beyond the BH sphere of influence. Of
course, the underlying assumption is that both # and i are
constant.
4.2.1. The Intrinsic Surface Brightness Distribution
of Emission Lines
The intrinsic surface brightness distribution of emission
lines # is an important ingredient in the model computa­
tions because it is the weight in the averaging of the
observed quantities. The ideal situation would be to have an
emission­line image with spatial resolution higher than the
spectra. Not having that, Barth et al. (2001) have used a
deconvolved HST image to model their HST data. This
helped them to match the microstructure of the velocity field
at large radii, but it is important to remember that the infor­
mation about the BH mass comes predominantly from
within the central unresolved source of unknown luminosity
distribution. Our approach is to attempt to match the emis­
sion­line fluxes as observed in the spectra (or in emission­
line images) together with a variety of synthetic realizations
of this unknown central emissivity distribution. We use this
to estimate its impact on the derived value of M BH . Thus, we
assume simple analytical forms of the line flux distribution
that well match the observed one along the slit, after folding
with the telescope and instrument responses as described
above.
It is clear from the right­hand panel of Figure 3 that the
flux distribution observed at the three slit positions is very
complex and cannot be described by a single radially sym­
metric component in the disk plane. Therefore, we choose
the approach of modeling directly the flux distribution on
the plane of the sky by exploring various simple analytical
forms of the line flux distribution. We take the intrinsic light
distribution at point x, y in the plane of the sky in the NUC
slit reference frame (x is the position across the NUC slit,
while y is the position along it) to be defined as
IÏx; y÷ ¼ X i
I 0i f i
r i
r 0i
# # ; Ï3÷
where f i is a circularly symmetric function of characteristic
radius r 0i and weighting factor I 0i . Each component function
f i is centered at (x 0i , y 0i ), and r i is the radial distance from
this location r i ¼ ½Ïx # x 0i ÷ 2
× Ïy # y 0i ÷
2# 0:5
. This is an
approach similar to that of Barth et al. (2001), except that
we use a synthetic realization of the line flux map. Although
each bright patch is assumed to be circularly symmetric, its
material is in circular rotation about the galactic nucleus.
Hence, patches are constantly shearing rather than moving
as coherent units.
The free parameters characterizing the flux distribution
are therefore (x 0i , y 0i ), I 0i , and r 0i and are chosen in order to
match the observed flux distribution along the slit after fold­
ing the model with telescope and instrument. For example,
one of the functional forms used is the combination of four
exponentials and a constant baseline, i.e., f Ïr i =r 0i ÷ ¼
expÏr i =r 0i ÷. Figure 12 shows the contours of this line flux
distribution, which matches the observed profile along the
slit after convolving with the instrumental response (see
Fig. 13).
4.2.2. The Standard Approach
The standard approach followed so far in gas kinematical
analysis is to assume that (1) gas disks around BHs are not
warped; i.e., they have the same line of nodes and inclina­
tions as the more extended components; and (2) the stellar
population has a constant mass­to­light ratio with radius
(e.g., van der Marel & van den Bosch 1998; Barth et al.
2001). In the present case, the blueshift of the inner disk
(x 3.1.2) indicates that the standard approach must be gen­
eralized to allow for the kinematical decoupling between
inner and large­scale disks.
We use the emission­line flux distribution derived in the
previous section and fix the inclination to i ¼ 20 # , i.e., the
inclination of the galactic disk. The free parameters of the fit
are then M BH , #, s 0 , b, h, V sys , and DV sys , the velocity shift
of the extended component with respect to the nuclear one.
We perform the fit using the mass density profiles derived
for the I and K band, in both the spherical and disk cases.
The results of the fit are shown in Table 2. Statistical errors
878 MARCONI ET AL. Vol. 586

at the 2 # level on log M BH and log # are #0.2 and #0.1,
respectively; thus, the derived BH mass is M BH ¼ Ï1 ×0:6
#0:7 ÷
#10 7 M # . Figure 14 shows the best­fit model (solid line)
obtained with the mass density distribution derived from
the K­band light profile with the assumption of spherical
symmetry. The dotted line is the best­fit model without a
BH. The left­hand panel of Figure 16 shows the fit of the
NUC data from the model with the mass distribution
derived from the I band. Figures 14 and 16 indicate that a
model with no BH cannot account for the observed nuclear
rotation curve, producing a velocity gradient that is shal­
lower than observed. Note that the P.A. of the line of nodes,
a free parameter of the fit, is the same as the one inferred
from the large­scale CO velocity map by Sakamoto et al.
(1999). This supports the assumption that the nuclear disk is
the continuation at small scales of the galactic disk. How­
ever, the fit confirms that the nuclear disk is blueshifted by
#10 km s #1 with respect to the extended component, and
this argues against the assumption in the previous sentence.
The fit is greatly improved if a velocity shift of 8 km s #1 is
allowed for the POS2 data. This velocity shift is the conse­
quence of an error on the absolute wavelength calibration
of the POS2 data and is well within the expected STIS per­
formances (Leitherer et al. 2001).
The value of # in the K band derived from the fit varies
between 0.2 (disk light distribution) and 0.5 (spherical light
distribution). This range of values is in good agreement with
the typical K­band mass­to­light ratios of spiral bulges
(Moriondo, Giovanardi, & Hunt 1998). Similarly, the value
of # in the I band ranges between 2.2 and 3.6, in agreement
Fig. 12.---Contours of the modeled emission­line flux distribution of the
nuclear gas disk with overplotted slit positions (vertical dashed lines), the
position of the kinematic center ( filled square), and the line of nodes derived
from fitting the kinematic data (solid straight line). The dotted line is the line
of nodes of the extended material. The ellipse indicates the nuclear gas disk
outer limit for an inclination of 20 # with respect to the line of sight.
Fig. 13.---Assumed flux distribution compared with the data after fold­
ing with the telescope and instrument responses. The model values are
connected by straight lines in order to guide the eye.
TABLE 2
Model Fit Results in the Standard Approach
Band a q b
log M BH
(M # ) # c
s 0
(arcsec)
b
(arcsec)
h
(deg)
V sys
(km s #1 )
DV sys
(km s #1 ) # 2
red
K ..................... 1.0 7.04 0.52 #0.04 0.02 #45 1212 11 2.32
K ..................... 1.0 0.0 d 1.00 #0.09 0.06 #58 1211 11 2.88
K ..................... 0.1 7.09 0.22 #0.07 0.08 #32 1213 11 2.13
I ...................... 1.0 6.86 2.28 #0.01 0.00 #42 1210 12 1.86
I ...................... 1.0 0.0 d 2.65 #0.03 0.03 #45 1210 12 2.09
I ...................... 0.1 7.14 1.65 #0.03 0.01 #47 1212 10 2.22
a Band from which the mass density profile was derived.
b Assumed axial ratio of the mass distribution.
c Mass­to­light ratio in used band.
d The parameter was fixed to this value.
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 879

with measurements within the inner 2 kpc of spiral galaxies
(Giovanelli &Haynes 2002).
4.2.3. Alternative Model: The Decoupled Nuclear Disk
To date, everyone who has determined a central BH mass
from gas kinematics assumed that the disks are unwarped,
i.e., coplanar at all radii, and that the stellar mass­to­light
ratio is constant with radius.
The high surface brightness and velocity o#set of the
nuclear disk with respect to the extended component pre­
sented in x 3 indicate that, at least for NGC 4041, this
assumption could be untrue. The nuclear and large­scale
disks might be characterized by di#erent geometrical prop­
erties like P.A. of the line of nodes, inclination, and systemic
velocity. Recently, Cappellari et al. (2002) have shown a
discrepancy between the BH mass in IC 1459 determined
from gas kinematical (Verdoes Klein et al. 2000) and stellar
dynamical models. The authors propose that a possible sol­
ution to this discrepancy could be an error on the assumed
gas disk inclination.
Therefore, in the alternative approach presented here, we
first fit the extended component data in order to determine
the mass­to­light ratio # of the stellar component to be used
in the fit of the nuclear data. For the extended component
we assume a constant line flux distribution and determine #
in the two extreme cases of the spherical and disklike light
distribution (i.e., with an axial ratio of q ¼ 0:1).
From the CO velocity map by Sakamoto et al. (1999) the
P.A. of the kinematic line of nodes is 86 # , i.e., almost exactly
east­west. Since the P.A. of the slit in the STIS observations
is 43 # , this means that the angle between the slit and the
disk line of nodes is # ¼ #43 # (for the definition of h see
Appendix B). We also assume that the inclination of the disk
is i ¼ 20 # , as specified in the previous section. The CO veloc­
ity map also clearly indicates that the large­scale velocity
field is circularly symmetric.
The fitted circular rotation curves of the extended gas are
shown in Figure 15. The derived # for the spherical or disk
geometry are shown in Table 3 for the mass density profiles
derived from both the K­ and I­band light. In our derivation
we allowed for a constant velocity shift between the data
points in the o#­nuclear slit positions and the nuclear one.
This is caused by the uncertainty in the absolute wavelength
calibration and is very small, #1 and 8 km s #1 for POS1 and
POS2, respectively.
Similarly to the results of the standard model, mass­to­
light ratios in the K and I band are consistent with typical
values for the central regions of spiral galaxies.
We now fit the nuclear data (r # 0>4) allowing for values
of # or i di#erent from those of the extended component.
Results of the fit are shown in Figure 15 and Table 3. The
right­hand panel of Figure 16 shows the fit of the NUC data
from the model with the mass distribution derived from the
I band.
In all these cases the BHmass has been set to zero and the
goodness of the fit suggests that no BH is needed to fit
the data with this approach. We estimate the upper limit on
the BH mass by a one­parameter variation; i.e., we assign a
Fig. 14.---Best­fit standard model of the observed rotation curves (solid line) compared with the data. The dotted line is the best­fit model without a BH. The
model values are connected by straight lines in order to guide the eye. Note that points from external and nuclear regions are not connected because they are
kinematically decoupled.
880 MARCONI ET AL. Vol. 586

fixed value of M BH and perform the fit. When we have
D# 2
# 1, 4, 9 we have reached the 1, 2, 3 # upper limit. As an
example, in Figure 17 we show the # 2 plot corresponding to
the K­ and I­band spherical cases. There the 3 # upper limit
is 10 6.8 M # .
From Table 3 it can be seen that the position of the
nucleus along and across the slit is independent of the
assumed mass distribution. The inclination varies but,
like #, is just a scaling factor. Indeed, varying # produ­
ces the same results, as shown in the same table. The
Fig. 15.---Best­fit alternative model of the observed rotation curves compared with the data. The solid line is the best­fit model obtained by fixing i and
varying #, while the dotted line represents the case in which i has been varied and # kept fixed. The model values are connected by straight lines in order to
guide the eye. Note that points from external and nuclear regions are not connected because they are kinematically decoupled. The right­hand panel is a
zoom on the nuclear disk region. The plotted model uses the mass density distribution derived from the K­band light profile with the assumption of spherical
symmetry.
TABLE 3
Model Fit Results in the Alternative Approach
Band a q b
log M BH
(M # ) # c
s 0
(arcsec)
b
(arcsec)
h
(deg)
i
(deg)
V sys
(km s #1 ) # 2
red
Extended Component
K ..................... 1.0 0.0 d 0.48 0.0 d 0.0 d
#43 d 20 d 1224 2.7
K ..................... 0.1 0.0 d 0.31 0.0 d 0.0 d
#43 d 20 d 1224 2.7
I ...................... 1.0 0.0 d 2.29 0.0 d 0.0 d
#43 d 20 d 1224 2.5
I ...................... 0.1 0.0 d 1.38 0.0 d 0.0 d
#43 d 20 d 1224 2.5
Nuclear Disk
K ..................... 1.0 <6.8 1.29 0.0 #0.05 #40 20 d 1204 0.66
K ..................... 1.0 <6.6 0.48 d
#0.01 #0.04 #39 35 1205 0.67
K ..................... 0.1 <6.8 0.88 #0.01 #0.05 #40 20 d 1205 0.69
I ...................... 1.0 <6.8 4.05 #0.01 #0.04 #41 20 d 1206 0.60
I ...................... 0.1 <6.8 3.77 0.00 #0.04 #41 20 d 1205 0.58
a Band from which the mass density profile was derived.
b Assumed axial ratio of the mass distribution.
c Mass­to­light ratio in used band.
d The parameter was fixed to this value.
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 881

angle between the slit and the line of nodes is very well
determined in the range #39 # to #41 # . This means that
the P.A. of the kinematical line of nodes is 43 # (the P.A.
of the slit) minus #40 # , i.e., #83 # , consistent with the
P.A. of the kinematical line of nodes on large scales. The
systemic velocity is 1204--1206 km s #1 and is definitely
blueshifted with respect to the extended rotation where it
is 1224 km s #1 .
As pointed out by Barth et al. (2001), the FWHM of the
lines might be a worrisome issue in the sense that FWHM
much larger than instrumental might indicate motions that
could deviate from circular (apart from the broadening of
the line due to unresolved rotation around a point mass). In
Figure 18 we plot the observed FWHMs and compare them
with the expected one in the framework of the alternative
model assuming an intrinsic velocity dispersion of only 35
km s #1 . In addition, the dashed line represents the case with
a 10 7 M # BH showing that, in this case, the FWHM cannot
pose good constraints on the BHmass.
5. DISCUSSION
5.1. How the BHMass Determination
is A#ected by Assumptions
In x 4.2 we have shown that the inclination of the nuclear
disk, i, and the mass­to­light ratio of the stellar population,
#, play a critical role in detecting the presence of a nuclear
BH. The standard model, in which i and # are the same for
Fig. 16.---Best­fit models of the observed rotation curves computed using the mass densities derived from the I­band light profiles. The model values are
connected by straight lines in order to guide the eye. The left­hand panel refers to the standard model (the dotted line is the model without a BH), while the
right­hand panel refers to the alternative model.
Fig. 17.---Statistical upper limits on the BH mass in the alternative
approach for the cases in which the mass density has been derived by the
I­ and K­band light profiles.
882 MARCONI ET AL. Vol. 586

the nuclear disk and the more extended galactic disk,
requires the presence of a BH with M BH ¼ Ï1 ×0:6
#0:7 ÷ # 10 7 M #
(x 4.2.2). Conversely, allowing either # or i to vary, the kine­
matical data can be fitted without a BH (the so­called alter­
native model; x 4.2.3) and M BH < 6 # 10 6 M # . These results
show a caveat in gas kinematical searches, where it is always
assumed that i and # are constant at all radii.
To distinguish between the standard and alternative mod­
els, it is necessary to find some extra physical constraints. In
particular, the value of mass­to­light ratio # determined
from the observations should be consistent with realistic
stellar populations (see x 4.2). Little can be said about i,
which, in contrast with the other geometrical parameter of
the disk, h, is poorly constrained by the rotation curves.
Usually one can infer an estimate of i from the morphology
of the disk observed in emission­line images, especially in E
and S0 galaxies (e.g., van der Marel & van den Bosch 1998;
Barth et al. 2001). The root of the problem lies in the cou­
pling between M BH , #, and i that arises because they are all
essentially scaling factors on the amplitude of the rotation
curves.
In the standard approach (Table 2) M BH does not depend
on the assumed mass density profile (either from the K or I
band). The angle h, which is a free parameter of the fit, is
quite well constrained around #40 # ; i.e., the PA of the kine­
matical line of nodes is the same as the one at large scales
derived from the two­dimensional CO map by Sakamoto et
al. (1999). The constancy of h, which is not an intrinsic
parameter of the system but only depends on the galaxy ori­
entation in the plane of the sky, suggests that i (Di between
the nuclear and galactic disk is an intrinsic parameter of the
system) should also be constant. No similar argument can
be placed on #, which can then be allowed to vary, and it is
then possible to fit the kinematical data without the pres­
ence of a BH. This requires that # in the nuclear region is a
factor of #2 larger than that in the extended disk.
To verify if this # is consistent with the observed colors,
we plot in Figure 19 I#K, # K , and # I as a function of time
from the burst for a single stellar population experiencing
an instantaneous burst of star formation. We have used the
models by Bruzual & Charlot (1993), updated in 1996; we
considered solar metallicities, Salpeter (1955) and Scalo
(1986) initial mass functions (IMFs; solid and dashed lines
in the figure, respectively); and we used the theoretical stel­
lar libraries. The solid line in the top panel marks the upper
limit we placed on the I#K color of the stellar cluster (x 4.1);
the dotted line marks the same upper limit but dereddened
for A V ¼ 1:2 mag (x 3.2.1). Finally, the dotted lines in the
middle and bottom panels shows the # required to fit the
data in the alternative model for both the spherical and
disklike distribution of stars. The #­values implied by the
analysis are thus consistent with a very old nuclear star clus­
ter, with age #10 Gyr. This is apparently in contrast with
the high H# equivalent width measured in the nuclear spec­
trum, '70 A š , which seems to require the presence of very
young stars. A possibility is that the H# emission is not due
to young stars but to a low­luminosity AGN, and indeed
the total H# luminosity in the nuclear region is
LÏH#÷ ¼ 4:9 # 10 4 L # , corresponding to the low end of the
L(H#) distribution for AGNs (Ho et al. 1997). To prove
Fig. 18.---Observed and model FWHMs at the NUC position. The solid
line is the case without a BH (alternative model), and the dashed line is with
a 10 7 M # BH (standard model). The assumed intrinsic velocity dispersion is
35 km s #1 .
Fig. 19.---I#K (top), log M=LK
Ï ÷ (middle), and log M=L I
Ï ÷ (bottom) vs.
time for a single stellar population (instantaneous burst) with solar metal­
licity (models by Bruzual & Charlot 1993, updated in 1996). The solid and
dashed lines are for a Salpeter (1955) and a Scalo (1986) IMF, respectively.
In the top panel, the solid horizontal line represents the observed upper
limit on the I#K color. The dotted line is the same upper limit after correc­
tion for AV ¼ 1:2 mag. In the middle and bottom panels, the horizontal
dotted lines indicate the range of mass­to­light ratios allowed by the
observations.
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 883

that the emission lines are excited by an AGN, one could
use the so­called diagnostic line ratios [N ii]/H# and [O iii]/
H# (e.g., Osterbrock 1989). Ground­based spectra with 1 00
resolution (D. J. Axon et al. 2003, in preparation) indicate
that this object is on the boundary between H ii regions and
LINERs with ½N
ii# =H# # 1 and ½O
iii# =H# < 0:5. Unfor­
tunately, the spectra coverage of our HST data is limited to
the 6500--7000 A š region, but the [N ii]/H# ratio at 0>1
resolution remains similar to that observed from the
ground, suggesting that the same could be true for the
[O iii]/H# ratio.
If, on the other hand, the star cluster is young, it can fully
account for the H# emission. Assuming case B recombina­
tion, the L(H#) luminosity implies an ionizing photon rate
of QÏH÷ # 2 # 10 50 s #1 (e.g., Osterbrock 1989), which cor­
responds to #2--20 O stars (Panagia 1973). When con­
fronted with the stellar population synthesis model
(Leitherer et al. 1999), the observed H# equivalent width
indicates an age of #10 7 yr. In this case, as shown in Figure
19, the I­band mass­to­light ratio is more than a factor of 10
lower and the fit of the kinematical data requires a BH even
in the alternative model. However, when using the stellar
mass density profile derived from the K band, the data can
still be explained by stars alone. The reason is that the con­
tribution of the star cluster to the total light density profile
is small, as can be seen from the rotation curves in Figure
11. Note that in order to explain the high [N ii]/H# ratio
with photoionization from stars one needs to allow for
higher than solar N abundances in the nucleus as has been
suggested many times in the past since the work by Burbidge
&Burbidge (1962).
In order to proceed further and firmly establish the age of
the star cluster, one needs HST/NICMOS infrared images
of the nuclear region. The optical--to--near­IR colors of the
star cluster constrain its age and mass­to­light ratio. More­
over, with HST/STIS blue spectra one can measure the
[O iii]/H# ratio, in order to estimate the AGN contribution,
and verify the existence of deep H Balmer lines in absorp­
tion expected if the cluster is young.
The other assumptions behind the modeling are the shape
of the stellar density profile (spherical vs. disklike, I vs. K
band, e#ects of reddening), the intrinsic flux distribution of
the emission lines, and, finally, the assumption of circular
motions.
The stellar density profile depends on q, the intrinsic axis
ratio of the extended mass distribution, and on the relative
importance between the nuclear star cluster and the
extended stellar component. As shown in the previous sec­
tion, the four cases examined there cover the most extreme
cases: spherical or disklike distribution of the extended com­
ponent (q ¼ 1 and 0.1, respectively), stellar cluster that gives
a small contribution to the total stellar density distribution
(K band) or that dominates over the nuclear disk size, i.e.,
r < 0>4 (I band). From Tables 2 and 3 we can conclude that
the value of the BH mass or its upper limit does not depend
on the assumed stellar density profile. From the same tables
it can be seen that the position of the nucleus along and
across the slit is also independent of the assumed mass dis­
tribution. In x 4.1 we decided not to apply a reddening cor­
rection to the light profiles because of the uncertainties
involved. This reddening correction would have a#ected the
values of the mass­to­light ratios derived from the rotation
curves (see x 4.2) whose accurate measurement is not the aim
of this paper. Roughly, assuming an average AV ¼ 1:2 mag
(x 3.2.1), the I light would increase by #70% and # I would
decrease by the same amount (#0.2 in log). The K light
would similarly increase only by #10% with a similar
decrease in # K (#0.04 in log). These corrections would have
a marginal e#ect for the conclusions on the cluster age
derived in Figure 19. Finally, the reddening correction
e#ects on the shape of the stellar light/mass density profiles
would be negligible for the estimate of the BH mass upper
limit. Indeed, the same value is obtained regardless of the
use of mass density profiles derived from the I­ or K­band
light, which are very di#erently a#ected by reddening.
Another important issue is the influence of the intrinsic
line flux distribution on M BH . We have tried several other
di#erent realizations of the observed flux distribution by
varying both the number of components and the functional
form of single components (exponential, Gaussian, power
law). We have found that accurately describing the observed
flux distribution is important. Models that did not do this
produced fits significantly worse than those described in the
previous section. When the fits of fluxes and velocities are
acceptable, M BH or its upper limit does not depend on the
assumed flux distribution.
The assumption of circular motions is probably one of
the most critical ones. Circular motions in gaseous disks are
expected as a result of the dissipative nature of the gas; how­
ever, anisotropies in the stellar potential rapidly lead to non­
circular streaming (e.g., Athanassoula 1992). The e#ect of
noncircular motions would be that of `` distorting '' the rota­
tion curves and increasing the intrinsic line width of the gas.
Therefore, obtaining very good fits (# 2
red < 1) of the rotation
curves in the alternative model strengthens our confidence
on the assumption of circular motions. This assumption is
also supported by the small intrinsic line width (# ' 35 km
s #1 ) required by the observations. Indeed, #=v circ can be
used to quantify the e#ects of turbulence or noncircular
motions on the BH mass determination (e.g., Verdoes Klein
et al. 2000; Barth et al. 2001). For the best­fitting models
#=v circ is always less than 0.4 (provided that there is a BH
with a mass close to the upper limit in the case of the alterna­
tive models). This is similar to what Barth et al. (2001) found
in the case of their best­fit model of NGC 3245. There the
e#ect of the asymmetric drift correction, a possible way to
include noncircular motions in the analysis (e.g., Binney &
Tremaine 1987), is to increase the estimate of the BH mass
by just 10%. Apart from these qualitative arguments, the
presence of significant noncircular motions cannot be
excluded, and it is not possible at the moment to quantify
their e#ects on the method we have adopted and described
here.
Noncircular motions are certainly present in the nuclear
region of NGC 4041 as indicated by the presence of the blue
wing, but these have been singled out in the deblending pro­
cedure as in Winge et al. (1999).
5.2. The Blueshift of the Nuclear Disk
A proposed explanation for the blueshift observed in the
nuclear disk is that the star cluster is oscillating across the
galactic plane. In this picture, the cluster is bound to the gal­
axy, is old, and has a large # ensuring that it is massive
enough and not subjected to tidal disruption.
Any velocity component perpendicular to the disk will
give a large contribution to the observed velocity because
the galaxy is almost face­on. If the cluster is oscillating
884 MARCONI ET AL. Vol. 586

across the galactic plane, the observed blueshift may be
translated into a velocity modulus of similar magnitude.
Thus, the cluster velocity (#10--20 km s #1 ) is smaller than
the rotational velocity (see Fig. 11) and the cluster is bound
to the galaxy. If the cluster is very massive, as in the alterna­
tive model, then the gravitational potential over the star
cluster size is dominated by the star cluster itself and is not
subjected to tidal disruption. Also in this picture the gaseous
disk is completely dominated by the gravitational potential
of the star cluster. The total cluster mass within r < 0>3 as
derived from the I­band mass density profile (the one in
which the cluster dominates) is #6 # 10 7 M # , and this gives
an average density of '200 M # pc #3 , not an unusually high
value, since it can be found in Galactic globular clusters
(Pryor &Meylan 1993).
In the picture in which the cluster is young (from the fit of
the K band we have I#K ¼ 0:9, i.e., #10 7 yr from Fig. 19),
the stability against tidal disruption could pose a problem,
since the mass in the nuclear region is completely dominated
by the bulge stars. However, this problem could be solved if
the lifetime of the cluster against tidal disruption is less than
or equal to the cluster age.
5.3. M BH ­Galaxy Correlations
It is now clear that a large fraction of hot spheroids con­
tain a BH and, moreover, it seems that the hole mass is pro­
portional to the mass (or luminosity) of the host spheroid.
Quantitatively, M BH =M sph # 0:001 (e.g., Merritt &
Ferrarese 2001a). This relation is still controversial, how­
ever, both because the sensitivity of published searches is
correlated with bulge luminosity and because there is sub­
stantial scatter in M BH at fixed M sph . Recently, Ferrarese &
Merritt (2000) and Gebhardt et al. (2000) have shown that a
tighter correlation holds between the BH mass and the
velocity dispersion of the bulge. The two groups, however,
find two di#erent slopes of the correlation, M BH / # 5
# and
M BH / # 4
# , respectively (Tremaine et al. 2002). Clearly, any
correlation of BH and spheroid properties would have
important implications for theories of galaxy formation in
general and bulge formation in particular. Indeed, several
authors have shown that the slope of the M BH ­# * correla­
tion yields information on the formation process. If the for­
mation process is self­regulated, i.e., if the BH growth is
limited by radiation pressure, M BH / # 5
# (Silk & Rees
1998). Conversely, if the growth is regulated by ambient
conditions, M BH / # 4
# (Adams, Gra#, & Richstone 2001;
Cavaliere & Vittorini 2002). The uncertainties on the slope
of the M BH ­# * correlation do not allow one to distinguish
between the two cases. To solve this problem, more BH
mass measurements are needed in the low­mass range (10 6 --
10 7 M # ) where spiral galaxies are expected to fit.
These correlations are also very important to estimate
BH masses quickly and easily instead from very complex
dynamical and kinematic measurements. Therefore, the
quantities involved in the correlations must be measured as
carefully as possible in order to reduce the scatter to its
intrinsic value. For instance, it has been shown that with
careful estimates of the bulge luminosity the M BH ­L sph cor­
relation has the same scatter as the M BH ­# * correlation, in
contrast with previous claims (McLure &Dunlop 2002).
Given the low value of the BH mass in NGC 4041,
M BH < 2 # 10 7 M # , it is worthwhile to verify the relation of
this galaxy with the proposed correlations. The B magni­
tude from the RC3 catalog is 11.9, becoming 11.8 after
extinction correction (see NED). The morphological type is
Sbc/Sc and T ¼ 4:0 # 0:3. From Simien & de Vaucouleurs
(1986) the bulge­to­total luminosity ratio is '0.16, resulting
in Dm ¼ 2. Thus, the bulge magnitude is 13.8. The adopted
distance, 19.5 Mpc, corresponds to a distance modulus of
31.5; thus, the absolute bulge magnitude is #17.7, corre­
sponding to 1:8 # 10 9 L B# . The best fit of the M BH ­L B;sph
correlation gives M BH ¼ 0:8 # 10 8
ÏL B;sph =10 10 L B# ÷ 1:08
(Kormendy &Gebhardt 2001). Thus, the expected BHmass
in NGC 4041 would be 1:2 # 10 7 M # , in agreement with the
BH mass estimate or upper limit. The best fit of the M BH ­# *
correlation gives M BH ¼ 1:3 # 10 8
Ï# # =200 km s #1 ÷ 4:0
(Tremaine et al. 2002) or M BH ¼ 1:4 # 10 8
Ï# # =200 km
s #1 ÷ 4:8
(Merritt & Ferrarese 2001a). Using INTEGRAL/
WYFFOS at the William Herschel Telescope (WHT), we
have recently measured the stellar velocity dispersion in the
central 2 00 of NGC 4041 (D. Batcheldor et al. 2003, in prepa­
ration). With ## ¼ 95 # 5 km s #1 , the expected BH masses
are then 7 # 10 6 and 4 # 10 6 M # , both consistent with the
result from this paper.
Since the main goal of our project is to determine whether
or not spirals do in fact follow the M BH ­L sph and M BH ­# *
relations, it is important to observe also objects that, a pri­
ori, are expected to have marginally detectable or non­
detectable BHs. Indeed, even if the present measurement is
only an upper limit, this is still useful in ruling out the pres­
ence of unusually massive central BHs in late­type spiral
galaxies.
6. CONCLUSIONS
We presented HST/STIS spectra of the Sbc spiral galaxy
NGC 4041, which were used to map the velocity field of the
gas in its nuclear region. We detected the presence of a com­
pact (r ' 0>4 ' 40 pc), high surface brightness, circularly
rotating nuclear disk cospatial with a nuclear star cluster.
This disk is characterized by a rotation curve with a peak­
to­peak amplitude of #40 km s #1 and is systematically
blueshifted by #10--20 km s #1 with respect to the galaxy
systemic velocity.
We have analyzed the kinematical data assuming that the
stellar mass­to­light ratio is constant with radius and that
the gaseous disk is not warped, having the same inclination
as the large­scale galactic disk. We have found that, in order
to reproduce the observed rotation curve, a dark point mass
of Ï1 ×0:6
#0:7 ÷ # 10 7 M # is needed, very likely a supermassive
BH.
However, the observed blueshift suggests the possibility
that the nuclear disk could be dynamically independent.
Following this line of reasoning, we have relaxed the stan­
dard assumptions varying the stellar mass­to­light ratio and
the disk inclination. We have found that the kinematical
data can be accounted for by the stellar mass provided that
either the mass­to­light ratio is increased by a factor of #2
or the inclination is allowed to vary. This model resulted in
a 3 # upper limit of 6 # 10 6 M # on the mass of any nuclear
BH.
Combining the results from the standard and alternative
models, the present data only allow us to set an upper limit
of 2 # 10 7 M # to the mass of the nuclear BH.
If this upper limit is taken in conjunction with an esti­
mated bulge B magnitude of #17.7 and with a central stellar
velocity dispersion of '95 km s #1 , the putative BH in NGC
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 885

4041 is not inconsistent with both the M BH ­L sph and M BH ­
# * correlations.
Support for proposal GO­8228 was provided by NASA
through a grant from the Space Telescope Science Institute,
which is operated by the Association of Universities for
Research in Astronomy, Inc., under NASA contract NAS
5­26555. This work was partially supported by the Italian
Space Agency (ASI) under grants I/R/35/00 and I/R/112/
01 and by the Italian Ministry for Instruction, University,
and Research (MIUR) under grants Cofin00­02­35 and
Cofin01­02­02. We thank Peter Erwin for useful discussions
and the referee, Aaron Barth, for careful reading of the
manuscript, suggestions, and comments, which improved
the paper. This publication makes use of the LEDA data­
base and data products from the Two Micron All Sky
Survey, which is a joint project of the University of Massa­
chusetts and the Infrared Processing and Analysis Center/
California Institute of Technology, funded by the National
Aeronautics and Space Administration and the National
Science Foundation.
APPENDIX A
DERIVING THE LUMINOSITY DENSITY OF THE STARS FROM THE OBSERVED SURFACE
BRIGHTNESS PROFILE
The stellar luminosity density can be inferred by inverting the observed surface brightness profiles. Following van der Marel
& van den Bosch (1998), we assume an oblate spheroidal density distribution, which we parameterize as
#Ïm÷ ¼ # 0
m
r b
# # ##
1 ×
m
r b
# # 2
" # ##
: ÏA1÷
Here m is defined as m 2
¼ x 2
× y 2
× z 2 =q 2 , where xyz is a reference system with the x­y plane corresponding to the principal
plane of the potential and q is the intrinsic axial ratio of the mass distribution. If # is the mass­to­light ratio, the observed
surface brightness distribution is given by 14
# ¼
1
4##
Z ×1
#1
# ds ; ÏA2÷
where the integration is performed along the line of sight. It can be shown that
# m 0
Ï ÷ ¼
1
4##
q
q
Z ×1
m 02
# m 2
Ï ÷dm 2
###################### m 2
#m 02
Ï ÷
p ; ÏA3÷
where m 02 ¼ x 02 × y 02 =q 02 and x 0 y 0 is a reference system on the sky, with the x 0 ­axis aligned along the apparent major axis; q 0 is
the observed axial ratio of the isophotes, which is related to the intrinsic axial ratio of the mass distribution by
q 02 ¼ cos 2 i × q 2 sin 2 i, where i is the inclination of the line of sight (i ¼ 90 # is the edge­on case). The observed surface
brightness results from the convolution of #with the PSF P of the system (i.e., telescope and optics) and the detector pixels:
# app ÏX ; Y÷ ¼ Z X×DX
X#DX
Z Y×DY
Y#DY
Z ×1
#1
Z ×1
#1
# true Ïx; y÷P x 0 # x; y 0 # y
Ï ÷dx dy
# # dx 0 dy 0
4DXDY
; ÏA4÷
where X, Y is the center of an aperture with size 2DX # 2DY . The integration on dx dy can be directly carried out on the PSF
P, which is described by a sum of Gaussians:
PÏx; y÷ ¼
1
N K
X N
i¼1
k exp #
x 2
× y 2
2# 2
i
# # ; ÏA5÷
with NK ¼ P N
i¼1 2## 2
i k i . Then # app is given by
# app ÏX ; Y÷ ¼ Z ×1
#1
Z ×1
#1
#Ïx; y÷ true K X # DX # x; X × DX # x; Y # DY # y; Y × DY # y
Ï ÷
dx dy
4DXDY
; ÏA6÷
where K is the convolution kernel,
KÏX 1 ; X 2 ; Y 1 ; Y 2 ÷ ¼
1
N K
X N
i¼1
2## 2
i
k i
4
E
X 1
# i ### 2
p
# # # E
X 2
# i ### 2
p
# #
# # E
Y 1
# i ### 2
p
# # # E
Y 2
# i ### 2
p
# #
# # ; ÏA7÷
where E is the complementary error function. The integration is then carried out numerically using the Gauss Legendre
approximation.
14 Note the 1/4# factor. This derives from the fact that #/V is a density (luminosity per unit volume) and # is a surface brightness (luminosity per unit area
per unit solid angle).
886 MARCONI ET AL. Vol. 586

The best­fitting model is determined by minimizing the reduced # 2
red , defined as
# 2
red ¼
1
N d
X N
i¼1
log # # log #m r i ; #r i ; p 1 ; . . . ; p m
Ï ÷
# log # i
# # 2
; ÏA8÷
where i ¼ 1; . . . ; N indicates a data point with surface brightness # i # ## i at radius r i . #m Ïr ; #r i ; p 1 ; . . . ; p m ÷ is the model
surface brightness, averaged over radii r # #r # r # r i × #r i , which is a function of m free parameters p 1 ; . . . ; pm . N d ¼ N #m
is the number of degrees of freedom. The # 2
red is minimized to determine the m free parameters using the downhill simplex
algorithm by Press et al. (1992).
When the mass density has been determined as described above, the circular velocity in the principal plane is given by the
relation (e.g., Binney & Tremaine 1987)
V 2
c Ïr÷ ¼ 4#Gq 2 Z r
0
#Ïm 2
÷m 2 dm
#################### r 2
#m 2 e 2
p ; ÏA9÷
where e is the eccentricity related to q by q 2
¼ 1 # e 2 . Similarly, the mass enclosed within the homoeoid defined by m < m 0 is
MÏm 0 ÷ ¼ 4#q Z m 0
0
# m 2
# # m 2 dm : ÏA10÷
APPENDIX B
THE ROTATION CURVE MODEL
The velocity field along the line of sight, # vv, can be easily computed in the case of a circularly rotating thin disk. Let XY be a
reference frame in the plane of the sky with Y­axis fixed along the direction given by the slit position (hereafter called the slit
reference frame). Consider a reference frame X d Y d in the sky with the X d ­axis aligned along the disk line of nodes and the origin
coincident with the disk center (hereafter called the `` disk reference frame ''; see Fig. 20). The origin of XY is chosen in such a
way that the disk center has coordinates x ¼ b and y ¼ 0 in the slit reference frame. Then a given disk point P with coordinates
(x, y) in the slit reference frame has coordinates (x d , y d ) in the disk reference frame given by
x d ¼ Ïx # b÷ sin # × y cos # ;
y d ¼ #Ïx # b÷ cos # × y sin # : ÏB1÷
If the disk has an inclination angle i (i ¼ 0 in the face­on case), then P is at the disk radius r given by
r 2
¼ x 2
d × # y d
cos i # 2
: ÏB2÷
Fig. 20.---Geometry of the disk
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 887

The circular velocity of P, in the case of a spherical mass distribution, is then given by
V c Ïr÷ ¼ r
d#
dr# # # # # # # #
1=2
¼
GMÏr÷
r
# # 1=2
; ÏB3÷
where M(r) is the enclosed mass at radius r, a constant value M BH in case of a point mass. In the case of an oblate spheroidal
mass distribution, V c (r) is given by equation (A9).
The velocity component along the line of sight is finally given by
# vv ¼ V sys # V c Ïr÷ sin i x d
r ¼ V sys # GM BH
Ï ÷ 0:5 sin i x d
r 1:5
; ÏB4÷
with the latter expression describing the simple case of a point mass M BH .
This velocity field # vv has to be convolved with the instrumental response in order to simulate the observed quantities. Below
we show how this can be done for any velocity field # vv.
Consider again the reference frame XY in the plane of the sky defined above. The light distribution can be written as
#Ïx 0 ; y 0 ; v÷ ¼ IÏx 0 ; y
÷# Ïv # # vv÷ ; ÏB5÷
where IÏx 0 ; y 0 ÷ is the total intensity at the point (x 0 , y 0 )
and# Ïv # # vv÷ is the intrinsic line profile centered at the velocity along the
line of sight # vvÏx 0 ; y 0 ÷. After passing through the telescope, the light gets convolved with the instrumental PSF PÏx # x 0 ; y # y ÷,
and at point (x, y) in the focal plane the light distribution is described by
#Ïx; y; v÷ ¼ Z Z ×1
#1
dx 0 dy I x 0 ; y 0
Ï
÷# Ïv # # vv÷P x # x 0 ; y # y 0
Ï ÷ : ÏB6÷
As already defined, the slit in the focal plane is aligned along the Y­axis, and let its center cross the X­axis at x 0 . Let the slit
width be 2Dx. Past the slit, the light falls into the detector plane, where the spatial coordinate x and the velocity v are combined
into one single detector coordinate w ¼ v × kÏx # x 0 ÷, which we identify with the observed velocity (we denote it w here in con­
trast to the intrinsic velocity denoted by v). Here the coe#cient k is given by lDw/Dy, where l, the anamorphic magnification,
accounts for the di#erent scales on the dispersion and slit directions. In the case of STIS, the scale along the slit is 0>05071 and
along dispersion is 0>05477; thus, l ¼ 1:080 (Bowers & Baum 1998). The light distribution in the detector plane, #Ïw; y÷, is
calculated by ousting v and integrating the light contribution across the slit
#Ïw; y÷ ¼ Z x 0 ×Dx
x 0 #Dx
dx # x; y; w # k x # x 0
Ï ÷
½# ¼ Z x 0 ×Dx
x 0 #Dx
dx Z Z ×1
#1
dx dy
0# w # # vv x 0 ; y
Ï ÷ # k x # x 0
Ï ÷
½# I x 0 ; y 0
Ï ÷P x # x 0 ; y # y
Ï ÷ : ÏB7÷
Note the component kÏx # x 0 ÷ in the velocity profile, which is the `` spurious '' velocity for light entering o# from the slit
center. Properties of thus defined two­dimensional light distribution with the spurious velocity shift were investigated in detail
by Maciejewski & Binney (2001).
The detector integrates over finite pixel sizes; therefore, we calculate the expected line fluxes, average velocities, and widths
for the line profile that is obtained by integrating the light distribution on the detector plane #Ïw; y÷ over the width 2Dy of the
jth pixel along the slit and convolving it with the shape of the pixel in the dispersion direction: a top hat of width 2Dw. Thus,
the expected line profile in the detector is
~
## Ïw÷ ¼ Z y ×Dy
y #Dy
dy Z w×Dw
w#Dw
dw# w 0 ; y
Ï ÷ ; ÏB8÷
where the jth pixel has the coordinate y j along the slit. Note that the observed intensities are measured at discrete values of w
corresponding to the pixel centers. In order to calculate the expected line fluxes, average velocities, and widths, one has to
evaluate moments of ~
## j Ïw÷:
Z ×1
#1
dww n ~
## j Ïw÷ ¼ Z y ×Dy
y j #Dy
dy Z x 0 ×Dx
x 0 #Dx
dx ZZ ×1
#1
dx 0 dy 0 P Z ×1
#1
dww n Z w×Dw
w#Dw
dw
0# w 0 # w 0
Ï ÷ ; ÏB9÷
where we used abbreviations P ¼ IÏx 0 ; y 0 ÷PÏx # x 0 ; y # y ÷ and w 0 ¼ # vvÏx ; y 0 ÷ × kÏx # x 0 ÷. The last two integrals can be
simplified by inverting the order of integration:
Z ×1
#1
dww n Z w×Dw
w#Dw
dw
0# w 0 # w 0
Ï ÷ ¼ ZZ ×1
#1
dw dw 0 w
n# w 0 # w 0
Ï ÷H w # w 0
Ï ÷ ¼ Z ×1
#1
dw# Ïw # w 0 ÷ Z w×Dw
w#Dw
dw 0 w 0n
¼ Z ×1
#1
dw# Ïw # w 0 ÷ Ïw × Dw÷ n×1
# Ïw # Dw÷ n×1
n × 1
;
888 MARCONI ET AL. Vol. 586

where HÏw # w 0 ÷ is 1 for w # w 0
j j < Dw and 0 otherwise. This leads to the following formulae for the moments of ~
## j Ïw÷:
Z ×1
#1
~
## j Ïw÷dw ¼ 2Dw Z x 0 ×Dx
x 0 #Dx
dx Z y j ×Dy
y #Dy
dy ZZ ×1
#1
dx 0 dy 0 P ; ÏB10÷
Z ×1
#1
w ~
## j Ïw÷dw ¼ 2Dw Z x 0 ×Dx
x 0 #Dx
dx Z y j ×Dy
y #Dy
dy ZZ ×1
#1
dx 0 dy 0 w 0 P ; ÏB11÷
Z ×1
#1
w 2 ~
## j Ïw÷dw ¼ 2Dw Z x 0 ×Dx
x 0 #Dx
dx Z y j ×Dy
y j #Dy
dy ZZ ×1
#1
dx 0 dy 0 w 2
0 ×
Dw
Ï ÷ 2
3 × # 2
" # P : ÏB12÷
To obtain them, we assumed that the intrinsic velocity
profile# Ïv÷ is bound and symmetric, i.e.,
Z ×1
#1#
Ïv÷dv ¼ 1 Z ×1
#1
v# Ïv÷dv ¼ 0 Z ×1
#1
v
2# Ïv÷dv ¼ # 2 ;
where # 2 is the intrinsic velocity dispersion.
Equations (B10), (B11), and (B12) can be used to compute the expected line fluxes, average velocities, and widths. For
example, in the most simple case of a constant line intensity I ¼ const and velocity field # vv ¼ const, one can write
hv j i ¼ R ×1
#1 w ~
## j Ïw÷dw
R ×1
#1
~
## Ïw÷dw ¼ # vv ;
hv 2
j i ¼ R ×1
#1 w 2 ~
## j Ïw÷dw
R ×1
#1
~
## Ïw÷dw ¼ # vv 2
× # 2
× ÏDw÷ 2
3 × ÏkDx÷ 2
3
: ÏB13÷
Then the expected velocity dispersion, given by # 2
j ¼ hv 2
j i # hv j i 2
, is, understandably, larger than the intrinsic velocity
dispersion #. This is due to the convolution with the pixel size and slit width, which add quadratically. However, note that
while Dw and # enter the integral given by equation (B12) as constants, w 0 is a linear combination of the spurious velocity shift
and the intrinsic velocity. These two contributions can cancel, resulting in the expected line profile being broadened by the
pixel size only, and not by the width of the slit. This implies that the wide slit can probe the disk on the scale of the pixel size
rather than the slit width. Maciejewski & Binney (2001) explore consequences of this finding.
In order to compute the model given by equations (B10), (B11), and (B12), we create a grid in x and y with sampling given
by # PSF /n. Here # PSF is the spatial rms of the PSF and n is the subsampling factor. We have verified that the optimal sub­
sampling factor used is n ¼ 3 since larger values do not produce any appreciable di#erences in the final results. The PSF used
is the one generated by TinyTim (V6.0; Krist &Hook 1999) at 6700 A š . Convolution with the PSF is done using the fast Fourier
transform algorithm (Press et al. 1992). Following Barth et al. (2001), we have also introduced the CCD scattering function
(Leitherer et al. 2001), but as already noticed by them, it does not have any appreciable e#ect in the final results.
We compare models to the observed spectrum, which is essentially an array of intensities # ij , after the observed line profile
is derived in the following way: to the sequence of intensities # ij for a given row j along the slit, we fit a baseline and a
continuous analytical function ~
## j Ïw÷ obs , which we interpret as the observed equivalent of the expected line profile ~
## j Ïw÷
(eq. [B8]). The best­fitting model is determined by minimizing the reduced # 2
red defined as
# 2
red ¼
1
N d
X 3
k¼1
X N k
j¼1
v kj # v j
# # k p 1 ; . . . ; pm
Ï ÷
#v kj
" # 2
×
W kj # W j
# # k Ïp 1 ; . . . ; pm ÷
#W kj
" # 2
8 < :
9 = ;
; ÏB14÷
where the index k ¼ 1, 3 indicates the slit position and j ¼ 1, N k counts pixels along the slit. Here the characteristics of the
model are as follows. The velocity in the jth row along the slit hv j i k and the FWHM of the velocity profile hW j i k are calculated
directly from equation (B13), now clearly for variable intensity and velocity field. They both are functions of m free parameters
p 1 ; . . . ; pm , which are determined by # 2
red minimization. The FWHM is calculated from the expected velocity dispersion # j after
assuming a Gaussian line profile. The observed velocities (v ki # #v ki ) and velocity dispersions (W ki # #W ki ) are also derived
from equation (B13), but after ~
## j Ïw÷ has been replaced by ~
## j Ïw÷ obs defined above. N d ¼ P 3
k¼1 2N k #m is the number of
degrees of freedom.
The # 2
red is minimized to determine the m free parameters using the downhill simplex algorithm by Press et al. (1992). In
order to apply statistical methods when # 2
red is much larger than 1, we follow Barth et al. (2001) and rescale errors as
#v 0 ki
Ï ÷ 2
¼ #v 2
ki × #V 2 ; ÏB15÷
where #V is a `` systematic '' error determined such that the resulting # 2
c of the `` best '' model is 1. The fits presented in this
paper have # 2
red # 1, and the error rescaling was not performed.
No. 2, 2003 BLACK HOLE AT CENTER OF NGC 4041 889

REFERENCES
Adams, F. C., Gra#, D. S., &Richstone, D. O. 2001, ApJ, 551, L31
Athanassola, E. 1992, MNRAS, 259, 345
Ba#a, C., et al. 2001, A&A, 378, 722
Bahcall, J. N., &Wolf, R. A. 1976, ApJ, 209, 214
Barth, A. J., Sarzi, M., Rix, H., Ho, L. C., Filippenko, A. V., & Sargent,
W. L. W. 2001, ApJ, 555, 685
Binney, J., &Mamon,G.A. 1982, MNRAS, 200, 361
Binney, J., & Merrifield, M. 1998, Galactic Astronomy (Princeton:
Princeton Univ. Press)
Binney, J., & Tremaine, S. 1987, Galactic Dynamics (Princeton: Princeton
Univ. Press)
Bo ˜ ker, T., Laine, S., van der Marel, R. P., Sarzi, M., Rix, H., Ho, L. C., &
Shields, J. C. 2002, AJ, 123, 1389
Bower, G. A., et al. 1998, ApJ, 492, L111
Bowers, C., & Baum, S. 1998, STIS Instrument Science Rep. 98­23
(Baltimore: STScI)
Braatz, J. A., Wilson, A. S., &Henkel, C. 1997, ApJS, 110, 321
Bruzual, A. G., &Charlot, S. 1993, ApJ, 405, 538
Burbidge, E. M., &Burbidge, G. R. 1962, ApJ, 135, 694
Cappellari, M., Verolme, E. K., van der Marel, R. P., Kleijn, G. A. V.,
Illingworth, G. D., Franx, M., Carollo, C. M., & de Zeeuw, P. T. 2002,
ApJ, 578, 787
Cardelli, J. A., Clayton, G. C., &Mathis, J. S. 1989, ApJ, 345, 245
Carollo, C. M., Franx, M., Illingworth, G. D., & Forbes, D. A. 1997, ApJ,
481, 710
Carollo, C. M., Stiavelli, M., &Mack, J. 1998, AJ, 116, 68
Carollo, C. M., Stiavelli, M., Seigar, M., de Zeeuw, P. T., & Dejonghe, H.
2002, AJ, 123, 159
Cavaliere, A., &Vittorini, V. 2002, ApJ, 570, 114
Chokshi, A., &Turner, E. L. 1992, MNRAS, 259, 421
de Jong, R. S. 1996, A&A, 313, 377
de Vaucouleurs, G., de Vaucouleurs, A., Corwin, J. R., Buta, R. J., Paturel,
G., & Fouque, P. 1991, Third Reference Catalogue of Bright Galaxies
(New York: Springer)
Devereux, N., Ford, H., & Jacoby, G. 1997, ApJ, 481, L71
Ferrarese, L., Ford, H. C., & Ja#e, W. 1996, ApJ, 470, 444
Ferrarese, L., &Merritt, D. 2000, ApJ, 539, L9
Ford, H. C., Tsvetanov, Z. I., Ferrarese, L., & Ja#e, W. 1998, in The
Central Regions of the Galaxy and Galaxies, ed. Y. Sofue (Dordrecht:
Kluwer), 377
Ford, H. C., et al. 1994, ApJ, 435, L27
Gebhardt, K., et al. 2000, ApJ, 539, L13
Genzel, R., Pichon, C., Eckart, A., Gerhard, O. E., & Ott, T. 2000,
MNRAS, 317, 348
Giovanelli, R., &Haynes, M. P. 2002, ApJ, 571, L107
Greenhill, L. J., Moran, J. M., &Herrnstein, J. R. 1997, ApJ, 481, L23
Harms, R. J., et al. 1994, ApJ, 435, L35
Ho, L. C., Filippenko, A. V., & Sargent, W. L. W. 1997, ApJ, 487, 568
Ho, L. C., Sarzi, M., Rix, H., Shields, J. C., Rudnick, G., Filippenko,
A. V., &Barth, A. J. 2002, PASP, 114, 137
Ja#e, W., Ford, H., Ferrarese, L., van den Bosch, F., & O'Connell, R. W.
1996, ApJ, 460, 214
Kinney, A. L., Calzetti, D., Bohlin, R. C., McQuade, K., Storchi­
Bergmann, T., & Schmitt, H. R. 1996, ApJ, 467, 38
Knapen, J. H., Beckman, J. E., Heller, C. H., Shlosman, I., & de Jong, R. S.
1995, ApJ, 454, 623
Kormendy, J., &Gebhardt, K. 2001, in 20th Texas Symposium on Relativ­
istic Astrophysics, ed. J. C. Wheeler &H.Martel (Melville: AIP), 363
Kormendy, J., &Richstone, D. 1995, ARA&A, 33, 581
Kormendy, J., et al. 1996, ApJ, 459, L57
Krist, J., &Hook, R. 1999, The Tiny Tim User's Guide (Baltimore: STScI)
Leitherer, C., et al. 1999, ApJS, 123, 3
---------. 2001, STIS Instrument Handbook, Version 5.1 (Baltimore: STScI)
Lynden­Bell, D. 1969, Nature, 223, 690
Macchetto, F., Marconi, A., Axon, D. J., Capetti, A., Sparks, W., &Crane,
P. 1997, ApJ, 489, 579
Maciejewski, W., & Binney, J. 2001, MNRAS, 323, 831
Maiolino, R., &Rieke, G. H. 1995, ApJ, 454, 95
Maraston, C. 1998, MNRAS, 300, 872
Marconi, A., Capetti, A., Axon, D. J., Koekemoer, A., Macchetto, D., &
Schreier, E. J. 2001, ApJ, 549, 915
Marconi, A., & Salvati M., N. 2002, in ASP Conf. Ser. 258, Issues in Unifi­
cation of Active Galactic Nuclei, ed. R. Maiolino, A. Marconi, &
N. Nagar (San Francisco: ASP), 217
Marconi, A., Schreier, E. J., Koekemoer, A., Capetti, A., Axon, D.,
Macchetto, D., &Caon, N. 2000, ApJ, 528, 276
McLure, R. J., &Dunlop, J. S. 2002, MNRAS, 331, 795
Merritt, D. 1997, AJ, 114, 228
Merritt, D., &Ferrarese, L. 2001a, ApJ, 547, 140
---------. 2001b, in ASP Conf. Ser. 249, The Central Kiloparsec of Starbursts
and AGN: The La Palma Connection, ed. J. H. Knapen, J. E. Beckman,
I. Shlosman, & T. J. Mahoney (San Francisco: ASP), 335
Miyoshi, M., Moran, J., Herrnstein, J., Greenhill, L., Nakai, N., Diamond,
P., & Inoue, M. 1995, Nature, 373, 127
Moriondo, G., Giovanardi, C., &Hunt, L. K. 1998, A&A, 339, 409
Osterbrock, D. E. 1989, Astrophysics of Gaseous Nebulae and Active
Galactic Nuclei (Mill Valley: University Science Books)
Panagia, N. 1973, AJ, 78, 929
Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 1992,
Numerical Recipes in C: The Art of Scientific Computing (2d ed.;
Cambridge: Cambridge Univ. Press)
Pryor, C., & Meylan, G. 1993, in ASP Conf. Ser. 50, Structure and
Dynamics of Globular Clusters, ed. S. G. Djorgovski & G. Meylan (San
Francisco: ASP), 357
Richstone, D. 1998, BAAS, 30, 1140
Sakamoto, K., Okumura, S. K., Ishizuki, S., & Scoville, N. Z. 1999, ApJS,
124, 403
Salpeter, E. E. 1955, ApJ, 121, 161
Scalo, J. M. 1986, Fundam. Cosmic Phys., 11, 1
Schreier, E. J., et al. 1998, ApJ, 499, L143
Shaw, M., Axon, D., Probst, R., &Gatley, I. 1995, MNRAS, 274, 369
Silk, J., &Rees, M. J. 1998, A&A, 331, L1
Simien, F., & de Vaucouleurs, G. 1986, ApJ, 302, 564
Sootan, A. 1982, MNRAS, 200, 115
Statler, T. S. 1987, ApJ, 321, 113
Tremaine, S., et al. 2002, ApJ, 574, 740
van der Marel, R. P., Cretton, N., de Zeeuw, P. T., & Rix, H. 1998, ApJ,
493, 613
van der Marel, R. P., & van den Bosch, F. C. 1998, AJ, 116, 2220
Verdoes Kleijn, G. A., van der Marel, R. P., Carollo, C. M., & de Zeeuw, P.
T. 2000, AJ, 120, 1221
Whitmore, B. 1995, in Calibrating Hubble Space Telescope: Post Servicing
Mission, ed. A. Koratkar &C. Leitherer (Baltimore: STScI), 269
Wilkins, T. W., & Axon, D. J. 1992, in ASP Conf. Ser. 25, Astronomical
Data Analysis, Software and Systems I, ed. D. M. Worral, C. Biemesder­
fer, & J. Barnes (San Francisco: ASP), 427
Winge, C., Axon, D. J., Macchetto, F. D., Capetti, A., &Marconi, A. 1999,
ApJ, 519, 134
890 MARCONI ET AL.