The Green Bank Ammonia Survey (GAS): First Results of NH3 mapping the Gould Belt

TL;DR

The Green Bank Ammonia Survey (GAS) provides detailed NH3 emission maps of Gould Belt star-forming regions, revealing correlations with dust, gas kinematics, and temperature variations linked to star formation activity.

Contribution

This study presents the first data release of GAS, mapping NH3 in multiple Gould Belt regions and analyzing its relation to dust emission, gas dynamics, and temperature, offering new insights into star formation environments.

Findings

NH3 emission closely matches dust emission in most regions.

NH3 is detected in over 60% of dense lines-of-sight, except B18.

Gas temperature and velocity dispersion increase with star formation activity.

Abstract

We present an overview of the first data release (DR1) and first-look science from the Green Bank Ammonia Survey (GAS). GAS is a Large Program at the Green Bank Telescope to map all Gould Belt star-forming regions with $A_V \gtrsim 7$ mag visible from the northern hemisphere in emission from NH$_3$ and other key molecular tracers. This first release includes the data for four regions in Gould Belt clouds: B18 in Taurus, NGC 1333 in Perseus, L1688 in Ophiuchus, and Orion A North in Orion. We compare the NH$_3$ emission to dust continuum emission from Herschel, and find that the two tracers correspond closely. NH$_3$ is present in over 60\% of lines-of-sight with $A_V \gtrsim 7$ mag in three of the four DR1 regions, in agreement with expectations from previous observations. The sole exception is B18, where NH$_3$ is detected toward ~ 40\% of lines-of-sight with $A_V \gtrsim 7$ mag. Moreover, we find that the NH$_3$ emission is generally extended beyond the typical 0.1 pc length scales of dense cores. We produce maps of the gas kinematics, temperature, and NH$_3$ column densities through forward modeling of the hyperfine structure of the NH$_3$ (1,1) and (2,2) lines. We show that the NH$_3$ velocity dispersion, ${\sigma}_v$, and gas kinetic temperature, $T_K$, vary systematically between the regions included in this release, with an increase in both the mean value and spread of ${\sigma}_v$ and $T_K$ with increasing star formation activity. The data presented in this paper are publicly available.

Figures

Tables

Table 1. Table 1: GAS DR1 region data

Cloud	Region	Distance	Number of
		(pc)	footprintsaaEach footprint is 10′$\times$ 10′.
Taurus	B18	$135\pm 20$bbSchlafly et al. (2014)	11
Ophiuchus	L1688	$137.3\pm 6$ccOrtiz-León et al. (2017)	10
Perseus	NGC 1333	$260\pm 26$bbSchlafly et al. (2014)	8
Orion	Orion A-North	$414\pm 7$ddMenten et al. (2007)	12

Table 2. Table 2: GAS VEGAS spectral line setup

Species	Transition	Rest frequency	Reference
		(MHz)
NH3	$(1,1)$	23694.4955
NH3	$(2,2)$	23722.6336
NH3	$(3,3)$	23870.1296
HC5N	$9-8$	23963.9010	LOVAS
HC7N	$21-20$	23687.8974(6)	CDMS
HC7N	$22-21$	24815.8772(6)	CDMS
C2SaaObserved in central KFPA beam only	$2_1-1_0$	22344.030(1)	CDMS

Table 3. Table 3: GAS relative beam calibration

Pol	Beams
	1	2	3	4	5	6	7
LL	0.894 (0.071)	0.865 (0.101)	0.847 (0.049)	0.772 (0.074)	0.925 (0.062)	0.677 (0.081)	0.851 (0.080)
RR	0.845 (0.060)	0.846 (0.072)	0.816 (0.074)	0.745 (0.077)	0.764 (0.057)	0.511 (0.046)	0.970 (0.081)

Equations

$$W(r)=\frac{J_1(\pi r/a)}{(\pi r/a)}{e}^{-r^2/b^2}\frac{1}{T_{\mathrm{sys}}^2}$$

$$T_{\mathrm{MB}}(\nu )=\left[J(T_{\mathrm{ex}})-J(T_{\mathrm{bg}})\right]\left[1-e^{-\tau (\nu )}\right].$$

$$N_{(2,2)}=N_{(1,1)}\frac{5}{3}\exp \left(-\frac{\Delta E}{k{T}_{\mathrm{rot}}}\right)$$

$$T_{\mathrm{rot}}=T_K\left\{1+\frac{T_K}{T_0}\ln \left[1+\frac{3}{5}\exp \left(-\frac{15.7\text{leavevmode}\mathrm{K}}{T_K}\right)\right]\right\}.$$

$$Z_{\mathrm{para}}$$

$$\frac{N_u}{N_l}=\frac{g_u}{g_l}\exp \left(-\frac{h{\nu }_0}{k{T}_{\mathrm{ex}}}\right),$$

$$N_l=\frac{N_{(1,1)}}{1+\frac{g_u}{g_l}\exp \left(-\frac{h{\nu }_0}{k{T}_{\mathrm{ex}}}\right)}\text{leavevmode}.$$

$${\alpha }_{\nu }=\frac{c^2}{8\pi }\frac{1}{{\nu }_0^2}\frac{g_u}{g_l}{n}_l{A}_{ul}\left[1-\exp \left(-\frac{h{\nu }_0}{k{T}_{\mathrm{ex}}}\right)\right]\varphi (\nu ),$$

$$\int \varphi (\nu )d\nu =1.$$

$$\tau$$

$${\tau }_{(1,1)}$$

$$A_{ul}=\frac{64{\pi }^4{\nu }_0^3}{3h{c}^3}\frac{1}{g_u}|{\mu }_{ul}{|}^2,$$

$$|{\mu }_{ul}{|}^2={\mu }^2\frac{K^2}{J(J+1)},$$

$$\varphi (\nu )=\sum _i\frac{w_i}{\sqrt{2\pi }{\sigma }_{\nu ,i}}\exp \left[-\frac{(\nu -{\nu }_0-\delta {\nu }_i{)}^2}{2{\sigma }_{\nu ,i}}\right].$$

$$\mathrm{argmin}$$

$$I_{\nu }=(1-e^{-{\tau }_{\nu }})B_{\nu }(T_d)$$

$${\tau }_{{\nu }_0}={\kappa }_{{\nu }_0}\mu m_{\mathrm{H}}N({\mathrm{H}}_2)$$

Full text

stars:formation — ISM:molecules — ISM:individual (OrionA molecular complex) — ISM:individual (NGC1333) — ISM:individual (B18) — ISM:individual (L1668)

††journal: ApJS††software: Astropy (Robitaille et al., 2013), Matplotlib (Hunter, 2007), pyspeckit (Ginsburg & Mirocha, 2011)

1 Introduction

The past several years have seen tremendous advancements in our ability to characterize the structure of nearby molecular clouds and the substructures in which dense star-forming cores are born. Within 500 pc of the Sun, nearly all the ongoing, predominantly low-mass star formation is contained within a ring of young stars and star-forming regions named the Gould Belt. These nearby star-forming clouds range in size and activity, from the nearly inactive Pipe Nebula, through the star-forming but largely quiescent Taurus molecular cloud, to stellar group and cluster-forming clouds like Perseus and Ophiuchus. The Orion molecular cloud is our nearest example of a high-mass star-forming region. This variation in star formation activity, coupled with their relatively nearby distances, make the Gould Belt clouds excellent survey targets to understand how star formation proceeds in different environments, and to constrain star formation theories and simulations through detailed observations.

Consequently, the Gould Belt star-forming regions have been surveyed by several Key Programs in the infrared and sub-millimeter continuum. The infrared surveys have characterized the embedded young stellar populations, identifying the loci of active star formation and constraining the lifetimes of the starless and embedded stages of star-forming cores (the Spitzer Gould Belt and c2d surveys; Enoch et al. 2008; Harvey et al. 2008; Kirk et al. 2009 and Evans et al. 2003, 2009; Dunham et al. 2014). At longer wavelengths, the emission from cold dust within dense, star-forming cores, and the larger filaments and clumps that form the cores and young stellar objects (YSOs), are being surveyed extensively through Gould Belt Legacy surveys with the Herschel Space Observatory (Herschel GBS; André et al., 2010) and the James Clerk Maxwell Telescope (JCMT GBLS; Ward-Thompson et al., 2007). The Herschel GBS has revealed the striking prevalence of filaments within star-forming regions as seen in the dust emission (André et al., 2010; Könyves et al., 2010). While the presence of filaments in star-forming regions was not a new discovery (e.g., Schneider & Elmegreen, 1979), their ubiquity at both low and high column densities suggests that they are integral in the build-up of the required mass to form stars. Furthermore, the Herschel GBS has found gravitationally unstable cores predominantly within (apparently) gravitationally unstable filaments. Previous studies have shown that dense molecular cores are mostly found above a visual extinction threshold of $A_{\mathrm{V}}>5\sim 7$ mag (N(\mathrm{H}_{2})$$\sim 5-7\times 10^{21} cm*-2*; e.g., Onishi et al. 1998; Johnstone et al. 2004; Kirk et al. 2006). With a suggested characteristic filament width of $\sim 0.1$ pc (Arzoumanian et al., 2011, but see Panopoulou et al. 2017), this extinction threshold agrees within a factor of two with the critical mass per unit length of an isothermal, thermally-supported cylinder ($M_{\mathrm{line,crit}}\sim 16$ M⊙ pc*-1* at 10 K; Inutsuka & Miyama 1997; Hocuk et al. 2016). Therefore, understanding filament instability might be central to understanding the formation of dense cores and, consequently, stars.

If the massive and long filaments found by the Herschel GBS are an important piece of the star formation puzzle (André et al., 2014), then we must understand what are their main formation mechanisms, evolutionary drivers, and fragmentation triggers. The power of the legacy surveys described above lies in the large areal coverage and consistency in observing strategies between nearby molecular clouds. A major gap in the present data is a comparable survey to characterize the dense gas properties. Continuum data show the dust column density structure, but do not provide information about the gas kinematics. Kinematics and gas temperatures are key to understanding the history and future fate of star-forming material.

Kinematic studies of filaments are relatively few (e.g., Pineda et al., 2011; Hacar & Tafalla, 2011; Hacar et al., 2013; Pineda et al., 2015; Hacar et al., 2016; Henshaw et al., 2016a). The CARMA Large Area Star Formation Survey (Storm et al., 2014) targeted several regions within nearby Gould Belt clouds, at higher resolution but on a smaller scale relative to the JCMT and Herschel continuum surveys. The JCMT GBLS has produced the widest survey thus far, using observations of CO (3-2) and its isotopologues to trace moderately dense gas ($10^{3}$ cm*-3*, e.g., Graves et al. 2010). CO, however, is frequently very optically thick and the gas phase abundance can be severely affected by depletion onto dust grains (e.g., Caselli et al. 1999; Tafalla et al. 2002; Christie et al. 2012) in star-forming filaments and cores, and therefore it is not a useful tracer of the dense gas.

The lower NH3 metastable inversion transitions primarily trace gas of density $n\geq 2\times 10^{3}$ cm*-3*(for gas at 10 K, Shirley 2015). Furthermore, NH3 does not typically suffer from depletion. Since NH3 observations allow us to identify actual high volume density features rather than column density peaks, NH3 is then an ideal tracer of the hierarchical nature of star-forming regions, and a powerful probe of the kinematics of star-forming filaments and cores.

\floattable\twocolumngrid

\floattable

Past NH3 observations of dense gas were primarily pencil-beam pointings or small maps around interesting targets (e.g., Jijina et al. 1999; Rosolowsky et al. 2008). The NH3-identified cores are generally cold and quiescent (Goodman et al., 1998), either in early stages of collapse or near-critical equilibrium (Benson & Myers, 1983), and often show significant velocity gradients (Goodman et al., 1993). The local environment influences key core properties, such as the ratio of thermal vs. non-thermal motions and temperature (Myers et al., 1991; Jijina et al., 1999; Foster et al., 2009), both important parameters for investigations of core stability and collapse. Larger, more sensitive maps are needed to relate the physical properties of the cores to their surrounding environment. These recent studies have highlighted the sharp transition between turbulent and thermal line widths in dense cores (Pineda et al., 2010; Seo et al., 2015), identified incongruities between structures identified via continuum and line emission (Friesen et al., 2009), shown that cores in more active environments, such as Orion, tend to be more unstable, despite having similar internal motions (Li et al., 2013), identified large-scale gravitational instability in a young filamentary cluster (Friesen et al., 2016), and revealed evidence for filamentary accretion in Serpens South (Friesen et al., 2013, also in N2H+ by Kirk et al. 2013).

Here, we describe the Green Bank Ammonia Survey (GAS), an ambitious Large Project to survey all the major, nearby ($120<d<500$ pc), northern Gould Belt star-forming regions with visual extinctions $A_{\mathrm{V}}\gtrsim 7$ (matching the extinction threshold for dense cores in continuum studies ) in emission from NH3, as well as the carbon-chain molecules C2S, HC5N, and HC7N. Within this distance range, the GBT beam at 23 GHz (32″ FWHM) subtends $0.02\ \mathrm{pc}-0.08\ \mathrm{pc}$, able to resolve both the Jeans length (0.12 pc for $n=10^{4}$ cm*-3* and $T=10$ K) and the $\sim 0.1$ pc typical filament width suggested by Arzoumanian et al. (2011), and is furthermore well-matched to the resolution of Herschel SPIRE at 500 µm ($\theta_{\mathrm{beam}}=36$″ FWHM; Griffin et al., 2010). This coordinated, large-scale survey of NH3  in the Gould Belt clouds will play a key role in understanding the evolution of dense gas in star-forming regions as a function of environment, expanding in a consistent way the NH3 studies described above, and will complement and enhance the Gould Belt continuum surveys.

In this paper, we present the first data release (hereafter referred to as DR1) of a subset of the survey targets, describing in detail the calibration, imaging, and analysis pipelines, and discussing some general trends in the data. The target regions presented are Barnard 18 in Taurus (hereafter B18), NGC 1333 in Perseus, L1688 in Ophiuchus, and Orion A (North). Future research will focus on the key scientific goals of GAS, including analyzing the structure and stability of the dense gas, the dissipation of turbulence from large to small scales, and the evolution of angular momentum in dense cores.

In § 2, we discuss the source selection, observations, data calibration and imaging. In § 3, we describe in detail the hyperfine modeling of the NH3 emission and the determination of key physical parameters of the emitting gas. In § 4, we present the moment maps and noise properties of the regions included in this first data release, as well as the results from the NH3 line modeling. We present a summary in § 5.

2 Data

2.1 Source selection and description

The goal of GAS was to map NH3 emission toward all areas within the nearby ($d<500$ pc), GBT-visible, star-forming molecular clouds of the Gould Belt with extinctions $A_{\mathrm{V}}>7$ mag [N(H2) $\gtrsim 6.7\times 10^{21}$ cm*-2*, assuming $N(\mathrm{H}_{2})=9.4\times 10^{20}$ cm*-2* ($A_{\mathrm{V}}$ mag*-1*)] . The initial map extents were identified in each cloud using either publicly released Herschel GBS 500 µm continuum data, unreleased JCMT GBLS 850 µm continuum data, or extinction maps derived from 2MASS data (Dobashi et al., 2005; Ridge et al., 2006). Further refinement of mapping targets was done as new data were publicly released (e.g., Herschel GBS).

For our initial observing season we focused on four regions that span a range of physical environments and star formation activity: B18, NGC 1333, L1688, and Orion A (North). These regions comprise the first data release (DR1) of GAS.

B18 is a filamentary structure located in the south-east of the Taurus molecular cloud, and runs roughly parallel in projection with the molecular gas structure containing Heiles Cloud 2 and the B213/L1495 filament. B18 has a projected length of $\sim 3\arcdeg$ ($\sim 7$ pc) with a total mass of $380-440$ M⊙ (Heyer et al., 1987; Mizuno et al., 1995), giving a mean gas number density $n\sim 300$ cm*-3*, and contains $\sim 9$ molecular clumps as traced by CO (Heyer, 1988). Based on the numbers of associated young stars, the star formation efficiency (SFE), where $SFE=\mathrm{current}\ M_{*}/M_{\mathrm{gas}}$, is $\sim 1-2$ % (Mizuno et al., 1995).

NGC 1333 is the most active star-forming region within the Perseus molecular cloud and contains approximately 450 M⊙ of molecular gas (Warin et al., 1996) and $\sim 150$ young stars (Gutermuth et al., 2008) in a region $\sim 1$ pc2. The mean gas density is therefore a few $\times 10^{3}$ cm*-3*. Rough estimates suggest the star formation efficiency is $\sim 13-15$ % in the region over the past $10^{6}$ years (Walsh et al., 2006; Jørgensen et al., 2008).

L1688 is the centrally-concentrated, dense hub of the Ophiuchus molecular cloud and is $\sim 1-2$ pc in radius and $\sim 2\times 10^{3}$ M⊙ in mass (Loren, 1989). Similarly to NGC 1333, the mean gas density is approximately a few $\times 10^{3}$ cm*-3*. L1688 contains regions of extremely high visual extinction, with $A_{\mathrm{V}}\sim 50-100$ mag (e.g., Wilking & Lada, 1983), and over 300 young stellar objects (Wilking et al., 2008). The star formation efficiency of the dense gas cores, identified in submillimeter continuum emission from dust, is $\sim 14$ % (Jørgensen et al., 2008).

Orion A (North), an integral-shaped filamentary structure of dense gas that extends both north and south of the Orion Nebula and its associated young stellar cluster, including the Trapezium stars, is our nearest example of ongoing high-mass star formation. The filament is characterized by varying physical conditions along the $\sim 12$ pc length mapped for DR1. A compact, narrow ridge north of the Orion Nebula, the filament becomes wider and less dense toward the south, and contains a total mass of $\sim 5\times 10^{3}$ M⊙ in the area mapped here (e.g., Bally et al., 1987; Bally, 2008). Overall, the star formation efficiency in Orion A is $\sim 3-5$ % in gas with extinction $A_{\mathrm{K}}>2$, with significantly higher values (up to several tens of percent) found toward young protostellar groups and clusters (Megeath et al., 2015).

Table 1 lists the assumed distances and number of footprints (areal coverage of 10′ $\times$ 10′) completed for the GAS DR1 regions.

\floattable

\twocolumngrid

2.2 Observations

For all regions, observations were performed using the seven-beam K-Band Focal Plane Array (KFPA) at the GBT, with the VErsatile GBT Astronomical Spectrometer (VEGAS) backend. GAS uses VEGAS configuration Mode 20, allowing eight separate spectral windows per KFPA beam, each with a bandwidth of 23.44 MHz and 4096 spectral channels, for a spectral resolution of 5.7 kHz, or $\sim 0.07$ km s*-1* at 23.7 GHz. All spectral lines observed are listed in Table 2. Due to limitations on the maximum separation in GHz between spectral lines in a single VEGAS bank, the GAS setup includes six spectral lines observed in each KFPA beam, plus the C2S $2_{1}-1_{0}$ line in a single central beam. Observations were performed using in-band frequency switching to maximize on-source time, with a frequency throw of 4.11 MHz.

GAS maps were generally observed in 10′ $\times$ 10′ footprints, scanning in Right Ascension (R.A.), with scan rows separated by 13″ in Declination (Dec.) to ensure Nyquist sampling. Each 10′ square submap was observed once. To best match the observed regions to known sub-millimetre continuum emission, some regions included submaps of 10′ $\times$ 5′, or 5′ $\times$ 10′ (where scans were in Dec., spaced by 13″ in R.A. to minimize overheads). For all maps, the telescope scan rate was 6$\farcs$2 s*-1*, with data dumped every 1.044 s. A fast frequency-switching rate of 0.348 s was used to meet the proposed sensitivity per map in the time allotted, but resulted in a small fraction (approximately a few percent) of the data being blanked due to synthesizer settling.

2.3 Flux Calibration

Flux calibration was performed using observations of the Moon and Jupiter. Because the Moon is large relative to the KFPA footprint, simple on-target and off-target observations allowed relative beam calibrations to be performed. In angular size, Jupiter was comparable to the 32″ GBT beam at 23.7 GHz. For these sources, nod observations were performed to alternate on- and off-source beams, cycling through all beams. Relative beam calibrations were then calculated in GBTIDL111http://gbtidl.nrao.edu/. Models for each solar system target, including the effects of solar illumination and apparent size, were used to determine the expected brightness temperature in the gain analysis222https://safe.nrao.edu/wiki/pub/GB/Knowledge/GBTMemos/GBTMemo273-11Feb25.pdf.

The Moon was used as the primary calibrator when available, as the planetary calibrations could vary substantially when the weather was not ideal, particularly when the on-site winds were near the 5 m s*-1* limits for KFPA operations at 23 GHz. This can be seen in Figure 1, where we show the beam- and polarization-averaged gains derived from observations of Jupiter and the Moon over the 15A observing semester. Jupiter and Moon calibrations agree well when the winds are low. The last three dates when Jupiter was observed had relatively high winds ($\gtrsim 5$ m s*-1*), lowering the pointing accuracy of the telescope . Jupiter’s angular size over the observation dates was $\sim 40$″, only slightly larger than the 32″ (FWHM) GBT beam. Small pointing errors introduced by high winds are thus able to impact greatly the derived gains. In contrast, we show in Figure 1 that the gains determined from the Moon were very stable over the 15A semester.

Figure 2 shows the gains as a function of beam and polarization, averaged over the observing semester. The final three Jupiter measurements that were impacted by winds have been omitted from the analysis. Although the gains vary between beams and polarizations, they remained consistent over the observing period with only small variation. The gains determined and used for calibration (next section) are listed in Table 3, along with the standard deviation of the measurements for each beam over the semester.

Two known targets with bright NH3 (1,1) lines (L1489_PPC and L43) were also observed frequently in a single beam. We found that the peak line brightness varied by $<10$ % (1-$\sigma$) for these sources over the semester, in agreement with the calibration uncertainty estimated from our directly-measured beam gain variations.

2.4 Data reduction and calibration

All GAS data reduction and imaging codes (next section) are publicly available on GitHub333https://GitHub.com/GBTAmmoniaSurvey/GAS and can be easily adapted for similar observations.

The data reduction was performed using a Python wrapper that called the GBT KFPA data reduction pipeline (Masters et al., 2011). The data were calibrated to main beam brightness temperature ($T_{\mathrm{MB}}$) units using the relative gain factors for each of the beams and polarizations discussed above and listed in Table 3. Based on the 1 $\sigma$ variation in the beam gain calibrations, we estimate our calibration uncertainty to be $\sim 10$ %.

2.5 Imaging

We grid the On-The-Fly (OTF) observations following the prescription by Mangum et al. (2007). For each region, our mapping algorithm establishes a world coordinate system and corresponding pixel grid that covers the footprint of the observations, using a pixel size of 1/3 of the beam FWHM ($\theta_{\mathrm{beam}}$) of the NH3 (1,1) line. We then calculate the value of each map pixel as the weighted average of spectra at the center of the map pixel, using a truncated, tapered Bessel function for the weighting scheme. Each integration is given a weight, $W$, based on the location of the observation relative to the center of the pixel in the resulting map and the observations rms:

[TABLE]

where $a=1.55(\theta_{\mathrm{beam}}/3)$, $b=2.52(\theta_{\mathrm{beam}}/3)$, $r$ is the angular distance between the observed datum and the map datum, and $T_{\mathrm{sys}}$ is the recorded system temperature. Integrations with $r>\theta_{\mathrm{beam}}$ are excluded from the average, as suggested by Mangum et al. (2007). For simplicity of comparison, we grid data for all lines to a common world coordinate system.

Online Doppler tracking is only applied for the NH3 (1,1) line. For all other lines, we apply a frequency shift in the Fourier domain to each individual integration to align all spectra before convolution (gridding), i.e. Doppler tracking is done in software. Baseline subtraction is done in both the individual integrations (linear baseline) and in the individual pixels in the final data cubes (Legendre polynomials of order 3). In the latter case, sections of the spectrum with significant emission or within 2 $\mbox{km s}^{-1}$ of the expected systemic velocity are excluded from the fit.

The KFPA on the GBT is a sparse array, with seven beams arranged in a hexagonal pattern, and beam centers separated by 94.88″ on the sky. The mapping method described in § 2.2 ensured Nyquist-spaced sampling over 10′ $\times$ 10′ footprints on the sky. At the edges of these footprints, coverage is variable, and the rms noise increases substantially as noted previously. For the final moment maps and parameters maps for each region, we mask pixels within 3 pixels of the map edges to remove the noisiest regions.

In Figure 3, we show histograms of the rms noise for all regions and lines observed. The rms noise map for the NH3 (1,1) observations toward each DR1 region are presented in Appendix A.

In all regions, the noise histograms have strong peaks that agree well for the NH3 (1,1), (2,2), (3,3), HC5N $9-8$, and HC7N $21-10$ observations, with typical 1 $\sigma$ noise values of $\lesssim 0.1$ K per velocity channel across the DR1 regions. Note that the typical noise in individual map footprints can vary slightly depending on the conditions in which the region was mapped. Figure A show that in general this variation is minimal. The C2S observations have greater rms noise values since C2S was observed with a single beam rather than all seven beams in the KFPA. The HC7N $22-21$ data are noisier than the data of other lines observed with all KFPA beams, a systematic difference seen over all observations.

3 NH3 line analysis

3.1 Line fitting

We calculate the properties of the column of ammonia molecules using a forward-modeling approach in the pyspeckit package (Ginsburg & Mirocha, 2011), which builds from the works of Rosolowsky et al. (2008) and Friesen et al. (2009) relying on the theoretical framework laid out in Mangum & Shirley (2015). Under this approach, we consider a beam-filling slab of ammonia gas with a set of physical properties including para-ammonia column density ($N_{\mathrm{p-NH}_{3}}\equiv N$ for notational compactness in this section), kinetic temperature ($T_{K}$), velocity dispersion $\sigma_{v}$, and line of sight velocity in the (kinematic) LSR frame ($v_{\mathrm{LSR}}$). While we observe the $(J,K)=(3,3)$ line of ortho-ammonia, we do not directly incorporate any information from this line into our fit because of unknown ortho-to-para ratios that might differ from LTE (e.g., Faure et al., 2013).

Given a total column density for the molecule and a velocity dispersion for the slab, we use standard relations to determine the fraction of molecules found in each of the states, namely $(J,K)$ rotational states of the species and the upper ($u$, antisymmetric inversion symmetry) and lower ($l$, symmetric inversion symmetry) within those states. Our model incorporates three distinct temperatures. The kinetic temperature $T_{K}$ describes the velocity distributions of the particles in the system, notably of the collider species H2 and He. The rotation temperature, $T_{\mathrm{rot}}$ is the temperature that characterizes the rotational level population. Our work usually focuses on the temperature that defines the population ratio between the (1,1) and (2,2) states. Finally, $T_{\mathrm{ex}}$ is the radiative excitation temperature of the observed spectral lines.

NH3 fitting usually proceeds under the assumption that the (1,1) and the (2,2) line have the same radiative excitation temperature and that the frequencies are approximately equal and $\ll kT_{\mathrm{ex}}/h$. The former assumption is only approximately true: using RADEX (van der Tak et al., 2007) for $n=10^{4}\leavevmode\nobreak\ \mathrm{cm}^{-3}$ and $T_{K}=15\ \mathrm{K}$ and $N(\mathrm{p-NH}_{3})=10^{14}\leavevmode\nobreak\ \mathrm{cm^{-2}}$ gives $T_{\mathrm{ex},(1,1)}=8.5\leavevmode\nobreak\ \mathrm{K}$ and $T_{\mathrm{ex},(2,2)}=6.9\leavevmode\nobreak\ \mathrm{K}$. Each of the inversion transitions consists of multiple, resolved hyperfine components and we further assume a constant $T_{\mathrm{ex}}$ for all these transitions, which precludes our model representing hyperfine anomalies (e.g., Stutzki et al., 1984). We see, however, no evidence in our survey that these effects are important.

Our observations consist of antenna temperature as a function of frequency, and so we generate a model of the spectral line given the input parameters ($N$, $T_{K}$, $T_{\mathrm{ex}}$,$\sigma_{v}$, $v_{\mathrm{LSR}}$). We assume a single radiative excitation temperature, $T_{ex}$, characterizes both the (1,1) and the (2,2) transitions. We generate a model of optical depth as a function of frequency $\tau(\nu)$ to produce a model spectrum of the main-beam radiation temperature (assuming the beam filling fraction $=1$):

[TABLE]

The optical depth, $\tau(\nu)$, is a function of the input parameters. To derive $\tau(\nu)$ from the input parameters, we first approximate the molecular population as having most of the p-NH3 in the (1,1) and (2,2) states so that we can define a $T_{\mathrm{rot}}$ by the ratio of the column densities in the (2,2) and (1,1) states:

[TABLE]

where the $5/3$ is the ratio of the rotational statistical weights of the states ($g_{J}=2J+1$), and $N(J,K)$ is the sum of the column densities of the two inversion levels. From here we have that $T_{\mathrm{rot}}$ can be related to the kinetic temperature by analytically solving the equations of detailed balance (Swift et al., 2005):

[TABLE]

Here $T_{0}=41.5\leavevmode\nobreak\ \mathrm{K}$ is the energy difference between the $(2,2)$ and the $(1,1)$ states. We then use this rotational temperature to calculate a partition function for the p-NH3 molecule, though the next metastable level is the $(4,4)$ line at 178 K above the ground state and is negligibly populated ($\tau<1$) for $T_{K}<65\mbox{ K}$. The partition function is given by

[TABLE]

where $B_{\mathrm{rot}}$ and $C_{\mathrm{rot}}$ are the rotational constants of the ammonia molecule, 298117 and 186726 MHz, respectively (Pickett et al., 1998).

Since the excitation temperature defines the population ratio of the lower to the upper state, we can calculate $N_{l}$ from the total column density in the $(1,1)$ state, $N_{(1,1)}$.

[TABLE]

so

[TABLE]

Given the column density in the lower state, we then calculate a total, integrated optical depth in the (1,1) and (2,2) lines using standard radiation relationships. In particular, the opacity per unit frequency,

[TABLE]

where $\nu_{0}$ is the line rest frequency, $n_{l}$ is the volume density in the lower state of the transition, $g_{u}$,$g_{l}$ are statistical weights of the upper and lower states of a transition and whose ratio is unity for each pair of inversion transitions, $A_{ul}$ is the Einstein $A$ coefficient, and $\phi(\nu)$ is the line profile function defined such that

[TABLE]

We integrate over line of sight and frequency to get the optical depth as a function of the column density in the lower transition of the spectral line:

[TABLE]

or, in terms of $N_{(1,1)}$, which is the total column of molecules in the symmetric and antisymmetric (upper and lower) states of the J=1 K=1 level,

[TABLE]

The optical depth for the (2,2) transition has the same form. The column densities are related to the total column $N$ through the partition function via $T_{\mathrm{rot}}$.

As per Mangum & Shirley (2015),

[TABLE]

where $\mu_{ul}$ is the dipole moment of a given transition, which is related to the dipole moment of the molecule $\mu^{2}$ by

[TABLE]

where $\mu=1.468\times 10^{-18}\mathrm{\leavevmode\nobreak\ esu\leavevmode\nobreak\ cm}$ for ammonia.

Finally, we calculate the line profile function $\phi(\nu)$ from the hyperfine structure of the different transitions. Since each hyperfine transition is assumed to have the same radiative excitation temperature, $\phi(\nu)$ is given as the weighted sum of Gaussian line profiles:

[TABLE]

The weights, $w_{i}$, and frequency offsets from the rest frequency, $\delta\nu_{i}$, are set by quantum mechanics and are tabulated in Mangum & Shirley (2015). The frequency widths of the lines are calculated in terms of the velocity widths from the Doppler formula: $\sigma_{\nu,i}=(\nu_{0}+\delta\nu_{i})\sigma_{v}c^{-1}$.

This model provides a complete spectrum $T_{\mathrm{MB}}$ as a function of input parameters. We then use a non-linear gradient descent algorithm (MPFIT; Markwardt, 2009) implemented in Python to calculate the optimal values of the input parameters and their uncertainties (e.g., Press et al., 2002).

[TABLE]

Non-linear least squares fitting can be unstable unless the algorithm is provided with a good set of initial conditions. In particular, the $v_{\mathrm{LSR}}$ value of the fit is critical, so we use guesses based on the velocity centroid of the line (or of the center hyperfine component for ammonia lines). The velocity dispersion is taken initially to be the intensity weighted second moment of the velocity around the centroid. We adopt $\log_{10}(N/\mathrm{cm}^{-2})=14.5$, $T_{K}=12$ K, and $T_{\mathrm{ex}}=3$ K as our initial guesses. The fit is performed serially for each position in the map, starting from the pixel with the highest integrated intensity and working down in reverse order of integrated intensity. If a given position has a neighboring pixel with a good fit, the parameters of the successful fit are used as initial guesses for the next optimization.

In this release, we do not attempt to fit multiple velocity components or non-Gaussian velocity profiles to our spectra. Such approaches require human supervision under current implementations (Henshaw et al., 2016b). Our fits produce good quality results for most of the studied regions, but some lines of sight in Orion A show good evidence for needing improved models for their underlying velocity components.

Recent work by Estalella (2017) presented the HfS hyperfine fitting tool. The assumptions embedded in our approach are essentially equivalent to those embedded in that tool, though Estalella (2017) used a relationship between $T_{\mathrm{rot}}$ and $T_{K}$ based on the collision coefficients in Maret et al. (2009), which differs from the Swift et al. (2005) expression by $>10\%$ at $T>50\mbox{ K}$. The optimization approach differs significantly, in that we use a gradient descent non-linear regression, whereas the HfS code uses a genetic algorithm. Our approach requires more accurate guesses for the initial conditions, but the strategy we have adopted leads to good convergence. The uncertainty estimates both follow standard estimates based on $\chi^{2}$ fitting.

3.2 $v_{\mathrm{LSR}}$-tracked integrated intensity maps

When calculating the integrated intensity map of a specific line, it is common practice to define a fixed velocity range that includes the emission of all hyperfine components for a single region. The maps presented here, however, cover large areas of individual star-forming clouds. In these regions, large velocity gradients make a fixed range sub-optimal.

To create the NH3 (1,1) and (2,2) integrated intensity maps, we therefore first produce a model emission cube given the best fit models over the entire map. We set a threshold value $>0.0125$ K (i.e., slightly greater than machine precision) to identify the spectral channels with model line emission. We then select all the voxels (3D “pixels”) with brightness greater than the threshold to calculate the integrated intensity map. In regions where no line fit was found, we calculate the integrated intensity within a spectral window defined by the mean cloud $v_{\mathrm{LSR}}$ and $\sigma_{v}$. We calculate the uncertainty at each pixel based on the channel range thus identified.

This method for calculating the integrated intensity while following the velocity gradients in a star-forming region depends on having a good model, which is usually the case for NH3 emission in the clouds presented here. There are locations, however, where the emission is not perfectly fitted with a single component. In those cases, this method may not include some of the flux in the final moment maps. A detailed study of the gas kinematics will be done in future papers focusing on individual regions.

Since the spatial and kinematic distribution of NH3 (3,3) and the carbon-chain molecules do not typically follow that of NH3 (1,1) and (2,2), and are significantly less extended, where detected we use a single velocity range to produce integrated intensity maps of these transitions.

3.3 Data masking

For all regions, NH3 (1,1) and (2,2) lines were fit where the NH3 (1,1) line had a $\mathrm{SNR}>3$ through the property map analysis described above. Here $\mathrm{SNR}=T_{\mathrm{peak}}/\mathrm{rms}$, where $T_{\mathrm{peak}}$ is the peak brightness temperature in the NH3 (1,1) cube at a given pixel. At this SNR, however, not all line parameters could be fit well over all pixels. We consequently mask the output parameter maps to remove poor fits to the NH3 lines. Since $v_{\mathrm{LSR}}$, $\sigma_{v}$, and $T_{\mathrm{ex}}$ require good fits to the NH3 (1,1) line only, we mask pixels with a $\mathrm{SNR}<3$ in NH3 (1,1) integrated intensity. The parameters $T_{\mathrm{K}}$ and $N$ require detections of the NH3 (2,2) line, and for these property maps we additionally mask any pixels with a $\mathrm{SNR}<3$ in the peak NH3 (2,2) line intensity. Finally, we apply a mask that removes any remaining isolated pixels that are unconnected with larger-scale features and are generally noisy datapoints. Additional masking based on the returned uncertainty in the fit parameters may be done for detailed analysis in future papers.

4 Results

4.1 NH3 moment maps

In Figures 4, 5, 6, and 7, we show the NH3 (1,1) integrated intensity maps of B18, L1688, NGC 1333, and Orion A (North). We further overlay the NH3 (1,1) integrated intensity contours on maps of the H2 column density, $N(\mathrm{H}_{2})$, in each region, derived from spectral energy distribution (SED) modeling of continuum emission from dust observed with the Herschel Space Observatory (A. Singh et al., in preparation). The SED modeling is described in more detail in § 4.3. The locations of Class 0/I and flat spectrum protostars are highlighted, based on infrared analysis of Taurus (Rebull et al., 2010), Perseus and Ophiuchus (Dunham et al., 2015), and Orion (Megeath et al., 2012).

In general, the NH3 (1,1) integrated intensity emission follows closely the H2 column density as traced by dust continuum emission. The GAS data highlight individual cores and filamentary structure, as expected, but at the achieved sensitivity level also reveal low-SNR emission that extends between the higher column density structures, and follows, by visual inspection, the lower H2 column density material. As an example, Figure 4 shows how NH3 (1,1) integrated intensity peaks highlight the dense cores in Taurus B18, while extended emission matches the shape of the surrounding material seen in $N(\mathrm{H}_{2})$. Toward L1688, the central, high column density cores are surrounded by material at lower $N(\mathrm{H}_{2})$ values that fade away toward the northeast in Figure 6. Low SNR NH3 (1,1) emission is similarly more extended in this direction. Extended, low SNR NH3 emission can also be seen toward the more distant NGC 1333 and Orion A regions.

We show in Appendix B the integrated intensity maps of NH3 (2,2) and NH3 (3,3) (where detected) toward the DR1 regions.

In appearance, the NH3 (2,2) integrated intensity maps match closely the NH3 (1,1) maps, but with lower SNR.

In NGC 1333, L1688, and Orion A, NH3 (3,3) emission is detected, indicative of warmer dense gas present in these regions. In NGC 1333, the NH3 (3,3) emission is compact and centered on two well-known sources, the bright Class 0/I object SVS 13 and Herbig-Haro object HH12. In L1688, the NH3 (3,3) emission is more diffuse, overlapping the NH3 (1,1) contours in the north-west but primarily visible on the western edge. This region is an interface between the cold, dense gas of the Ophiuchus A core and an edge-on photodissociation region associated with the B2V star HD 147889 (Loren & Wootten, 1986; Liseau et al., 1999). The NH3 (3,3) emission in Orion A peaks at the BN/KL region and traces the associated filamentary structure extending north. The NH3 (3,3) emission is also prominent toward the Orion bar, highlighting this feature significantly better than NH3 (1,1).

4.2 Detections of other lines

In the previous sections, we have focused primarily on emission from the NH3  transitions toward the DR1 regions. We additionally mapped the rotational transitions of several carbon chains, which are listed in Table 2. Not all lines were detected in all regions, and the HC7N transitions were not detected above the rms noise or in an integrated intensity map in any region. We show in Appendix C integrated intensity maps of the additional detected lines in each region. In Figure 8, we show the spectrum of each observed line averaged over each DR1 region after masking of the edge pixels as described in § 2.5.

Toward B18, both HC5N 9-8 and C2S $2_{1}-1_{0}$ are detected. For both molecules, the distribution of emission is significantly offset from NH3 (1,1), generally peaking away from the NH3 emission peaks. While HC5N is seen interior to some NH3 (1,1) contours, C2S, where seen, appears to surround the NH3 emission. Similar distributions of NH3 and C2S have been identified toward the L1495-B218 filaments in Taurus (Seo et al., in preparation). This offset between emission from NH3 and carbon-chain molecules is common in nearby star-forming regions and is generally explained by the depletion of C-bearing molecules onto dust grains at high densities, while NH3 remains in the gas phase. Asymmetric offsets like those seen in B18 might be indicative of the contraction of dense gas (Tafalla et al., 2004) or accretion of lower density material onto the dense cores (Friesen et al., 2013).

HC5N is also detected toward NGC 1333 in a single compact region (slightly greater in extent than the beam) located to the west of the central structure traced by NH3 (1,1). While C2S is not visible in the integrated intensity map of NGC 1333 (not shown), the averaged C2S spectrum in Figure 8 shows a significant detection. Similarly, C2S is clearly detected toward L1688 in the averaged spectrum, while the integrated intensity map shows only a faint, SNR $\sim 3-4$ feature near the map center. In future analysis we will investigate extended, low-level emission of NH3 and the other observed lines through spectral stacking techniques.

4.3 Relationship of NH3 emission to $N(\mathrm{H}_{2})$

4.3.1 Column density of H2

The H2 column density, $N(\mathrm{H}_{2})$, used here and in the maps above was determined in two steps (A. Singh et al., in preparation). First, the optical depth of dust was obtained by fitting a modified blackbody to the spectral energy distribution (SED) of the thermal continuum emission from cold dust using Herschel maps at 160 µm, 250 µm, 350 µm, and 500 µm. Second, this is converted to gas column density by a simple scale factor. An additional scale factor can be used to convert to $A_{\mathrm{V}}$.

To implement the first step, we note that for dust temperature $T_{d}$ the intensity of thermal dust emission is

[TABLE]

and we adopt $\tau_{\nu}=\tau_{\nu_{0}}(\nu/\nu_{0})^{\beta}$. The dust opacity spectral index was taken to be $\beta=1.62$, the mean value found by Planck over the entire sky (Planck Collaboration XI, 2014). All maps were first convolved to the 36$\farcs$3 resolution of the 500 µm data and put on the same grid. Zero offsets were determined for each wavelength by comparison with Planck emission maps. Minor color corrections were applied during the SED fitting. The SED parameters $\tau_{\nu_{0}}$ and $T_{d}$ were determined at each pixel using $\chi$-squared minimization.

The required scale factor to convert to $N(\mathrm{H}_{2})$ can be seen via the relation

[TABLE]

where $\mu=2.8$, $m_{\mathrm{H}}$ is the mass of a Hydrogen atom, and $\kappa_{\nu_{0}}$ is the opacity (incorporating a gas-to-dust mass ratio of 100), which is taken to be 0.1 cm2 g*-1* at $\nu_{0}=1000$ GHz (Hildebrand, 1983). The considerable systematic uncertainty in the opacity propagates uniformly throughout the map into a systematic uncertainty in the scale of $N(\mathrm{H}_{2})$.

Ideally, one could scale the dust-based $\tau_{\nu_{0}}$ directly into $A_{\mathrm{V}}$, which is also proportional to dust column density. This scaling too is uncertain (Planck Collaboration XI, 2014). Here we simply use $N(\mathrm{H}_{2})=9.4\times 10^{20}$ cm*-2* ($A_{\mathrm{V}}$/mag)), which is based on data for the diffuse interstellar medium, not molecular clouds, and so the two-figure precision given belies the actual uncertainty.

Regional variations in dust spectral properties could impact the values determined through SED modeling (and thus the threshold for NH3 emission discussed below), although this is likely a small effect compared to the uncertainty in the scaling: for example, a change in $\beta$ of 0.05 produces a change in $N(\mathrm{H}_{2})$ of only $\sim 5$ %.

4.3.2 $N(\mathrm{H}_{2})$ threshold for NH3 emission

Figure 9 (top) shows a smooth increase in the fraction of pixels with NH3 emission as $N(\mathrm{H}_{2})$ increases, with some variation between regions in the rate of rise in the correlation, and in the H2 column density where NH3 is found universally. The survey initially targeted regions with extinctions $A_{\mathrm{V}}>7$, based on previous observations of NH3 in nearby star-forming regions. That this was a good choice is confirmed by the DR1 data, where above this extinction value (highlighted by the gray dashed line in the Figure, assuming $N(\mathrm{H}_{2})=9.4\times 10^{20}$ cm*-2* ($A_{\mathrm{V}}$/mag)), $\gtrsim 60$ % of pixels have detectable NH3 emission above the survey noise levels in all DR1 regions, with the exception of B18, where the threshold appears to be slightly higher. Furthermore, the survey footprints were rectangular, not following the high column density contours precisely, and so lines of sight corresponding to lower column density were surveyed too.

4.3.3 Comparison to simulations

We compared the GAS data to results calculated from a synthetic observation of NH3 emission produced by a 3D hydrodynamic simulation. The simulated emission map is based on astrochemical modeling performed by Offner et al. (2014). The abundances were computed by post-processing the simulation with 3d-pdr, a three-dimensional photodissociation region code (Bisbas et al., 2012), assuming the simulation is irradiated with a uniform 1 Draine field. We next computed the emission using radmc-3d444http://www.ita.uni-heidelberg.de/~dullemond/software/radmc-3d. Finally, to model instrumental affects we convolved the emission map with a 32″ beam for a distance of 250 pc and added Gaussian noise with $\sigma_{\rm rms}=$ 0.15 K km s*-1*, a typical value from the GAS NH3 (1,1) integrated intensity maps.

Figure 9 (middle) shows the emission fraction for the total simulation domain ($n_{\rm H_{2}}=900$ cm*-3*) and for two sub-regions ($n_{\rm H_{2}}=1,700$ cm*-3* and $n_{\rm H_{2}}=830$ cm*-3*). These distributions agree well with the GAS distributions and demonstrate that variation in the underlying gas densities can account for some of the spread between regions. Differences between regions may also result from variations in the local UV radiation field and mean column density, which we will explore in future work.

Lastly, Figure 9 (bottom) shows the cumulative distribution of $N(\mathrm{H}_{2})$ in the observed footprints of the GAS DR1 regions. A more gradual increase is seen relative to the top panel, where the presence of NH3 rises steeply above $\mathrm{log}_{10}N(\mathrm{H}_{2})/[\mathrm{cm}^{-2}]\sim 21.5$.

Figure 9 shows that observations of NH3 inversion transitions become an excellent indicator of where there is gas at higher column densities, and thus probably higher volume volume densities, but miss material at lower column density that can be detected via the thermal dust continuum emission.

4.3.4 The appearance of filaments

The ability of NH3 to highlight higher density structures enables analyses of structure sizes, masses and concentrations that are complementary to similar analyses using dust continuum emission. For example, we noted previously that filaments are prevalent within star-forming regions, with a typical filament width of $\sim 0.1$ pc suggested by profile-fitting of the mid- and far-infrared dust emission (Arzoumanian et al., 2011). The authors define the filament ‘width’ as the FWHM of a Gaussian fit to the innermost section of the filament’s radial profile. In Figure 10, we show a filament in NGC 1333 in a) NH3 (1,1) integrated intensity, and b) $N(\mathrm{H}_{2})$. In panels c) and d), we show the filament radial profile in NH3 and $N(\mathrm{H}_{2})$, respectively, averaged along the long axis shown by the straight line in panels a) and b). We then fit a Gaussian, plus constant offset, to the innermost 0.1 pc of each radial profile, matching the analysis of Arzoumanian et al., with the results shown in red. Similar to the authors’ result, we find a filament FWHM of 0.086 pc in $N(\mathrm{H}_{2})$, but show that the filament is narrower in NH3 emission by a factor $\sim 1.4$.

Locating the range of $N(\mathrm{H}_{2})$ values from panel d) in the top panel of Figure 9 shows that the NH3 emission does not measure the full column density, only the denser core. Furthermore, in $N(\mathrm{H}_{2})$ the filament is embedded within lower column density material and crossed by additional, narrow features, while the NH3 emission highlights the filament only.

Future work using the GAS dataset will identify filamentary structures systematically using algorithms such as filfinder (Koch & Rosolowsky, 2015) and DisPerSE (Sousbie, 2011). This will probe the properties of filamentary structures at different spatial scales and determine what dominates the kinematics within the filamentary structures at each scale, thus providing insight into their origin.

4.4 Property Maps

We present the property maps resulting from the NH3 line fitting described in §  3. In Figures 11, 12, 13, and 14, we show maps of $v_{\mathrm{LSR}}$ and $\sigma_{v}$ for each of the DR1 regions. The gas temperature $T_{\mathrm{K}}$ is presented for all regions in Figures 15, 16, 17, and 18. Maps of the NH3 column density, N(NH3), are shown for all regions in Appendix D.

While detailed analyses of the dense gas properties traced by NH3 in the Gould Belt will be presented in future papers, we discuss here some general results and trends in the data.

Figures 11 - 14 show significant variations in the kinematic properties traced by NH3 between the DR1 regions. In both B18 and L1688, gas velocities are concentrated around the mean cloud value with only a small spread. In B18, at the survey sensitivity, we detect NH3 toward only a small fraction of the extended gas surrounding the NH3 peaks, and the $v_{\mathrm{LSR}}$ values shown represent mainly the $v_{\mathrm{LSR}}$ of the densest NH3 structures in the region. In the other DR1 regions, we detect more of the extended NH3 emission. The $v_{\mathrm{LSR}}$ map of L1688 in Figure 13 shows that distinct variations in $v_{\mathrm{LSR}}$ are seen between some of the NH3 structures and the surrounding dense gas, but the velocity changes are small ($\lesssim 1$ km s*-1*). Toward NGC 1333, a broader distribution in $v_{\mathrm{LSR}}$ is visible in Figure 12. Unlike the variation in $v_{\mathrm{LSR}}$ between compact and extended emission in L1688, gradients in $v_{\mathrm{LSR}}$ are seen along the filamentary structures in NGC 1333. A large gradient in $v_{\mathrm{LSR}}$ is also clearly visible extending along the Orion A filament from north to south, with a smaller velocity gradient perpendicular to the filament in the south.

We show in Figure 19 the distributions of $\sigma_{v}$, $T_{\mathrm{K}}$, N(NH3), and $X(\mathrm{NH}_{3})$ over the DR1 regions. For each parameter we show the probability density function normalized such that the integral over the range is one for each region. We include only those pixels with small uncertainties in the fitted values, as described in the Figure captions.

Figure 19 (a) reveals a systematic increase in the spread of the NH3 velocity dispersion in the DR1 regions that reflects the increase in the level of star-formation activity from the least active region, B18, to the most active region, Orion A. In B18, NGC 1333, and L1688, the $\sigma_{v}$ distribution shows a narrow peak at small velocity dispersions, and a ‘tail’ of greater velocity dispersions that increases in both the distribution amplitude relative to the narrow $\sigma_{v}$ peak, as well as the maximum $\sigma_{v}$ present in the region. In Orion A, the narrow peak might still exist, but the population of larger $\sigma_{v}$ values is substantially greater. The lower limit of measured $\sigma_{v}$ in all regions is greater than the velocity resolution of the NH3 data, and will be explored further in a future paper (Pineda et al., in preparation).

Similarly to the distributions of $\sigma_{v}$, both the mean gas kinetic temperature and the spread of temperatures measured varies systematically between the DR1 regions with increasing star-formation activity. Figure 19 (b) shows that B18 contains the coldest dense gas of the DR1 regions, with a mean T_{\mathrm{K}}$$\lesssim 10 K and only a small variation about that value. In both NGC 1333 and L1688, a moderate fraction of the dense gas is similarly cold ($T_{\mathrm{K}}\sim 10$ K), but with increasing spread to higher temperatures. Previous pointed NH3 observations of dense cores in Perseus found a typical temperature of 11 K for these objects, in agreement with our data (Rosolowsky et al., 2008), whereas in NGC 1333 warm gas is concentrated around known young stellar objects. In L1688, dense gas is more likely to be $\sim 14-15$ K, as seen previously in smaller NH3 maps of some of the cores in this region (Friesen et al., 2009). Toward Orion A, only a small fraction of the dense gas is as cold as 10 K, with typical temperatures $T_{\mathrm{K}}\sim 15-20$ K, up to very large values near the Orion KL region (seen in Figure 18). At these temperatures ($>$30 K), we expect significant populations of higher-order NH3 (J,K) levels, and our temperature estimates from the (1,1) and (2,2) lines become inaccurate.

Figure 19 (c) shows the distribution of NH3 column densities over the DR1 regions. For all regions, we find a similar increase in $N(p-\mathrm{NH}_{3})$555In most works the total column density of NH3 (including ortho- and para-states) is reported, $N(\mathrm{NH}_{3})$. A simple way to estimate total column density of NH3 is: $N(\mathrm{NH}_{3})\approx 2\times N(p-\mathrm{NH}_{3})$, if the ortho-to-para ratio is the LTE value of 1. at low column densities, and a similar peak at $\mathrm{log}\ N(p-\mathrm{NH}_{3})\sim 14$. B18 shows a relatively sharp cutoff in the maximum column density present. $N(p-\mathrm{NH}_{3})$ extends to higher values in the other regions, however, with the largest column densities found in Orion A.

The resulting distribution of NH3 abundances, $X(\mathrm{NH}_{3})=N(\mathrm{NH}_{3})/N(\mathrm{H}_{2})$, shown in Figure 19 (d), have very similar peak values of $\mathrm{log}X(\mathrm{NH}_{3})\sim-8.5\ \mathrm{to}\ -8.0$. A tail to very high NH3 abundances is seen in Orion A, but might result from difficulties in accurately fitting the submillimetre continuum SED in this active region. In NGC 1333, the abundance distribution appears bimodal. We will present more detailed analysis of the dust properties and NH3 abundance distributions in an upcoming paper.

5 Summary

We have presented data from the first release (DR1) for the Green Bank Ammonia Survey, including cubes, moment maps, and property maps, toward B18 in Taurus, L1688 in Ophiuchus, NGC 1333 in Perseus, and Orion A North in the Orion molecular cloud. We furthermore describe in detail the DR1 calibration, imaging, and line fitting pipelines. The four Gould Belt clouds observed span a range of H2 column density and star formation activity. The extensive and sensitive NH3, HC5N, HC7N, and C2S observations of these regions greatly increase the total contiguous areal coverage of dense molecular gas tracers in these well-studied star-forming regions.

All cubes and maps are publicly available , as are the imaging and analysis pipelines.

Some highlights from the sensitive GAS data include:

Extended NH3 (1,1) emission is detected toward all DR1 regions, revealing the physical properties of both the dense gas associated with prestellar and star-forming cores and filaments, as well as their moderately dense molecular envelopes. 2. 2.

NH3 (1,1) emission is a particularly good tracer of material at higher $N(\mathrm{H}_{2})$, tracking well the H2 column density as traced by dust continuum emission. In three of four DR1 regions, NH3 is present above our map rms noise levels in $\gtrsim 60$ % of pixels at $N(\mathrm{H}_{2})\gtrsim 6\times 10^{21}$ cm*-2*, or $A_{\mathrm{V}}\gtrsim 7$. 3. 3.

We observe some variation between the regions in the fraction of pixels with NH3 as H2 column density increases. We compare the GAS data to results calculated from a synthetic observation of NH3 emission produced by a 3D hydrodynamic simulation. The simulated distributions agree well with the GAS distributions and demonstrate that variation in the underlying gas densities can account for some of the spread between regions. Differences between regions may also result from variations in the local UV radiation field and mean column density, which we will explore in future work. 4. 4.

Previously unmapped NH3 (3,3) emission is detected toward NGC 1333, L1688, and Orion A. Highlighting warmer, dense gas, NH3 (3,3) emission is coincident with YSOs and PDR edges. 5. 5.

The carbon-chain molecules C2S and HC5N are strongly detected toward Taurus B18, but are faint or undetected in the integrated intensity maps of the other DR1 regions. In Taurus, the spatial distributions of C2S and HC5N are offset from the NH3 emission. We detect C2S toward NGC 1333 and L1688 only in spectra that were averaged over the entire observation extent. The observed HC7N lines are not detected toward any DR1 target. 6. 6.

The distributions of gas temperature and velocity dispersion from the NH3 line fits vary systematically between DR1 regions, in step with increasing star formation activity. Lower $T_{\mathrm{K}}$ is found overall in B18, increasing through NGC 1333 and L1688, and is greatest in Orion A. Similarly, $\sigma_{v}$ reaches greater values in L1688 and Orion A relative to B18 and NGC 1333.

RKF is a Dunlap Fellow at the Dunlap Institute for Astronomy & Astrophysics. The Dunlap Institute is funded through an endowment established by the David Dunlap family and the University of Toronto. JEP, AP, ACT, and PC acknowledge the financial support of the European Research Council (ERC; project PALs 320620). EWR, PGM and CDM are supported by Discovery Grants from NSERC of Canada. SSRO acknowledges support from NSF grant AST-1510021. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. Green Bank Telescope

Appendix A Maps of the rms noise in the NH3 (1,1) observations for all DR1 regions

Fig. SetA. NH3 (1,1) rms maps

Appendix B NH3 (2,2) and (3,3) integrated intensity maps, where detected

Fig. SetB. NH3 (2,2) Integrated Intensity maps

Fig. SetC. NH3 (3,3) Integrated Intensity rms maps for L1688, and Orion A

Appendix C HC5N $9-8$ and C2S $2_{1}-1_{0}$ integrated intensity maps, where detected

Appendix D NH3 column density maps for all DR1 regions