The Dark Matter Distributions in Low-mass Disk Galaxies. I. Hα Observations Using the Palomar Cosmic Web Imager

Explicit details and criteria used to select our targets can be found in Truong et al. (2017) but are summarized here. We first selected low-mass galaxies by considering those with H i line widths W₂₀ ≲ 200 km s⁻¹ and absolute magnitude M_B ≳ −18. Second, the galaxies were required to have an inclination between 30° and 70° and a diameter larger than 1'. Third, we considered those that were late type and generally symmetric and regular in shape, with no obvious central bar or distortion. We also limited our search to galaxies within 30 Mpc. Given our expected resolution of about 2'', the innermost radius we could hope to observe is around 1'', which at the maximum distance of 30 Mpc corresponds to 145 pc, and to 78 pc at our sample's median distance of 16 Mpc.

Our final sample consists of 26 galaxies with typical luminosities of 10⁹ L_☉ in the r band and rotation velocities of around 50–140 km s⁻¹. Two galaxies, UGC 3371 and UGC 11891, are classified as low surface brightness (LSB) galaxies. Twenty-two of our galaxies overlap with the 26 targeted for CO observations by Truong et al. (2017), though of their 26 only 14 had sufficiently bright and extended CO to produce viable rotation curves, and of those we overlap with 11 galaxies. We compare the Hα and CO velocity fields in Section 5 and will examine the inner dark matter slopes in Paper II.

The names and properties of our 26 Hα galaxies can be found in Table 1. To determine distances, we found that 13 of our galaxies have Tully–Fisher distances from Tully et al. (2013). Of the remaining galaxies, we used Tully–Fisher distances from Gavazzi et al. (1999) for NGC 4376 and NGC 4396 and from Theureau et al. (2007) for UGC 8516. For the rest, we used the Hubble flow distance with H₀ = 73.0 and corrected for large-scale flows toward the Virgo Cluster, the Shapley Supercluster, and the Great Attractor, as listed in the NASA/IPAC Extragalactic Database (NED).

Table 1.  Sample Galaxy Properties

Name	R.A.	Decl.	Distance	Vsys	$\mathrm{log}({L}_{r})$	$\mathrm{log}({L}_{\mathrm{Ch}2})$
			(Mpc)	(km s−1)
NGC 746	29.46294	44.91834	12 (b)	702	8.97	⋯
NGC 853*	32.92178	−9.30524	20 (b)	1518	9.48	⋯
NGC 949*	37.70259	37.13668	10 (a)	605	9.28	⋯
NGC 959	38.09967	35.49463	10 (a)	593	9.01	9.49
NGC 1012*	39.81184	30.15114	14 (b)	969	9.43	⋯
NGC 1035*	39.87136	−8.13283	16 (a)	1240	9.49	10.26
NGC 2644	130.38355	4.98104	30 (b)	1932	9.75	10.33
NGC 2976	146.81385	67.91668	4 (a)	1	8.35	9.58
NGC 3622	170.05169	67.24187	23 (a)	1317	9.46	9.9
NGC 4376	186.32526	5.74134	24 (c)	1132	9.39	9.82
NGC 4396	186.49512	15.67166	16 (c)	−120	9.31	9.89
NGC 4451*	187.16880	9.25915	26 (a)	849	9.74	10.31
NGC 4632*	190.63331	−0.08242	14 (a)	1702	9.48	10.15
NGC 5303*	206.93757	38.30463	28 (b)	1421	9.66	10.34
NGC 5692*	219.57549	3.41031	27 (b)	1585	9.57	⋯
NGC 5949*	232.00273	64.76313	13 (a)	434	9.35	9.94
NGC 6106*	244.69664	7.41082	24 (a)	1443	9.81	10.38
NGC 6207*	250.76550	36.83213	16 (a)	843	9.67	10.23
NGC 6503	267.36004	70.14434	6 (a)	29	9.47	10.1
NGC 7320	339.01423	33.94826	14 (b)	771	9.14	9.63
UGC 01104	23.17725	18.31687	10 (b)	680	8.42	8.76
UGC 3371	89.15165	75.31691	15 (b)	814	8.88	⋯
UGC 4169	120.63655	61.38821	30 (a)	1603	9.62	10.04
UGC 8516	202.96847	20.00127	22 (d)	1021	9.27	9.77
UGC 11891	330.89073	43.74893	8 (a)	463	9.11	⋯
UGC 12009	335.66767	37.97735	20 (b)	1221	9.26	⋯

Note. Galaxies marked with an asterisk are those that overlap the viable CO sample from Truong et al. (2017). Distances are marked to indicate sources as follows (see Section 2.1): a—Tully et al. (2013); b—Hubble flow distance; c—Gavazzi et al. (1999); d—Theureau et al. (2007). Values of V_sys are taken from our analysis (see Section 4) with typical random uncertainties of 1–2 km s⁻¹. The remaining parameters were derived from photometry, as described in Section 3. For galaxies with a tabulated IRAC Channel 2 luminosity, we will preferentially use infrared photometry in our analysis.

Download table as:  ASCIITypeset image

Since we ultimately would like to probe the inner dark matter distributions of these galaxies, our observations require an instrument with a sufficiently high spatial and spectral resolution to measure gas kinematics near their centers. We chose to use the PCWI, an integral field spectrograph at Palomar Observatory. It is mounted at the Cassegrain focus of the Hale 5 m telescope and has a resolving power of R ∼ 5000 (Matuszewski et al. 2010). An image slicer divides the approximately 60'' × 40'' field of view into 24 slices, each approximately 40'' × 26. The spatial resolution is limited by the seeing on the 40'' axis (along each slice) and by the 26 slice width along the perpendicular axis, although this limit can be improved by dithering between exposures, as described below.

The 26 galaxies in our sample were observed during 17 nights from 2013 to 2015. To locate each galaxy for the exposures, a nearby bright star was centered in the field of view, and then the telescope was offset to the position of the galaxy. For the pointing covering the center of each galaxy, where the highest possible spatial resolution is necessary, we obtained consecutive exposures with the telescope dithered by 13, half the slice width, in the direction perpendicular to the slices. By combining these interlaced exposures as described in the following section, we aim to recover some of the spatial information that would otherwise be lost due to the slices undersampling the seeing. The mean seeing for each observation date ranges from 13 to 27, with an overall sample mean of 18.

The majority of our galaxies were observed using the nod and shuffle technique. This method, developed by Sembach & Tonry (1996), allows for high-precision sky subtraction by regularly taking sky exposures close in time and position to the target, nodding the telescope without reading out the detector. A mask is used to block the top and bottom portions of the detector so that only the middle third is illuminated. For each galaxy, a blank sky position is selected and a nearby star is chosen for guiding in both the target and sky regions. The telescope is nodded between the target and sky fields periodically, with a coordinated shuffle of charge on the CCD so that the illuminated portion of the detector (the middle third) contains the appropriate target or sky exposure. The nod and shuffle technique greatly reduces read noise while allowing us to sample the temporal variations in the sky background. The process ultimately produces two spectra, one of the target with sky and one of only the sky, which can easily be subtracted to obtain the target spectrum.

When using this technique, we nodded between sky and target every 2–3 minutes for 10–20 minutes of total exposure time. We also obtained calibration exposures using the internal quartz and arc lamps regularly throughout the night. Since nod and shuffle requires the use of a nearby guide star, for targets where no such star was present or conditions prevented us from seeing the star, we could not use the technique. This was the case for five of our galaxies, and for these we instead took separate sky exposures to use for subtraction.

In most cases, the galaxy was larger than the PCWI field of view, so exposures had to be taken at several positions to cover it. These groups of exposures, referred to as tiles, are later stitched together to create a single image (see Section 2.3).

Details about the observations of each galaxy can be found in Table 2.

Individual exposures were reduced by using the PCWI data reduction pipeline (Martin et al. 2014). The pipeline performs bias subtraction, sky subtraction for exposures using nod and shuffle, flat-fielding, derivation of the wavelength solution, and rectification. This process produces a data cube for each exposure.

The data cubes are then further processed through a custom set of IDL routines designed to produce maps of the Hα intensity, velocity, and velocity uncertainty. A zero-point correction to the wavelength solution is determined for each image using night-sky emission lines. This correction accounts for instrumental flexure between the science observations and the calibrations. For observations without the nod and shuffle technique, the intensity of skylines and flexure changes over time, which precludes a straightforward subtraction of the sky exposure. Instead, we shift the sky spectra slightly in wavelength and rescale the intensity, adjusting until the residuals are minimized. During the reduction, the PCWI slices are binned along their length to create spaxels that are 13 × 26. The Hα line in each spaxel is then fit with a Gaussian profile, which we found to adequately describe the emission.

Our procedure for acquiring the galaxies during observation was not accurate to arcsecond precision, so the coordinate map produced by the PCWI pipeline must be adjusted. To do this, we compare the distribution of Hα in each data cube to a narrowband image. More details about the imaging observations are given in Section 3. The Hα intensity map produced from the PCWI data cube was cross-correlated in 2D with the narrowband Hα image to derive the absolute coordinates.

In cases where it took more than one pointing to cover the galaxy, the 2D flux and velocity fields are stitched together to form a single field. The process used is similar to the DRIZZLE method, introduced by Fruchter & Hook (2002). The separate exposures are drizzled onto a grid of 132 pixels, which is approximately half of a slice. The dithering method used in the central pointing similarly steps by half of a slice, which allows for improved resolution in the center.

Various steps in the data reduction process for a sample galaxy (NGC 5949) are displayed in Figure 1, while an r-band image, the final Hα flux map, and the velocity field for the same galaxy can be found in Figure 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Due to the angular resolution and sensitivity of PCWI, we generally detect Hα uniformly throughout the disk, ensuring highly detailed velocity fields. The observations sample almost the full range of azimuth at most locations, which is crucial for measuring noncircular motions. As we average over annuli, these patches will not have a significant effect on the derived rotation curves. The Hα flux maps and velocity fields for all galaxies can be found in the bottom middle and bottom right panels of the Appendix figures.

For each galaxy we require an r-band image and a narrowband image around the Hα line. Continuum emission is removed from the narrowband image to produce the maps of Hα that are used to register the PCWI exposures to absolute coordinates. We observed 21 of our galaxies using SPICAM on the 3.5 m telescope at Apache Point Observatory (APO) over three nights: 2013 October 5, 2014 March 5, and 2014 November 27. For four of the remaining galaxies (NGC 2976, NGC 6503, UGC 3371, UGC 11891), the r-band and Hα data were taken from the Local Volume Legacy Project. These observations were obtained with the Steward Observatory Bok 2.3 m telescope on Kitt Peak, as described in Kennicutt et al. (2008). For the remaining galaxy (NGC 2644), we used images taken from the 1.0 m Jacobus Kapteyn Telescope on La Palma as part of the survey done by James et al. (2004). We masked other galaxies and stars in the field of view by using Source Extractor (Bertin & Arnouts 1996) to locate them and place elliptical masks. The r-band images can be seen in the bottom left panels of the Appendix figures, while contours from the narrowband images are overlaid in the bottom middle panels.

For 18 of the galaxies, we have additional images from Channel 2 (4.5 μm) of Spitzer's Infrared Array Camera (IRAC; Fazio et al. 2004), obtained from the Spitzer Heritage Archive and the Spitzer Survey of Stellar Structure in Galaxies (Sheth et al. 2010). For the analysis described below, if a galaxy had both infrared and r-band observations, we preferentially used the infrared images to define the parameters, as they more closely trace the stellar mass distribution (Korsaga et al. 2018, 2019).

We use two independent procedures to estimate the center, ellipticity, and PA of our galaxies. Comparison of the results both helps ensure accuracy by verifying agreement and better addresses special cases where one method may be better suited. We first use a program called find_galaxy that is included in the multi-Gaussian expansion (MGE) package presented by Cappellari (2002). We then use the ELLIPSE task in IRAF to compare and verify the MGE values. If the galaxy has an obvious nucleus, we chose the coordinates of the nucleus to define the center rather than using a fitting procedure.

The find_galaxy method finds the center, ellipticity, and PA by computing a best-fit ellipse using the moments of the pixels above a given surface brightness threshold (the results were not very sensitive to the threshold). The IRAF ELLIPSE procedure instead fits ellipses to a sequence of isophotes. This allows us to test the variation of the ellipticity with radius and is much less susceptible to irregularities in the central region, so it is more appropriate for galaxies with bright or chaotic centers. In the cases of NGC 746, NGC 853, and NGC 1012, we chose to use ELLIPSE values over find_galaxy values because of their irregular central regions. This is generally a consequence of using r-band data, which is not as smooth as the 4.5 μm data, rather than a reflection of abnormal behavior at the center. For all other cases, the ELLIPSE and find_galaxy values agreed to ∼0.02 in ellipticity and 1° in PA. These estimates are needed to extract rotation curves from the PCWI data (see Section 4). In Paper II, we will analyze the surface brightness profiles obtained from these images in order to estimate the contribution of stellar mass to the rotation curve.

DiskFit requires, as input, the velocity field and its associated error. The formal velocity measurement uncertainties associated with fitting a Gaussian to an observed emission line are very small. Given small-scale interstellar medium (ISM) motions that may vary randomly across a galaxy, they are also not representative of the true uncertainty in the rotation velocity as a function of radius. To account for random motions within the ISM, a value for Δ_ISM is added to the velocity errors in quadrature when DiskFit is run. We chose Δ_ISM = 6 km s⁻¹, as this was the value that produced a median reduced χ² of 1 for the sample. This value is in agreement with that adopted by Spekkens & Sellwood (2007) and Kuzio de Naray et al. (2012).

DiskFit fits velocities at specified radii. We chose to start at an inner radius of r = 1 pixel and then increase using a bin size of 2 pixels until we covered the entire Hα field. Since 1 pixel is 132, r = 1 approximately marks the smallest rescaled radius at which we can obtain measurements, given our typical seeing of 18. The cutoff point for the outermost radius varies by galaxy and was determined to be where the number of points contained in the annulus fell below 4.

As mentioned above, DiskFit also requires initial estimates of the shape, position, and systemic velocity of the galaxy, and these are either kept fixed or allowed to vary in the fitting procedure. For all parameters except V_sys, we used values that were derived from our photometry, as described in Section 3. For V_sys, our initial estimate was the value given by SIMBAD and was allowed to vary during fitting. Initially we allowed the disk geometry (center, ellipticity, and PA) to vary in the kinematic fit for all of the galaxies in our sample. For most galaxies this produces estimates very close to the photometric values; however, for about a third of our sample DiskFit ultimately fit values that were significantly different. Table 3 shows a comparison of the photometric and kinematic estimates of the ellipticity, PA, and center from DiskFit for each galaxy, while comparison plots for each parameter are provided in Figure 3. The most common issue was that the kinematic inclination was far too low. The corresponding rotation curves reflected these discrepancies, displaying highly nonmonotonic behavior and large radial motions or dramatic shifts with increased errors. Because of this, we decided to fix all three geometric parameters to their photometric values for galaxies exhibiting this behavior. For comparison, we also ran DiskFit on the remaining galaxies with fixed geometry parameters.

Zoom In Zoom Out Reset image size

Figure 4 shows the tangential and radial rotation curves derived with all parameters fixed (Fixed Geom.) for NGC 5949 in blue, as well as those with all parameters allowed to vary (Fit Geom.) in green. Note that for the galaxies that converged in both cases, there is close agreement between the two methods of fitting. This can be seen in the top panels of the Appendix figures, which show the rotation curves for all galaxies in our sample.

Zoom In Zoom Out Reset image size

Two galaxies, UGC 01104 and UGC 3371, required larger bin sizes of 5 pixels, due to their sparse Hα distributions. UGC 3371 additionally displayed unphysical behavior from its innermost point, due to sparse Hα emission, leading us to omit the first bin from the rotation curve. For our two LSB galaxies, UGC 3371 and UGC 11891, there were not adequate data to robustly determine the radial components $\bar{{V}_{r}}$. We therefore did not fit for radial velocities for these two galaxies. A similar problem occurred in just the outer part of NGC 6207, where the Hα field becomes sparse. In this case, we fit for both tangential and radial velocities up to a radius of 4.2 kpc and set $\bar{{V}_{r}}=0$ for the remainder of the galaxy. Finally, our photometric PA for UGC 01104 did not agree with the kinematic PA found with DiskFit. Upon inspection, it became clear that the photometric PA is not correct. We therefore chose to allow DiskFit to fit the PA for this galaxy, initializing it with the kinematic value.

Several galaxies had a few isolated extreme outliers from the DiskFit model. We created normalized residual fields by subtracting the DiskFit model from the velocity field and dividing by the error. Bad pixels were defined as those with a normalized residual larger than 6 and were removed from the velocity field before refitting. Four galaxies required minor masking (an average of 2 pixels): NGC 1012, NGC 4396, NGC 4451, and UGC 11891. This had no material effect on the rotation curve. UGC 3371 also required masking, but for a larger region of 23 pixels. NGC 2644 had 15 pixels that met the criteria for masking; however, they are located toward the center of the galaxy, so their removal would mean losing a significant part of our data. Upon inspection, we chose to leave them in.

In any image slicing IFU, the measured velocity will depend at some level on the distribution of Hα within the slice. When using PCWI, the seeing is often smaller than the slice width, so the emission does not illuminate the entire slice. This is potentially problematic because we would measure a different velocity for a point source at the bottom of the slice versus the top. PCWI images have a pixel scale of 058 and a dispersion of 0.32 Å per pixel, so we expect the difference in Hα velocity across a 26 slice to be ΔV = 26/058 × 0.32 Å/6563 Å × c ≈ 67 km s⁻¹ for a perfect point source.

In practice, the discrepancy is not as severe as the ΔV given above. In addition to the diffuse nature of Hα reducing the magnitude of this effect, the seeing helps blur out the target over a larger portion of the slice. Furthermore, the effects are still less pronounced in a rotation curve, where velocities are averaged over annuli. Even so, the effect is potentially significant, particularly for our LSB galaxies, whose H ii regions are more point-like, so we chose to model and correct for this effect.

To correct for this effect, we first need to determine where the emission arrived within each slice. We convolve the narrowband Hα images to match the (generally worse) seeing we had with PCWI, and then we project the PCWI slice array onto the image.

In each spaxel, we compute the centroid of the Hα emission within the 26 slice width. We expect this position to correlate with a velocity shift, which we can test by examining the correlation with residuals from the DiskFit model. To determine the correct offsets, we took the DiskFit models made using the uncorrected velocity field and plotted the residuals against slice location. The top row of Figure 5 shows the original data and the steps of this process for NGC 5949. As expected, there is a clear correlation between the position in the slice and the residual.

Zoom In Zoom Out Reset image size

We use the slope of this correlation to estimate the correction in velocity as a function of position. We repeated the data reduction process using this correction to shift the velocities. The resulting corrected velocity field, DiskFit model, residual field, and residual plot can be seen in the bottom row of Figure 5. The differences in the velocity fields are subtle, but the residual correlation plot shows that the correction did remove the dependence on slice location.

This procedure was done for all of the galaxies in our sample, and in all cases using the corrected field improved the χ² value of the DiskFit output. As mentioned, the change to the velocity field is subtle, and consequently the correction generally did not make a significant change to the rotation curve, though it often decreased the errors. The maximum difference in velocity ranges from less than 1 to 15 km s⁻¹, with a median of 3.4 km s⁻¹. The typical velocity correction is thus comparable to the motions of the ISM (6 km s⁻¹) included in our DiskFit model.

The top panels in the Appendix figures show the rotation curves from the corrected fields and their errors in bold, while the rotation curves from the uncorrected fields (No Illum. Corr.) are shown with reduced opacity. For the remainder of our work, we will use the fields corrected for slit illumination.

While our intent was to select galaxies with no obvious central bar, it is still possible that some of our galaxies might have noncircular motions induced by a bar that are dynamically significant. In particular, NGC 2976 (Menéndez-Delmestre et al. 2007; Spekkens & Sellwood 2007; Valenzuela et al. 2014) and NGC 6503 (Kuzio de Naray et al. 2012) have previously been suggested to contain bars. To investigate this possibility, we used DiskFit to fit for bisymmetric noncircular motions in each galaxy (see Equation (2)) and compared the results to the rotation curves previously derived from the radial flow model. In most cases, there was no appreciable difference between the two. For six galaxies, the fit did not converge and the new rotation curves were oscillatory and had enormous errors. If a galaxy fell into either of these cases, we discarded the nonaxisymmetric model. The remaining four galaxies (NGC 949, NGC 2976, NGC 3622, and NGC 4376) showed reasonable variation and could plausibly contain a dynamically significant bar.

The top panels in the Appendix figures show the bisymmetric rotation curves in purple for the four galaxies listed. NGC 949 shows good agreement in the inner part but shows a drop in the outer parts of the curve that is not seen in the radial flow model. While the shape and magnitude of the different curves for NGC 2976 are nearly the same, the behavior in the center is different. As this region is important to our mass modeling, we chose to continue to consider the bisymmetric case. NGC 3622 and NGC 4376 both show larger velocities in the outer regions of the curve for the bisymmetric fit, and NGC 3622 also disagrees in the center.

Upon closer inspection, for NGC 3622 the residuals between the Hα velocity field and the radial flow model show a distinct central distortion. This leads us to conclude that NGC 3622 is a clear case for containing a central bar, so we will use the bisymmetric model as our fiducial case for this galaxy. We were not able to come to a clear conclusion for the remaining three galaxies. For this reason, for all four galaxies discussed above, we will consider both the radial flow and bisymmetric rotation curves in our analysis and compare results in the mass-model stage of Paper II.

DiskFit uses a bootstrap technique to measure the scatter in the model parameters about their optimal values. A bootstrap sample is generated by adding randomly drawn residuals (from the distribution with the minimum χ² value) to the optimal velocity field model (Spekkens & Sellwood 2007). The minimization is done again using the bootstrap velocities, and this process is repeated 500 times. The standard deviation of each parameter from all bootstrap realizations is then given as the uncertainty.

DiskFit provides a log containing a list of the best-fitting model parameters determined at each bootstrap realization. For the N radial bins in each galaxy, we will use these logged data to compute the N × N covariance matrix for use in the mass modeling contained in Paper II. Correlations between radial bins could easily arise in our data, due to dependence on parameters like inclination or PA. Using the covariance matrix is preferable because it takes such correlations into account.