The Enigmatic (Almost) Dark Galaxy Coma P: The Atomic Interstellar Medium

The ALFALFA blind H i extragalactic survey has now produced the largest and most diverse catalog of H i-detected systems ever assembled (Giovanelli et al. 2005; Haynes et al. 2011). The final source catalog contains more than 31,500 extragalactic sources of high significance covering over 7 dex in H i mass. The ALFALFA catalog provides unique perspectives on the population of gas-rich extragalactic sources.

Our group has been actively investigating the physical properties of optically faint galaxies as the ALFALFA catalog has matured and neared completion. These efforts have included both interferometric H i spectral line imaging as well as deep ground-based imaging of selected sources. We have targeted both systems that are candidate Local Group/Volume objects and systems that are almost certainly well outside the Local Group.

This first line of inquiry has led to the characterization of extreme regions of extragalactic parameter space. In the "Survey of H i in Extremely Low-mass Dwarfs" (hereafter "SHIELD"; Cannon et al. 2011a) we are studying those ALFALFA-detected galaxies with H i mass reservoirs smaller than ${10}^{7.2}\,{M}_{\odot }$. Using Hubble Space Telescope (HST) imaging, McQuinn et al. (2015a) show that these galaxies are characterized by fluctuating, non-deterministic, and inefficient star formation. The nature of the star formation process in these galaxies is studied in further detail in Teich et al. (2016), where it is shown that stochasticity plays an important role in governing the relationships between recent star formation and neutral hydrogen gas. A companion paper (McNichols et al. 2016) explores the rotational dynamics of these sources, highlighting that galaxies in this mass range probe the transition regime from rotational to pressure support.

This line of inquiry has also led to the ALFALFA discovery of some of the most extreme Local Volume (D ≲ 11 Mpc) galaxies known. A series of manuscripts describes Leo P, an extremely low-mass galaxy just outside the Local Group with a single H ii region and an extremely low oxygen abundance (Giovanelli et al. 2013; McQuinn et al. 2013, 2015b, 2015c; Rhode et al. 2013; Skillman et al. 2013; Bernstein-Cooper et al. 2014; Warren et al. 2015). Subsequently, the even more metal-deficient Leoncino Dwarf (AGC 198691) was discovered by ALFALFA (Hirschauer et al. 2016).

Moreover, follow-up optical and H i observations have also produced multiple detections of resolved stellar populations in "ultra-compact high-velocity clouds" (UCHVCs). These are extended H i systems that have properties such as compact, low-mass halos if they are indeed nearby (within ∼1–2 Mpc). We refer the reader to Adams et al. (2013, 2015, 2016) and Janesh et al. (2015, 2017) for details about the population of these nearby systems.

The second line of inquiry targets optically faint objects that appear to be outside the Local Volume. This includes a very small subset of ALFALFA H i detections (<1%) that have clear optical counterparts, that have velocities that indicate they are clearly extragalactic, and that are not obviously tidal debris. We have labeled these objects as "(Almost) Dark" galaxies. These enigmatic sources push the boundaries of what might typically be labeled as a "galaxy"; they thus offer unique and powerful insights into the extrema of various fundamental scaling relations.

In Cannon et al. (2015) we presented the basic selection criteria that are used to identify clearly extragalactic (i.e., with recessional velocities sufficient to cleanly separate them from Milky Way H i gas) sources that are (Almost) Dark galaxy candidates. In short, such systems lack an obvious optical counterpart in Sloan Digital Sky Survey (SDSS) imaging (which can be a result of intrinsically low surface brightness or of spatial coincidence) and are unlikely to be tidal in origin. These systems naturally have elevated ${M}_{{\rm{H}}{\rm{I}}}/{L}_{B}$ ratios (Cannon et al. 2015 note that empirically the sources usually have ${M}_{{\rm{H}}{\rm{I}}}/{L}_{B}\gtrsim 10$). As discussed in Leisman (2017) and Leisman et al. (2017), roughly 1% of the sources in the ALFALFA catalog meet these criteria and are worthy of detailed further exploration to determine their physical properties.

Characterizing these (Almost) Dark galaxies requires intensive follow-up observations, including H i synthesis observations and deep ground-based optical imaging. Previous work by our group has shown that, to date, all of the candidate (Almost) Dark objects turn out to not be truly dark galaxies. Follow-up observations reveal complex tidal systems (Leisman et al. 2016), peculiar sources with large pointing offsets between the ALFALFA and the interferometric centroids (Cannon et al. 2015; Singer et al. 2017), galaxies with low surface brightness similar to ultradiffuse galaxies (Leisman et al. 2017), and sources with other interpretive issues (see Leisman 2017).

The source that is the focus of this article lies at the intersection of these two efforts. Drawn from the population of more distant (Almost) Dark sources, recent revisions to its distance show that it is actually quite nearby—an extreme Local Volume source that still pushes typical definitions of "galaxy." This particular source attracted attention as being completely undetected in SDSS, but detected with a signal-to-noise ratio S/N > 99 in ALFALFA data products. Examining archival GALEX imaging of the field suggested an optical counterpart. Subsequent detailed H i and optical follow-up observations presented in Janowiecki et al. (2015) revealed a stellar population of extremely low surface brightness with extended H i gas. This source, AGC 229385, we hereafter will refer to as "Coma P."

Figure 1 shows the same comparison of the optical and H i morphologies of Coma P as shown in Janowiecki et al. (2015). This galaxy is part of a trio of H i line sources with similar recessional velocities that were identified by ALFALFA and confirmed in the Westerbork Synthesis Radio Telescope (WSRT¹⁴ ) H i images presented in Janowiecki et al. (2015). Very deep optical images acquired with the partially populated One Degree Imager (pODI) on the WIYN¹⁵ 3.5 m telescope at Kitt Peak Observatory¹⁶ revealed an extended stellar component of low surface brightness associated with the brightest H i source (Coma P). Stellar counterparts of the other two H i sources (AGC 229384 and AGC 229383) remain undetected, with upper limits on the surface brightness of 27.9 and 27.8 mag arcsec⁻² in the ${g}^{{\rm{{\prime} }}}$ filter, respectively. These sources have similar systemic velocities to Coma P (${V}_{{\rm{H}}{\rm{I}}}=1309\pm 1\,\mathrm{km}\,{{\rm{s}}}^{-1}$ and 1282 ± 4 $\mathrm{km}\,{{\rm{s}}}^{-1}$ for AGC 229384 and AGC 229383, respectively) and are offset by ∼97 and ∼15' from Coma P (see Figure 1). Janowiecki et al. (2015) discuss the physical properties of these sources, based on their adopted distance of 25 Mpc.

Zoom In Zoom Out Reset image size

Although the recessional velocity of Coma P (${V}_{{\rm{H}}{\rm{I}}}=1348\pm 1\,\mathrm{km}\,{{\rm{s}}}^{-1}$; Janowiecki et al. 2015) and Virgocentric flow model (Masters 2005) for this region of the sky suggest that it is located well beyond the Virgo Cluster, recent HST images of the source have demonstrated that it is located significantly closer to us. As discussed in detail in S. Brunker et al. (2017, in preparation), the red giant branch of the HST color–magnitude diagram is best fit by a distance of only 5.50 ± 0.28 Mpc. The physical properties of Coma P, which are summarized in Table 1, become increasingly extreme as a result of this change in distance. We note that the HST images do not cover either of the other H i sources that are detected in the ALFALFA and WSRT data as shown in Janowiecki et al. (2015) and Figure 1. While their proximity in both space and velocity strongly suggests association, we do not know whether these objects are at the same distance as Coma P or not.

Table 1.  Physical Characteristics of Coma P

Parameter	Value
R.A. (J2000)	12h32m103
Decl. (J2000)	+20°25'23''
Adopted distance (Mpc)	5.50 ± 0.28a
${M}_{B}$ (mag)	−10.71 ± 0.14a
${M}_{\star }$ (${M}_{\odot }$)	(1.0 ± 0.3) × 106a
${M}_{{\rm{H}}{\rm{I}}}$ (${M}_{\odot }$)	(3.48 ± 0.35) × 107b
${D}_{{\rm{H}}{\rm{I}}}$ (kpc)c	4.0 ± 0.2
${\mathrm{SFR}}_{\mathrm{FUV}}$ (${M}_{\odot }$ yr−1)	(3.1 ± 1.8) × 10−4

Notes.

^aS. Brunker et al. (2017, in preparation). ^bUsing ${S}_{{\rm{H}}{\rm{I}}}=4.87\pm 0.04$ Jy km s⁻¹ from Haynes et al. (2011). ^cMeasured at the 1 ${M}_{\odot }$ pc⁻² level.

Download table as:  ASCIITypeset image

In this work, we present new H i line synthesis observations of the neutral interstellar medium of Coma P constructed by combining new higher resolution Very Large Array (VLA¹⁷ ) H i line imaging with the WSRT H i line data set presented in Janowiecki et al. (2015). The resultant spectral grids are used to study the morphology and dynamics of the H i gas at a variety of angular resolutions. The higher resolution of this new combined data set allows detailed kinematic and morphological analysis of the H i gas in Coma P, which we present below. We also leverage the spatially resolved nature of the stellar population from the HST images in order to investigate the star formation process in this enigmatic galaxy.

This paper is organized as follows. In Section 2 we discuss the H i data handling and combination. Section 3 contains detailed discussion about the H i morphology and kinematics. In Section 4 we compare the H i properties with those of the stellar populations. The H i dynamics of Coma P are studied and modeled in Section 5. We contextualize Coma P in Section 6 and summarize our conclusions in Section 7.

2.1. H i Data Products

2.2. Combination Imaging

In order to leverage both the surface brightness sensitivity of the WSRT data and the angular resolution of the VLA data, we use both sets of calibrated visibilities to create three-dimensional data cubes of H i line emission. While the combination of single-dish and interferometric data is common (the so-called "zero-spacing correction"), the combination of data sets from different interferometric observatories is less frequently employed. In addition to the different synthesized beam shapes of the individual observatories, they also typically employ different weighting schemes for calibrated visibilities. We refer the reader to Carignan & Purton (1998) for an early discussion of techniques.

There are effectively two strategies for the combination of interferometric data sets. In the first, appropriately weighted uv visibilities are combined in the image plane (e.g., Sorgho et al. 2017). In the second, the visibilities are weighted, combined in the uv plane, and imaged by deconvolving with a beam that combines the shape of the beams for the input data sets (e.g., Bernstein-Cooper et al. 2014). We favor the second strategy for image combination, since the CASA CLEAN algorithm is specifically designed to be able to image multiple data sets on the fly.

The relative weights of the individual data sets were calculated by measuring the rms as a function of time and of baseline in line-free channels using the CASA task statwt. We then performed rigorous checking for matching H i fluxes per channel between data sets by imaging the individual, calibrated databases separately and by enforcing appropriate tapering in the uv domain to achieve similar beam sizes. As expected, the full-resolution VLA images recover less flux than the WSRT images, and both recover less than the ALFALFA flux (see below). After tapering to common resolutions we find agreement between the fluxes per channel in all data sets.

With the individual uv databases appropriately weighted, we then imaged all of them simultaneously in the CASA CLEAN algorithm. We employ the "Briggs" weighting scheme and vary the "robust" weighting parameter (which ranges from −2 to +2 for uniform to natural weighting) and the uv-domain tapering to produce final cubes at high, medium, and low angular resolutions. Table 2 summarizes basic parameters of each final data set. The final data cubes were produced using interactive clean boxes and were cleaned to the 0.5σ level (determined in line-free channels of cubes with very small numbers of clean iterations). Residual flux rescaling was employed on each resulting cube following Jörsäter & van Moorsel (1995). We present the H i data cubes at low, medium, and high resolution in channel map format in Figures 2, 3, and 4, respectively.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2.  Details of Final H i Mapping Data Products

Data Product	Beam Size	rms	Peak ${N}_{{\rm{H}}{\rm{I}}}$ a	${S}_{{\rm{H}}{\rm{I}}}$
		(mJy beam−1)	(cm−2)	(Jy km s−1)
VLA data only	690 × 496	0.77	17.25 × 1020	3.02 ± 0.30
Combined data, High resolution	75 × 75	0.69	12.8 × 1020	2.92 ± 0.29
Combined data, Medium resolution	125 × 125	0.56	9.3 × 1020	2.96 ± 0.30
Combined data, Low resolution	17'' × 17''	0.66	8.5 × 1020	3.56 ± 0.36
WSRT data only	625 × 178	0.58	5.89 × 1020	4.43 ± 0.44b

Notes.

^aH i column density as defined assumes a uniform source that fills the beam, so it should always increase with increasing angular resolution if there is structure. ^bJanowiecki et al. (2015).

Download table as:  ASCIITypeset image

Each resulting data cube was threshold-blanked at the 2σ level and then closely examined by hand for spectral and spatial continuity of features. The global H i line emission in these blanked cubes is plotted in Figure 5. The systemic velocity of Coma P in each cube is assigned to the intensity-weighted mean. The resulting blanked data cubes were then collapsed to create two-dimensional images of H i mass surface density, H i intensity-weighted velocity, and H i intensity-weighted velocity dispersion. Note by examining Table 2 that the integrated H i flux from Coma P increases as the angular resolution decreases, and that the products with the lowest angular resolution (WSRT data only, highest surface brightness sensitivity) recover ${S}_{{\rm{H}}{\rm{I}}}$ = 4.43  ±  0.44 Jy km s⁻¹. This agrees with the ALFALFA H i flux integral (${S}_{{\rm{H}}{\rm{I}}}=4.87\pm 0.04$ Jy km s⁻¹) within errors.

Zoom In Zoom Out Reset image size

Our combined VLA+WSRT images reveal the properties of the H i gas in Coma P in unprecedented detail. Our low-resolution images are in excellent agreement with those presented in Janowiecki et al. (2015). We refer the reader to that work for a discussion of the gas with extended low surface brightness in Coma P, as well as for details about the other H i line sources that are in close angular and velocity proximity.

3.1. H i Morphology

In Figure 6 we show two-dimensional representations of the H i column density (calibrated in units of 10²⁰ cm⁻²), the intensity-weighted H i velocity, and the intensity-weighted H i velocity dispersion. Each of these images is shown at low (17'' beam), medium (125 beam), and high (75 beam) angular resolutions. The corresponding elements of physical resolution are 453 pc, 333 pc, and 200 pc, respectively, assuming the distance of 5.5 Mpc from S. Brunker et al. (2017, in preparation).

Zoom In Zoom Out Reset image size

In panels (a), (d), and (g) of Figure 6, the H i column density distributions are shown in grayscale and with uniformly spaced overlying contours. At ∼450 pc resolution the H i is smoothly distributed along an axis running roughly southeast to northwest, while H i gas of lower surface brightness extends beyond the outermost contour shown. The overall H i morphology at low angular resolution agrees well with the WSRT-only images as presented in Janowiecki et al. (2015).

At low and medium resolutions, the H i morphology is highlighted by moderate H i mass surface densities throughout most of the galaxy. At low (Figure 6(a)) and medium (Figure 6(d)) angular resolutions, H i columns exceed ${N}_{{\rm{H}}{\rm{I}}}\gtrsim 5\times {10}^{20}$ cm⁻² ($4.0\,{M}_{\odot }$ pc⁻²). Two maxima in H i column density are located on either end of the H i distribution. The densest H i gas (${N}_{{\rm{H}}{\rm{I}}}\gtrsim 9\times {10}^{20}$ cm⁻² or $7.2\,{M}_{\odot }$ pc⁻²) is located in the southeastern end of the H i distribution, while the secondary H i maximum peaks at only slightly lower H i column densities. In between these maxima is the morphological centroid of the H i gas, where the H i column densities are slightly lower.

At 200 pc resolution (Figure 6(g)) we resolve only the highest peaks in mass surface density in the H i gas. Interestingly, the H i maxima seen in the lower resolution images resolve into multiple individual H i peaks. At a few locations the H i column density exceeds the "canonical" threshold of ${N}_{{\rm{H}}{\rm{I}}}={10}^{21}$ cm⁻² that was empirically associated with ongoing star formation in dwarf irregular galaxies (Skillman 1987). The southeastern H i maximum seen at low and medium resolutions achieves the highest H i mass surface density in the entire galaxy (${N}_{{\rm{H}}{\rm{I}}}=1.28\times {10}^{21}$ cm⁻² or $10.3\,{M}_{\odot }$ pc⁻²).

As we discuss in detail below, the H i gas in Coma P is more complicated than a simple thin disk such as those found in many other dwarf galaxies. With this caveat in mind, we define the morphological center of Coma P by examining the medium-resolution H i data. Using the 3 × 10²⁰ cm⁻² contour as a guide, we find that an ellipse with an axial ratio of 2.1:1, oriented along a major axis rotated 155° east of north, provides a suitable representation of the H i morphology at this surface density level. The H i centroid derived using this approach is (α, δ) = (12^h32^m103, +20°25' 23''). The corresponding inclination is 63°, which agrees with the i = 63° ± 4° derived by Janowiecki et al. (2015).

The H i distribution is resolved along both the major and minor axes. We thus use the parameters derived above to create the radial profiles of the H i gas in Coma P that we present in Figure 7. These profiles (created using the ELLINT task in the GIPSY¹⁹ software package; van der Hulst et al. 1992) are explicitly corrected for beam smearing using the Gaussian approximation discussed in Wang et al. (2016). Specifically, the corrected and uncorrected H i sizes (${D}_{{\rm{H}}{\rm{I}}}$ and ${D}_{{\rm{H}}{\rm{I}},0}$, respectively) are related to the lengths of the major and minor axes of the beam (${B}_{\mathrm{maj}}$ and ${B}_{\min }$, respectively) via ${D}_{{\rm{H}}{\rm{I}}}=\sqrt{{D}_{{\rm{H}}\,{\rm{I}},0}^{2}-{B}_{\mathrm{maj}}\cdot {B}_{\min }}$. Beam smearing is negligible when the observed H i size is much larger than the angular size of the beam; this is the case for our H i observations of Coma P.

Zoom In Zoom Out Reset image size

The radial profiles allow us to quantify both the H i size and the H i scale length of Coma P. The average (radially integrated) H i mass surface density in the inner region is nearly flat within errors at a value of $\sim 4.25\,{M}_{\odot }$ pc⁻² (at 125 resolution). There is a small increase in mass surface density at ∼0.75 kpc in both profiles, corresponding to where the radial integration passes through both the primary and the secondary H i maxima. Moving outward the profile falls off smoothly. The horizontal dotted lines in Figure 7 demonstrate that the H i diameter is 4.0 ± 0.2 kpc measured at the $1\,{M}_{\odot }$ pc⁻² level. The H i scale length of Coma P is measured by fitting an exponential function to the six data points at radii larger than 1 kpc in the low-resolution radial profile (Figure 7(b)). We find a scale length of 0.8 ± 0.1 kpc, where the uncertainty is estimated by including and then excluding the next inner radial profile point in the exponential fit.

3.2. H i Kinematics

To understand the relationship of the H i gas to the underlying stellar populations in Coma P, in Figure 8 we overlay the same sets of H i column density contours as shown in panels (a), (d), and (g) of Figure 6 on a color HST image (created using the F606W filter as blue, the F814W filter as red, and a linear average of the two filters as green; see S. Brunker et al. 2017, in preparation). Figure 8(a) shows the same HST color image as in S. Brunker et al. (2017, in preparation) with no H i contours so as to allow detailed inspection of the stellar population and to identify foreground stars and background galaxies. As described in detail in Janowiecki et al. (2015) and S. Brunker et al. (2017, in preparation), the stellar population of Coma P has a remarkably low surface brightness (peak surface brightnesses in the ${g}^{{\rm{{\prime} }}}$, ${r}^{{\rm{{\prime} }}}$, and i' bands of 26.4 ± 0.1 mag arcsec⁻², 26.5 ± 0.1 mag arcsec⁻², and 26.1 ± 0.1 mag arcsec⁻², respectively). With a morphology similar to that of the H i images shown in Figure 6, the stellar population arcs from southeast to northwest. Nearly all of the compact and luminous sources seen in this field are foreground stars or background galaxies. The individual stars in Coma P are faint but spatially resolved, resulting in the diffuse distribution of stars that is contained within the H i contours. A few stellar clusters are apparent, and these are discussed in further detail below.

Zoom In Zoom Out Reset image size

In panels (b)–(d) of Figure 8, the H i contours as shown in Figure 6 are color-coded for ease of interpretation. Comparing the low- and medium-resolution H i contours to the HST images, we find excellent spatial agreement between the H i gas of high column density and the locations of the faint but resolved stars in Coma P. Effectively all of the stellar population of Coma P is located within H i gas with ${N}_{{\rm{H}}{\rm{I}}}\gtrsim 5\times {10}^{20}$ cm⁻² ($4.0\,{M}_{\odot }$ pc⁻²). The highest H i mass surface densities in the southeastern region of the galaxy are cospatial with what appear to be individual stellar clusters or extended star formation regions.

To further investigate the relationship between the H i gas and the resolved stars in Coma P, Figure 9 shows the locations of individual red and blue stars on the HST color image. Panel (a) shows the color–magnitude diagram as presented in S. Brunker et al. (2017, in preparation). Therein, blue and red stars are color-coded. Blue stars are selected as having F814W magnitudes brighter than 27.0 and F606W−F814W colors bluer than 0.4, while red stars are the individual red giants used for the distance determination as discussed in S. Brunker et al. (2017, in preparation). These red and blue stars are then plotted on the color HST image in panel (b) as cyan and red asterisks. H i column density contours from the medium-resolution (125 beam) image are shown at levels of (7, 8, 9) × 10²⁰ cm⁻². As expected for an old stellar population that is well mixed, the red stars are uniformly distributed across the galaxy. In contrast, the blue stars are decidedly concentrated within the regions of highest H i mass surface density.

Zoom In Zoom Out Reset image size

As discussed in detail in Janowiecki et al. (2015), the star formation properties of Coma P are extreme. The source is not detected in deep ground-based Hα imaging, implying that the current rate of massive star formation is zero. Curiously, Coma P is detected by GALEX in both the near-UV (NUV) and far-UV (FUV) bands (although at low S/N). We compare the low-resolution H i image with the GALEX FUV and NUV images in Figure 10. The ultraviolet emission detected by GALEX is cospatial with the H i gas detected in Coma P. By comparison with the HST images in Figure 8 it is evident that the flux from the resolved stellar populations is dominated by stars that are capable of producing UV photons but not of photoionizing H ii regions (S. Brunker et al. 2017, in preparation). Note that the southern maximum in H i column density is cospatial with the maxima in both GALEX images. This is the same location identified above in the HST images as harboring stellar clusters. We note, however, that there is also a significant background spiral galaxy nearby that makes this interpretation difficult (see the highest H i contour in Figure 8(c) and in both panels of Figure 10).

Zoom In Zoom Out Reset image size

We estimate the FUV-based star formation rate using the calibration derived from resolved stars in McQuinn et al. (2015d). The resulting global far-ultraviolet star formation rate is very low: ${\mathrm{SFR}}_{\mathrm{FUV}}=(3.1\pm 1.8)\times {10}^{-4}\,{M}_{\odot }$ yr⁻¹. This is equal to the lowest far-ultraviolet star formation rate among the galaxies in the SHIELD sample (Teich et al. 2016). However, the FUV flux in Coma P is extended over a much larger area than in the (physically smaller) SHIELD galaxies. A spatially resolved analysis of the FUV star formation rate density versus the H i mass surface density (e.g., as presented in Teich et al. 2016) is not possible at the low S/N of the GALEX FUV images. Deeper UV imaging of Coma P would be especially useful for studying the star formation on timescales of 100–200 Myr.

As discussed in Section 3.2, the H i kinematics of Coma P are not amenable to standard modeling techniques. The two- and three-dimensional data products (see Figures 2, 3, 4, 6) demonstrate this clearly. In an effort to understand the complex neutral gas dynamics of Coma P, we thus perform two types of analysis on the three-dimensional data cubes.

5.1. Spatially Resolved Position–Velocity Analysis

In the first line of investigation we apply the spatially resolved position–velocity analysis developed in Cannon et al. (2011b) and employed in McNichols et al. (2016) and Cannon et al. (2016). This technique allows us to estimate the projected rotation velocity of a source when its kinematics are not well described by a simple tilted ring. The approach relies on identifying a projected H i major axis, which is typically that axis along which the largest velocity gradient is present. In some cases this H i major axis aligns with the major axis of the stellar component. However, this is not necessarily the case for all dwarf galaxies, especially for those systems with extremely unusual gas kinematics (see application of the approach to the post-starburst galaxy DDO 165 in Cannon et al. 2011b). The H i minor axis is assumed to be offset from the major axis by 90°.

The spatially resolved position analysis then examines a collection of position–velocity slices for the major and minor axes. Each slice is the width of, and is offset from the neighboring slice by, the H i beam width. The series of major and minor axis slices allow the identification of the largest projected rotation velocity, both along the major axis and also as a change in velocity of the H i centroids of the minor axis slices. Further, this approach allows us to isolate kinematically distinct features (e.g., Cannon et al. 2011b) in the gas.

To balance angular resolution with sensitivity, in this analysis we examine the medium-resolution (125 beam) data cube (see Figure 3). We identify the H i major axis interactively by examining position–velocity slices at a range of position angles. The resulting major axis is oriented at a position angle of 335° measured east of north. This angle is in agreement with the direction of the maximum projected velocity gradient in the moment one images at both low and medium resolutions (see Figure 6), and it is also in rough agreement with the position angle of the stellar component (see Figure 1). The H i minor axis is 90° away (at 65°), again rotated east through north. The H i dynamical center is assumed to be cospatial with the adopted morphological center (see discussion in Section 3.1): (α, δ) = (12^h32^m103, +20°25'23''). Note that since the slices are the width of the H i beam, the absolute position of the dynamical center is not critical for interpretation.

Based on the 125 (333 pc) beam, we then symmetrically slice (i.e., use the same number of slices at positive and negative angular offsets from the major and minor axes) the data cube along both the major and the minor axes as shown in Figure 11. Five major axis slices are required to cover the H i distribution of Coma P, while 13 minor axis slices are required. Such an axial ratio is expected based on the inclination derived above. On both the major and the minor axes, the central slice is numbered "0." For the major axis slices, those with negative numbers are located east of the central slice and those with positive numbers are located west of it. Likewise, for the minor axis slices, those with negative numbers are located north of the central slice and those with positive numbers are located south of it.

Zoom In Zoom Out Reset image size

The resulting spatially resolved position–velocity slices for the major and minor axes are presented in Figures 12 and 13, respectively. In each figure the axes and color scalings are identical. Contours are overlaid at the 3σ, 6σ, 9σ, and 12σ levels in each panel. These figures contain a wealth of information about the kinematics of the H i gas in Coma P.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Beginning with the major axis slices shown in Figure 12, all slices show two opposing velocity gradients: a crescent-shaped structure in the plane of velocity versus offset, the shape and dynamic range of which are the same as shown in Janowiecki et al. (2015). This signature of two opposing velocity gradients is most easily interpreted as two separate velocity components within the H i distribution of Coma P. Moving from the central slice to either positive or negative offset reveals that the H i velocity structure on the major axis does not change as a function of offset from the center. This confirms the identification of the major axis of one of the H i components. Further, it demonstrates that the turnover in velocity as a function of offset, likely caused by the presence of a second H i kinematic component, extends throughout the entire H i distribution of Coma P.

The minor axis slices confirm the interpretation of two kinematic components in the H i gas. The H i velocity structure along the minor axis does indeed show a moving velocity centroid from one side of the galaxy to the center (e.g., examine the highest S/N contours in the series of slices with negative numbers). However, importantly, around the central slice this signature turns over in velocity space and then continues in the opposite direction on the other side of the galaxy (examine the highest S/N contours in the series of slices with positive numbers). Such a turnover in velocity space is very unusual. A galaxy with normal H i kinematics (that are well described by a simple tilted ring model) would show the signature of changing centroid position across all of the minor axis cuts. Note that even if the inclination were in fact very low (unlikely for Coma P based on the optical and H i morphologies) and the signature of projected rotation very subtle, it would still not result in a turnover in velocity space as a function of position.

Using the spatially resolved position–velocity analysis as a guide, we can coarsely separate the two H i gas components in the observed low- and medium-resolution data cubes. We stress that this separation is not clean, either spatially or spectrally. H i gas connects these two kinematic components in the data cubes and in the position–velocity slices through them. With this caveat in mind, we estimate that the relative H i mass ratio of the two H i components is in the range of (0.6:1) to (0.8:1); varying the masking and threshold levels used to identify the two components changes the resulting mass ratio. While the absolute value of the mass ratio is uncertain, the masses of the two components are of the same order.

It is difficult to interpret the velocity gradients seen in the position–velocity slices of Figures 12 and 13 as coherent rotation. If we take the projected velocity gradients at face value, then we can make coarse estimates of the dynamical masses of the two kinematic components. We use the turnover in projected velocity along the central major axis slice as the boundary between the two components. The more massive component (extending to positive angular offsets in Figure 12) has a maximum observed projected velocity gradient of Δv = 35  ±  10 $\mathrm{km}\,{{\rm{s}}}^{-1}$ that spans ∼60''. Likewise, the less massive component also spans an angular offset of ∼60'' with a maximum observed projected velocity gradient of v = 30 ± 10 $\mathrm{km}\,{{\rm{s}}}^{-1}$. If we interpret each of these maximum projected velocity gradients as the signature of rotation of each kinematic component, then the inclination-corrected rotational velocities of the more and less massive components are ∼20 ± 10 $\mathrm{km}\,{{\rm{s}}}^{-1}$ and ∼17 ± 10 $\mathrm{km}\,{{\rm{s}}}^{-1}$, respectively; each subtends ∼800 pc (30''). We stress that this interpretation is highly uncertain, since the dynamics of Coma P are certainly a mixture of rotation and non-circular motions.

Using the estimates of angular offset and rotational velocity above, a simple estimate of the dynamical masses (${M}_{{\rm{dyn}}}={\textstyle \tfrac{{V}^{2}R}{G}}$) of the two kinematic components of Coma P results in ${M}_{\mathrm{dyn}}\,=7.2\times {10}^{7}\,{M}_{\odot }$ and $5.3\times {10}^{7}\,{M}_{\odot }$. We conservatively estimate a 50% uncertainty to highlight the aforementioned difficulties in determining the physical properties of Coma P. Correcting the H i gas mass by a factor of 1.35 for helium, metals, and molecules, and assuming ${M}_{\star }=(1.0\pm 0.3)\times {10}^{6}\,{M}_{\odot }$ (see S. Brunker et al. 2017, in preparation), we find a total baryonic mass ${M}_{{\rm{bary}}}=4.8\times {10}^{7}\,{M}_{\odot }$. Summing the dynamical masses of the two components (${M}_{\mathrm{dyn}}\simeq 1.2\times {10}^{8}\,{M}_{\odot }$) within the extent of the H i shows that Coma P is dominated by dark matter at the level of roughly 2.5 to one.

5.2. Three-dimensional Modeling of the H i Distribution

The second line of investigation attempts to fit a three-dimensional model to the H i data cube using the "Tilted Ring Fitting Code" ("TiRiFiC") software package (Józsa et al. 2007). By working with the full three-dimensional data cube, TiRiFiC eliminates the loss of information when using only the velocity field for tilted ring fitting. Further, TiRiFiC allows for modeling of multiple H i disks and for inhomogeneities in the surface brightness distribution. In situations where complicated gas kinematics are present (as in Coma P), this allows a direct comparison of the resulting model with the observed data cube.

The TiRiFiC package was used in interactive mode to fit the medium-resolution (125 beam) data cube. The resulting best-fit models contain two H i components: a more massive component (located in the southeastern region of the galaxy) and a less massive segment that is undergoing either infall or outflow and is moving with respect to the center of the galaxy. No simple model with a single inclined disk was able to accurately describe the H i kinematics of Coma P.

Figures 14 and 15 present this best-fit model as contours overlaid on the same low-resolution and medium-resolution channel maps as shown in Figures 2 and 3, respectively. The model reproduces the complicated three-dimensional gas distribution and kinematics very well. Moving from low to high velocities, one can easily differentiate the two H i components, both spatially and spectrally. The opposing velocity gradients of the two components are also clearly visible. The two model kinematic components have a relative mass ratio of 0.64:1, with the more massive component in the southeastern region of the disk. This compares favorably with the H i and dynamical mass ratios of the two components estimated in Section 5.1.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The interpretation of the best-fit model from TiRiFiC is very similar to the interpretation of the spatially resolved position–velocity slices discussed above. Specifically, the best-fit model contains two kinematically distinct H i components. The kinematic discontinuity in the central region of the galaxy (highest velocity gas in Figures 14 and 15) is the result of the presence of two H i components of differing masses. Whether this event is best described by two low-mass galaxy disks or by the infall of H i gas into Coma P is not certain. Given the current and recent quiescence of the galaxy in terms of star formation rate (see discussion above and in S. Brunker et al. 2017, in preparation), it is unlikely that the kinematics are described as an outflow driven by star formation. The extended gas of low surface brightness would then be the result of this ongoing interaction or infall event.

Various observed properties of Coma P make it an extremely unusual dwarf galaxy when compared to other objects in the Local Volume. The proximity of the galaxy is unusual given its large recessional velocity. The source is currently quiescent in terms of massive star formation but yet has significant UV flux indicative of star formation on timescales of 100–200 Myr. The optical surface brightness is extremely low, and the H i properties are unusual for a low-mass halo. We now attempt to contextualize Coma P with respect to other Local Volume galaxies.

Comparing Coma P to the low-mass galaxies in the SHIELD sample (Cannon et al. 2011a; McNichols et al. 2016; Teich et al. 2016), Coma P is slightly more H i-massive than all but one of these galaxies. Given this H i mass, Coma P is large in physical size. As discussed in Section 3.1, the H i diameter is measured to be 4.0  ±  0.2 kpc at the $1\,{M}_{\odot }$ pc⁻² level. This allows us to place Coma P on the ${M}_{{\rm{H}}{\rm{I}}}\mbox{--}{D}_{{\rm{H}}{\rm{I}}}$ relation—the empirical relationship between H i size and total H i mass (see first characterization in Broeils & Rhee 1997). While the origins of this relation are likely multi-faceted, Wang et al. (2016) use data for hundreds of nearby galaxies to show that the relationship is linear over four orders of magnitude in H i mass ($7\lesssim \mathrm{log}({M}_{{\rm{H}}{\rm{I}}}/{M}_{\odot })\lesssim 11$). The relation appears to continue to lower masses, although the small number of sources makes this extension uncertain at the present time.

In Figure 16 we show the ${M}_{{\rm{H}}{\rm{I}}}\mbox{--}{D}_{{\rm{H}}{\rm{I}}}$ relation including Coma P and a subset of the data from the sample of Wang et al. (2016). Coma P is just consistent with the 3σ fit derived by Wang et al. (2016) (note that the parameterizations of the ${M}_{{\rm{H}}{\rm{I}}}\mbox{--}{D}_{{\rm{H}}{\rm{I}}}$ relation by Broeils & Rhee (1997) and Wang et al. (2016) are effectively identical). The galaxy can be considered to have a large H i distribution given its H i mass, a smaller total H i mass given its physical size, or an artificially extended H i size due to the presence of two H i components.

Zoom In Zoom Out Reset image size

It is important to stress that the ${M}_{{\rm{H}}{\rm{I}}}\mbox{--}{D}_{{\rm{H}}{\rm{I}}}$ relation is meaningful for isolated galaxy disks where the physical size of the H i component is directly related to the gravitational potential of the parent halo. A complicated system such as Coma P, where the H i kinematics favor a model with multiple H i components, may not be appropriately interpreted in the context of the M${}_{{\rm{H}}{\rm{I}}}\mbox{--}{D}_{{\rm{H}}{\rm{I}}}$ relation. With this caveat in mind, taken at face value, the H i is considerably extended in Coma P compared to expectations based on the total H i mass of the system.

In addition to the large H i size of Coma P, its total H i mass is large compared to its stellar mass. Janowiecki et al. (2015) estimate a very large ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$ ratio using stellar masses derived from standard mass-to-light scaling relations. While ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$ is independent of distance, S. Brunker et al. (2017, in preparation) suggest a significantly revised stellar mass from HST imaging of ${M}_{\star }=(1.0\pm 0.3)\times {10}^{6}\,{M}_{\odot }$. Regardless of the exact value of ${M}_{\star }$, Coma P is one of the most H i-rich galaxies known relative to its stellar mass, especially among sources with accurately measured ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$. The estimated change in stellar mass moves Coma P within the scatter of the sources with stellar masses near 10⁶ ${M}_{\odot }$ in Huang et al. (2012) and Janowiecki et al. (2015). However, we note that this is not a completely fair comparison, since these ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$ ratios are less well measured (deep optical observations tend to recover more flux, thus reducing the estimated ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$ ratio). The revised ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$ value also places Coma P slightly above the extrapolation of the best-fit trend line for galaxies of low surface brightness found in McGaugh et al. (2017).

Additionally, we note that, especially at low mass, other individual systems have been identified with higher ${M}_{{\rm{H}}{\rm{I}}}/{M}_{\star }$ ratios than the one that we find in Coma P. For example, Matsuoka et al. (2012) find that the northeastern component of the system HI 1225+01 is more extreme than Coma P, although Coma P is of lower surface brightness. Similarly, some of the recently discovered resolved stellar populations of the UCHVC population discovered by ALFALFA (Adams et al. 2013, 2016) have ${M}_{{\rm{H}}{\rm{I}}}$/M${}_{\star }$ ratios in excess of 100 (Janesh et al. 2015, 2017).

Coma P is shown on the baryonic Tully–Fisher relation (Tully & Fisher 1977) in Figure 17. The derivations of rotational velocity and baryonic mass use the same techniques for Coma P (Section 5.1) and the low-mass comparison samples (Bernstein-Cooper et al. 2014; McNichols et al. 2016). However, importantly, we stress that the dynamics of Coma P favor a multiple-disk scenario, and therefore it can only be shown on this plot with the estimate of rotation velocity from one of the two kinematic components. We use the inclination-corrected velocity of the more massive component (v $\sim \,20\pm 10\,\mathrm{km}\,{{\rm{s}}}^{-1}$) and the total baryonic mass of the entire Coma P system (which is well represented by the gas mass because the stellar mass is so low). The position of Coma P is highly uncertain and is meant to be demonstrative only.

Zoom In Zoom Out Reset image size

Compared to the SHIELD galaxies and other systems studied in McNichols et al. (2016), Coma P is slightly offset from the best-fit relation in the plane of ${M}_{{\rm{bary}}}$ versus ${V}_{{\rm{c}}}$, although in agreement within the errors on the measurement of the rotational velocity. However, we stress that the error bars on the estimates of the rotational velocity of Coma P make its placement representative only. Interestingly, while Coma P lies just slightly above the empirical threshold differentiating rotationally-supported and pressure-supported systems that was identified by McNichols et al. (2016), its interpretation is nonetheless very challenging due to its complicated H i gas dynamics.

We have presented new H i spectral line imaging of the curious Local Volume dwarf irregular galaxy Coma P (AGC 229385). This source is unique in many respects. It has an extremely low optical surface brightness, complicated H i gas dynamics, and very intriguing patterns of recent star formation. An understanding of the physical properties of this galaxy is important for galaxy evolution in the low-mass regime.

This galaxy is one of three H i line sources detected by ALFALFA at recessional velocities between 1282 and 1343 km s⁻¹ and in close angular proximity to each other. Remarkably, the new HST imaging presented in S. Brunker et al. (2017, in preparation) has shown that Coma P is much closer than the original Virgocentric flow model of Masters (2005) predicted (D ≃ 25 Mpc). At D = 5.50 ± 0.28 Mpc, the total H i mass is ${M}_{{\rm{H}}{\rm{I}}}=(3.48\pm 0.35)\times {10}^{7}\,{M}_{\odot }$, reshaping and complicating our interpretive framework for this already enigmatic source. At this time the physical association of Coma P with the other H i line sources (AGC 229384 and AGC 229383) remains unconfirmed, although based on their angular and velocity proximity it is likely that they are all located at roughly the same distance.

Using a combination of deep WSRT and VLA observations, we study the H i properties of Coma P at physical resolutions of ∼450, 333, and 200 pc. The H i morphology is characterized by two maxima in mass surface density near either end of an inclined H i distribution. The H i disk of Coma P is physically large (diameter 4.0 ± 0.2 kpc measured at the $1\,{M}_{\odot }$ pc⁻² level) compared to galaxies of comparable H i mass. Its globally averaged H isurface density is slightly lower than expected based on the known scaling between H i mass and H i diameter. At the same time, it has an exceptionally high ratio of H i mass (${M}_{{\rm{H}}{\rm{I}}}=(3.48\pm 0.35)\times {10}^{7}\,{M}_{\odot }$) to stellar mass (${M}_{\star }=(1.0\,\pm 0.3)\times {10}^{6}\,{M}_{\odot }$).

An intrinsically faint dwarf, Coma P's extreme optical surface brightness renders it completely undetectable in SDSS. Coma P is currently not forming massive stars (as traced by nebular Hα emission), even though the H i mass surface densities exceed 10²¹ cm⁻² in localized regions. The galaxy has formed stars during the most recent few hundred million years. Far-ultraviolet emission as traced by GALEX images shows good spatial agreement with the high column density of H igas in the inner regions of the galaxy.

The H i kinematics of Coma P are extremely complex and are not well described by a simple tilted ring model. There are highly irregular gas motions and macroscopic velocity gradients at multiple position angles. Using both three-dimensional modeling and spatially resolved position–velocity analysis, we find that the H i properties are best described either by two colliding H i disks or by a significant infall event. A coarse separation of these H i components shows that they are of comparable mass.

Using the most obvious signatures of rotation in the three-dimensional data, we estimate the total dynamical mass of the Coma P system to be (1.2 ± 0.6) × 10⁸ ${M}_{\odot }$. It is important to note that this is the sum of the estimates of the dynamical masses of the two kinematically distinct H i components. Coma P is not unusually dominated by dark matter (${M}_{{\rm{dyn}}}/{M}_{{\rm{bary}}}\simeq 2.5$) compared to other Local Volume dwarf galaxies with similar dynamical masses. In the context of the fundamental ${M}_{{\rm{H}}{\rm{I}}}\mbox{--}{D}_{{\rm{H}}{\rm{I}}}$ and baryonic Tully–Fisher relations, Coma P is slightly offset from, but in general agreement with, the trends seen for other Local Volume galaxies. We stress that the presence of two H i components complicates the interpretation of these fundamental scaling relations.

The origins of the exceptional physical properties of Coma P may result from one or more infall and/or interaction events. There are H i clouds with similar velocities to that of Coma P in close angular proximity. While we do not detect stellar components of these clouds in deep ground-based images (Janowiecki et al. 2015), if these are in fact physically associated with Coma P then they may suggest a past significant merger event. We refer the reader to the discussions in S. Brunker et al. (2017, in preparation) of the local environment and large-scale structure surrounding Coma P.

Assuming that the irregular kinematics of the H i gas in Coma P do result from a merger of two low-mass disks, this offers a unique opportunity to compare observations with predictions from simulations. Hydrodynamical simulations with realistic feedback prescriptions are just now able to model halos in the mass range of Coma P (see, for example, the discussion of the "Feedback in Realistic Environments" simulations that are discussed in S. Brunker et al. 2017, in preparation). The unusual gas kinematics in Coma P are qualitatively similar to those shown in the simulations of gas-rich dwarf galaxies and dark companions by Starkenburg et al. (2016). Further, those simulations naturally predict off-center star formation in the resulting merged system—similar to the properties of H i gas and star formation seen in Coma P.

Regardless of their origins, the physical characteristics of Coma P make it an exceptional galaxy. Its physical properties are among the most extreme of any known source in the Local Volume. Coma P will serve as a critical benchmark for our understanding of low-mass galaxy evolution.

J.M.C. and C.B. acknowledge support from past NSF grant AST-1211683. The authors are grateful to Macalester College for supporting this project. The ALFALFA team at Cornell is supported by NSF grants AST-0607007 and AST-1107390 to R.G. and M.P.H. and by grants from the Brinson Foundation. S.J. acknowledges support from the Australian Research Council's Discovery Project funding scheme (DP150101734). K.L.R. and W.F.J. are supported by NSF grant AST-1615483. M.G.J. acknowledges support from the grant AYA2015-65973-C3-1-R (MINECO/FEDER, UE).

Facilities: VLA - Very Large Array, HST - , GALEX - , WIYN 3.5 m. -