CARMA LARGE AREA STAR FORMATION SURVEY: DENSE GAS IN THE YOUNG L1451 REGION OF PERSEUS

Figure 3 shows the integrated intensity maps for HCO⁺, HCN, and N₂H⁺ ($J=1\to 0$) (∼8'' angular resolution), along with a Herschel  250 μm image (181 angular resolution) for a visual comparison between the dense gas and dust emission. The line maps were integrated over all channels with identifiable emission. The locations of the four Bolocam 1.1 mm sources (Enoch et al. 2006) and the one compact continuum core in L1451 are marked on each image of Figure 3. Four of the five sources are located near peaks of dense gas and dust emission. While the molecules and dust are tracing similar features around those sources, the exact morphological details vary. Below, we describe the qualitative emission features, and refer to the colored rectangles in Figure 3 for reference.

Zoom In Zoom Out Reset image size

All tracers show a curved structure surrounding L1451-mm and the two nearby Bolocam sources (see the red rectangle in Figure 3), with a peak of integrated emission at L1451-mm. The southwestern edge of the curved structure has a stream of emission that extends further to the southwest (see the dark blue rectangle in upper-left panel of Figure 3), which can be seen in all the maps, although it extends the farthest in the HCO⁺, HCN, and dust maps. The two other Bolocam sources to the far east of L1451-mm are surrounded by significant molecular gas structure (see the green rectangle in the lower-left panel). The HCO⁺ emission has the largest spatial extent in this region.

The integrated emission in the three lines is less similar across the western half of L1451 compared with the eastern half. There is a strong, condensed N₂H⁺ source (see the orange rectangle in the lower-right panel) that does not appear strongly in the HCN or HCO⁺ maps, but does correspond to a peak of emission in the dust map (see Section 8.3 for more details on this source). The strongest HCO⁺ feature in the western half of the map has a weaker counterpart in the HCN map (see the purple rectangle), which appears even weaker in N₂H⁺. Finally, there is HCO⁺ emission to the northwest of the curved structure (see the cyan rectangle) that closely mimics dust emission in that region; this emission is weakly detected in HCN but not in N₂H⁺. Because the $J=1\to 0$ transition of HCO⁺ traces densities about an order of magnitude lower than the other two molecules (Shirley 2015), the regions with strong HCO⁺ and weak HCN and N₂H⁺ are likely at a lower density compared with regions where all the molecules have strong emission.

Figure 4 shows channel maps of HCO⁺, HCN, and N₂H⁺ highlighting the bulk of the emission, which occurs from ∼4.8 to 3.5 km s⁻¹, with two-channel spacing (e.g., we skip the 4.66 km s⁻¹ channel between the 4.82 and 4.50 km s⁻¹ channels of HCO⁺). In general, it is clear from the channel maps in Figure 4 that strong HCO⁺ emission is more widespread compared with HCN emission, and particularly N₂H⁺ emission (as was also evident in Figure 3). We label qualitative features in the HCO⁺ channel maps, from A through I, in the order that they appear in velocity space (with eastern sources being labeled before western sources). We then place those same labels on the HCN and N₂H⁺ maps to aid the qualitative comparison of dense-gas features given below.

Zoom In Zoom Out Reset image size

Features A and C in the eastern half of the map are traced with all molecules, with varying strength. In Figure 4, Feature A appears strongly in HCO⁺ at 4.82 km s⁻¹, faintly at 4.82 km s⁻¹ in HCN, and not until 4.47 km s⁻¹ in N₂H⁺. Even when it finally appears as N₂H⁺ emission, the emission is much more spatially concentrated than what we observe in HCO⁺ and HCN. This spatial concentration is likely because the N₂H⁺ ($J=1\to 0$) line traces higher-density, colder material compared with HCO⁺ and HCN ($J=1\to 0$) lines (Shirley 2015). A 1.1 mm Bolocam core is located within Feature A, and corresponds with a peak in the N₂H⁺ emission. Feature A moves northeast to southwest from 4.82 km s⁻¹ to lower-velocity channels as Feature C appears to its southwest. Features A and C are possibly part of the same larger-scale structure, which will be explored in the next section when we analyze the kinematics of this region. Like Feature A, Feature C is more spatially concentrated in N₂H⁺, and contains a 1.1 mm Bolocam source at a peak of N₂H⁺ emission.

Feature B contains the L1451-mm compact continuum source. It appears strongly in all molecules, although it contains a prominent ridge of emission at 4.82 km s⁻¹ in HCO⁺ and HCN that does not appear in N₂H⁺ near that velocity.

Features D, E, F, G, H, and I are all identified based on the HCO⁺ emission. HCN emission appears weakly toward all features seen with HCO⁺, while N₂H⁺ only shows faint emission in one channel for Features G and H at 3.84 km s⁻¹. The descriptions below are based on HCO⁺. Feature D appears to the west of Features A and C, and to the south of Feature B, where it is seen at 4.50 km s⁻¹ as an elliptical feature. Feature E appears as a prominent, round emission feature at 4.50 km s⁻¹, with more extended emission in channels surrounding the 4.50 km s⁻¹ peak of emission. Feature F is emission that starts just to the northwest of Feature B at 4.17 km s⁻¹, peaks at 3.83 km s⁻¹, and then appears to get fainter while extending to the southwest at the unshown 3.67 km s⁻¹ channel, while then brightening to the southwest at 3.51 km s⁻¹. Feature G emission peaks at 3.83 km s⁻¹, and appears as a stream of emission to the east of Feature E. Feature H is a streamer to the southwest of Feature B and south of Feature F, which first appears at 3.83 km s⁻¹. It persists at 3.51 km s⁻¹ and more faintly extends into two lower-velocity channels not shown in the figure. Feature I first appears at 3.51 km s⁻¹, brightens in two lower-velocity channels not shown in the figure, and is not detectable at velocities lower than that.

Feature J, referred to as L1451-west in the rest of the paper, is relatively round emission that only appears strongly in the N₂H⁺ data. It first appears at 4.79 km s⁻¹ in Figure 4 and is visible across a total of four velocity channels per hyperfine component—the structure can be seen repeating for another hyperfine component, starting in the 3.84 km s⁻¹ channel. We discuss details of the structure and kinematics of L1451-west in Section 8.3.

We detect no HCN or HCO⁺ outflow emission in any channel, which suggests that L1451 is a young region with little to no protostellar activity. Figure 5 shows an example spectrum for each molecule from the location of L1451-mm within a single synthesized beam. The conversions from Jy beam⁻¹ to K for these data are 2.47 K/Jy beam⁻¹, 2.46 K/Jy beam⁻¹, and 2.42 K/Jy beam⁻¹ for HCO⁺, HCN, and N₂H⁺, respectively.

Zoom In Zoom Out Reset image size

A small fraction of HCO⁺ and HCN spectra across L1451 shows double peaks with a variety of line-shape characteristics—we estimate that ∼3% of cloud locations show double-peak features. We cannot determine the absolute cause of the double peaks without H¹³CN and H¹³CO⁺ ($J=1\to 0$) observations at each cloud location, so analyzing these features is beyond the scope of this paper. However, the most likely scenario is self-absorption from a foreground screen of lower-density gas with a significant HCO⁺ and HCN ($J=1\to 0$) population. An infall scenario can be ruled out in many locations because infall predicts a stronger blue peak (Evans 1999) while we often observe stronger red peaks; a scenario with two-components along the same line of sight can likely be ruled out in many locations where the HCO⁺ and HCN spectra do not both show double peaks.

We qualitatively described the dense-gas morphology of L1451 in Section 4. Here we quantitatively identify dense-gas structures from the three position–position–velocity (PPV) cubes and study the hierarchical nature of L1451 with the non-binary dendrogram algorithm described in Paper I.

A dendrogram algorithm identifies the emission peaks in a data set and keeps track of how those peaks merge together at lower emission levels. This method of identifying and tracking emission structures is advantageous compared to a watershed object-identification algorithm, such as CLUMPFIND (Williams et al. 1994), when the goals include understanding how the morphology and kinematics of star-forming gas connects from large to small scales. A full discussion of the most widely used dendrogram algorithm applied to astronomical data can be found in Rosolowsky et al. (2008b); details of dendrograms and our non-binary version of the dendrogram algorithm are found in the appendices of Paper I, along with a comparison with the results from using a more standard clump-finding algorithm on our CLASSy data.

We ran our non-binary dendrogram algorithm on the HCO⁺ emission, the emission from the strongest HCN hyperfine line, and the emission from the strongest N₂H⁺ hyperfine line. For other CLASSy regions, we limited our analysis to N₂H⁺ emission because the HCO⁺ and HCN lines were complicated by protostellar outflows and severe self-absorption effects not seen in L1451. Also, for other CLASSy regions, we ran the algorithm on the N₂H⁺ hyperfine component most isolated in velocity space because the strongest hyperfine component was often not spectrally resolved from the adjacent components. The N₂H⁺ hyperfine components in L1451 were resolved in all locations, letting us use the strongest component for our dendrogram analysis. A caveat is that some bluer emission from the higher-velocity, adjacent hyperfine component, and some redder emission from the lower-velocity, adjacent hyperfine component appear in the channels of the strongest component. Because there is no blending of hyperfine emission at the same location within the cloud, we masked the emission from the adjacent hyperfine components in the individual channels of the strongest hyperfine component in which they appeared. As an example, L1451-west (Feature J) appears in much bluer channels compared to Feature A, so L1451-west emission from the higher-velocity, adjacent hyperfine component also appears in the reddest channel of the strongest hyperfine component that shows Feature A. For this example, we masked out L1451-west emission from this red channel of the strongest component.

We ran the algorithm with similar parameters to those used for the Barnard 1 analysis described in Paper I, while following the prescription presented in the Appendix for comparing data cubes with different noise levels. The critical algorithm inputs and parameters were (1) a masked input data cube with all pixels greater than or equal to 4σ intensity, along with pixels adjacent to the initial selection that are at least 2.5σ intensity, where σ is the rms level of the given data cube; (2) a set of local maxima, where an individual local maximum is greater than or equal to all its neighbors in 10'' by 10'' by three channel (0.94 km s⁻¹) spatial velocity pixels, to act as potential dendrogram leaves; (3) a requirement that a local maximum peak at least 2σ_n above the intensity where it merges with another local maximum for a structure to be considered a leaf (referred to as the "minheight" parameter below and in the Appendix), where σ_n is the rms level of the noisiest data cube (N₂H⁺ at ∼0.14 Jy beam⁻¹, in this case); and (4) a requirement of at least three synthesized beams of spatial-velocity pixels for a structure to be considered a leaf (referred to as the "minpixel" parameter below and in the Appendix). The minheight and minpixel parameters act to prevent noise features from being identified as dendrogram leaves. Branching steps are restricted to integer values of the 1σ_n sensitivity of the data (referred to as the "stepsize" parameter below and in the Appendix) for our non-binary dendrograms when comparing data sets with different noise levels. The Appendix shows that using uniform minheight and stepsize allows a comparison of dendrogram properties that minimizes the impact of noise-level differences between data cubes.

Figures 7–9 show HCO⁺, HCN, and N₂H⁺ non-binary dendrograms for L1451, respectively. The vertical axis of the dendrograms represent the intensity range of the pixels belonging to a leaf or branch. The horizontal axis is arranged with the major features identified in Figure 4 progressing from east to west; we label certain branches that are associated with major features, and provide the numeric label for the structures referred to in the upcoming discussion. The isolated leaves are presented in numerical order. The horizontal dotted line in each dendrogram represents an intensity cut at 2.5σ_n that aids in the cross-comparison of dendrogram statistics (see the Appendix) and is discussed in Section 6.3.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The HCO⁺ dendrogram contains the largest number of structures: 86 leaves and 27 branches. The HCN dendrogram contains 33 leaves and 13 branches, while the N₂H⁺ dendrogram only contains 16 leaves and 6 branches. Leaves that peak at least 6σ_n in intensity above the branch they merge directly into are colored green and referred to as high-contrast leaves.¹⁹ The strongest leaf for every molecule is at or near the location of L1451-mm in Feature B; Leaf 66 is the strongest structure in the HCO⁺ dendrogram, with a peak intensity of 2.78 Jy beam⁻¹, leaf 30 is the strongest structure in the HCN dendrogram, with a peak intensity of 2.58 Jy beam⁻¹, and leaf 10 is the strongest structure in the N₂H⁺ dendrogram, with a peak intensity of 1.62 Jy beam⁻¹.

The leaves and branches of each dendrogram represent molecular structures. We fit for the spatial properties of each structure, as we did in Paper I, using the regionprops program in MATLAB. This program fits an ellipse to the integrated intensity footprint of each dendrogram structure to determine its R.A. centroid, decl. centroid, major axis, minor axis, and position angle. Columns 2–6 in Tables 4–6 list the spatial properties of each structure. To quantify the shape of each structure, we use the axis ratio and filling factor of the fit (Columns 7 and 8 of Tables 4–6). We then define the structure size (Column 9 of Tables 4–6) as the geometric mean of the major and minor axes, assuming a distance of 235 pc when converting to parsec units.

Table 4.  HCO⁺ Dendrogram Leaf and Branch Properties

No.	R.A.	Decl.	Maj. Axis	Min. Axis	PA	Axis	Filling	Size	$\langle {V}_{\mathrm{lsr}}\rangle$	ΔVlsr	$\langle \sigma \rangle$	Δσ	Pk. Int.	Contrast	Level
	(h:m:s)	(°:':'')	('')	('')	(°)	Ratio	Factor	(pc)	(km s−1)	(km s−1)	(km s−1)	(km s−1)	(Jy beam−1)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)	(12)	(13)	(14)	(15)	(16)
Leaves
13	03:25:31.6	+30:21:34.6	18.1	14.7	105.5	0.81	0.86	0.019	4.72(2)	0.04(2)	0.21(2)	0.04(2)	1.29	2.2	4
14	03:25:01.2	+30:21:13.7	140.5	65.4	139.4	0.47	0.62	0.109	4.63(2)	0.09(2)	0.21(2)	0.07(1)	0.99	4.8	0
15	03:25:26.1	+30:21:54.4	56.7	42.0	134.6	0.74	0.71	0.056	4.55(1)	0.05(1)	0.30(1)	0.05(1)	2.12	6.0	5
16	03:25:31.1	+30:22:56.5	117.4	53.9	121.4	0.46	0.75	0.091	4.76(1)	0.06(0)	0.17(0)	0.04(0)	1.74	7.4	2
17	03:24:55.8	+30:24:31.8	22.8	13.1	81.8	0.57	0.86	0.020	4.46(3)	0.06(2)	0.41(1)	0.03(1)	1.09	3.3	3
Branches
86	03:25:09.8	+30:23:52.9	47.6	42.9	118.3	0.90	0.73	0.051	3.89(1)	0.06(1)	0.19(1)	0.03(0)	2.35	...	6
87	03:25:03.4	+30:24:39.6	50.3	38.0	54.4	0.75	0.65	0.050	4.32(2)	0.07(1)	0.29(1)	0.04(1)	1.91	...	6
88	03:25:02.7	+30:24:17.9	110.8	59.8	49.9	0.54	0.68	0.093	4.20(2)	0.10(1)	0.26(1)	0.05(1)	1.77	...	5
89	03:25:09.3	+30:24:03.8	89.4	53.3	161.4	0.60	0.65	0.079	4.03(2)	0.11(2)	0.21(1)	0.03(1)	1.77	...	5
90	03:25:05.0	+30:24:13.9	186.3	102.2	93.1	0.55	0.75	0.157	4.18(1)	0.16(1)	0.24(0)	0.04(0)	1.62	...	4

Notes. (2)–(6) The position, major axis, minor axis, and position angle were determined from regionprops in MATLAB. We do not report formal uncertainties of these values because the spatial properties of irregularly shaped objects are dependent on the chosen method.

(7) Axis ratio, defined as the ratio of the minor axis to the major axis.

(8) Filling factor, defined as the area of the leaf or branch inscribed within the fitted ellipse, divided by the area of the fitted ellipse.

(9) Size, defined as the geometric mean of the major and minor axes, for an assumed distance of 235 pc.

(10) The weighted mean V_lsr of all fitted values within a leaf or branch. Weights are determined from the statistical uncertainties reported by the IDL MPFIT program. The error in the mean is reported in parentheses as the uncertainty in the last digit. It was computed as the standard error of the mean, ${\rm{\Delta }}{V}_{\mathrm{lsr}}/\sqrt{N}$, where ΔV_lsr is the value in column 11 and N is the number of beams worth of pixels within a given object. We report kinematic properties only for objects that have at least three synthesized beams worth of kinematic pixels.

(11) The weighted standard deviation of all fitted V_lsr values within a leaf or branch. The error was computed as the standard error of the standard deviation, ${\rm{\Delta }}{V}_{\mathrm{lsr}}/\sqrt{2(N-1)}$, assuming the sample of beams was drawn from a larger sample with a Gaussian velocity distribution.

(12) The weighted mean velocity dispersion of all fitted values within a leaf or branch. The error was computed as the standard error of the mean, ${\rm{\Delta }}\sigma /\sqrt{N}$.

(13) The weighted standard deviation of all fitted velocity dispersion values within a leaf or branch. The error was computed as the standard error of the standard deviation, ${\rm{\Delta }}\sigma /\sqrt{2(N-1)}$.

(14) For a leaf, this is the peak intensity measured in a single channel in the dendrogram analysis. For a branch, this is the intensity level where the leaves above it merge together.

(15) "Contrast" is defined as the difference between the peak intensity of a leaf and the height of its closest branch in the dendrogram, divided by 1σ_n.

(16) The branching level in the dendrogram. For example, the base of the tree is level 0, so an isolated leaf that grows directly from the base is considered to be at level 0. A leaf that grows from a branch one level above the base will be at level 1, etc.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

Table 5.  HCN Dendrogram Leaf and Branch Properties

No.	R.A.	Decl.	Maj. Axis	Min. Axis	PA	Axis	Filling	Size	$\langle {V}_{\mathrm{lsr}}\rangle$	ΔVlsr	$\langle \sigma \rangle$	Δσ	Pk. Int.	Contrast	Level
	(h:m:s)	(°:':'')	('')	('')	(°)	Ratio	Factor	(pc)	(km s−1)	(km s−1)	(km s−1)	(km s−1)	(Jy beam−1)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)	(12)	(13)	(14)	(15)	(16)
Leaves
5	03:25:24.1	+30:21:03.2	71.4	41.7	162.3	0.58	0.59	0.062	4.49(1)	0.05(1)	0.13(1)	0.03(0)	1.29	3.5	4
6	03:25:25.4	+30:21:45.4	27.1	19.7	82.7	0.73	0.86	0.026	4.56(2)	0.06(1)	0.15(1)	0.04(1)	1.30	3.6	4
7	03:25:26.3	+30:22:10.6	37.0	20.7	112.3	0.56	0.80	0.031	4.60(2)	0.06(1)	0.13(1)	0.06(1)	1.51	6.0	3
8	03:25:31.6	+30:22:59.0	53.1	29.1	101.8	0.55	0.72	0.045	4.73(1)	0.04(1)	0.11(1)	0.04(1)	0.89	2.7	2
9	03:25:28.4	+30:23:09.9	46.9	26.4	18.9	0.56	0.57	0.040	...	...	...	...	0.68	2.3	1
Branches
33	03:25:03.0	+30:24:41.8	58.4	46.1	37.9	0.79	0.67	0.059	4.30(2)	0.08(2)	0.15(1)	0.05(1)	1.14	...	3
34	03:25:09.5	+30:24:03.4	92.4	47.2	158.3	0.51	0.78	0.075	4.00(2)	0.10(1)	0.11(0)	0.03(0)	1.14	...	2
35	03:25:02.7	+30:24:32.7	112.1	44.4	37.4	0.40	0.68	0.080	4.28(2)	0.11(2)	0.15(1)	0.05(1)	0.99	...	2
36	03:25:05.6	+30:24:16.7	180.1	100.9	103.4	0.56	0.75	0.154	4.17(1)	0.14(1)	0.12(1)	0.06(0)	0.85	...	1
37	03:25:24.6	+30:21:11.0	98.8	64.4	21.4	0.65	0.65	0.091	4.50(1)	0.06(1)	0.14(1)	0.04(0)	0.78	...	3

Note. Columns are the same as Table 4.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

Table 6.  N₂H⁺ Dendrogram Leaf and Branch Properties

No.	R.A.	Decl.	Maj. Axis	Min. Axis	PA	Axis	Filling	Size	$\langle {V}_{\mathrm{lsr}}\rangle$	ΔVlsr	$\langle \sigma \rangle$	Δσ	Pk. Int.	Contrast	Level
	(h:m:s)	(°:':'')	('')	('')	(°)	Ratio	Factor	(pc)	(km s−1)	(km s−1)	(km s−1)	(km s−1)	(Jy beam−1)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)	(12)	(13)	(14)	(15)	(16)
Leaves
6	03:25:26.0	+30:21:44.0	51.1	36.0	30.7	0.70	0.66	0.049	4.57(1)	0.04(1)	0.13(1)	0.03(1)	1.25	4.0	3
7	03:25:03.9	+30:25:00.3	36.8	19.6	89.7	0.53	0.68	0.031	4.45(1)	0.03(1)	0.09(1)	0.02(0)	1.40	3.5	1
8	03:25:16.1	+30:19:47.4	35.7	28.9	54.8	0.81	0.68	0.037	...	...	...	...	0.89	3.7	0
9	03:24:55.8	+30:23:23.3	19.3	11.6	61.8	0.60	0.77	0.017	...	...	...	...	0.66	2.1	0
10	03:25:07.3	+30:24:32.1	38.5	21.1	96.6	0.55	0.71	0.032	4.24(1)	0.03(1)	0.08(0)	0.01(0)	1.62	3.0	2
Branches
16	03:25:08.5	+30:24:11.8	101.3	38.3	136.7	0.38	0.71	0.071	4.11(2)	0.12(1)	0.10(0)	0.02(0)	1.18	...	1
17	03:25:04.8	+30:24:18.4	196.1	108.9	95.0	0.56	0.62	0.166	4.26(2)	0.17(1)	0.11(0)	0.03(0)	0.90	...	0
18	03:25:17.7	+30:18:54.8	141.7	42.6	101.8	0.30	0.65	0.088	4.08(2)	0.07(2)	0.13(1)	0.04(1)	0.69	...	0
19	03:25:25.6	+30:21:34.2	96.5	65.4	41.1	0.68	0.62	0.091	4.58(1)	0.04(1)	0.10(1)	0.02(0)	0.67	...	2
20	03:25:27.2	+30:21:45.9	200.6	100.1	53.0	0.50	0.56	0.161	...	...	...	...	0.53	...	1

Note. Columns are the same as Table 4.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

Histograms of the size, filling factor, and axis ratio for all leaves and branches are plotted in the top rows of Figures 10–12. All the HCO⁺ and HCN leaves are smaller than 0.12 pc, while all the N₂H⁺ leaves are smaller than 0.06 pc. HCO⁺ branches are the largest of all molecules, peaking at 0.46 pc, while HCN branches peak at 0.24 pc and N₂H⁺ branches peak at 0.17 pc. This is expected because HCO⁺ shows the most widespread emission in the channel and integrated emission maps. The filling factor for all molecular structures is between 0.45 and 0.95. The structures with a filling factor closest to unity are leaves, indicating that leaves are more likely to be elliptically shaped objects, while branches are more likely to be irregularly shaped objects. The axis ratio for all structures is between 0.19 and 0.95, showing that there is a distribution from elongated to round structures without strong differences between leaves and branches. There are no obvious differences between the high-contrast leaves and the rest of the leaves.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We use the integrated intensity footprints of the dendrogram structures in combination with the centroid velocity and velocity dispersion maps to determine the kinematic properties of the leaves and branches. The four properties present in Tables 4–6 are the mean and rms centroid velocity ($\langle {V}_{\mathrm{lsr}}\rangle$ and ${\rm{\Delta }}{V}_{\mathrm{lsr}}$, respectively), and mean and rms velocity dispersion ($\langle \sigma \rangle$ and Δσ, respectively). We illustrated how to derive these properties for leaves and branches in Section 6 of Paper I.

Histograms of $\langle {V}_{\mathrm{lsr}}\rangle$, ${\rm{\Delta }}{V}_{\mathrm{lsr}}$, and $\langle \sigma \rangle$ are plotted in the bottom rows of Figures 10–12. HCO⁺ traces a larger variation in $\langle {V}_{\mathrm{lsr}}\rangle$ compared with HCN and N₂H⁺. We attribute this to HCO⁺ being sensitive to more widespread emission away from the densest regions of L1451, and therefore tracing more widespread centroid velocities relative to the systemic velocity of the most spatially compact gas. For all molecules, ΔV_lsr of the leaves generally extends from low to moderate velocities, while it is distributed from moderate to high velocities for the branches. This indicates a trend where ΔV_lsr is generally lower for smaller structures than larger structures. This trend was also seen for Barnard 1 gas structures, and we discussed explanations in Paper I; the primary reason was the scale-dependent nature of turbulence, which causes gas parcels separated by smaller distances to have smaller rms velocity differences between them (McKee & Ostriker 2007).

The distribution of $\langle \sigma \rangle$ is similar for leaves and branches, with neither distribution showing a preference to peak at high or low velocities. This trend was also seen in Paper I for Barnard 1, indicating that $\langle \sigma \rangle$ does not strongly depend on the projected size of a structure. The peak $\langle \sigma \rangle$ for HCO⁺ is higher compared to HCN and N₂H⁺: 0.42 km s⁻¹, 0.18, and 0.13 km s⁻¹, respectively. Because HCO⁺ traces effective excitation densities of an order of magnitude lower than HCN and N₂H⁺ (Shirley 2015), we are likely observing emission from lower-density gas that is more extended along the line of sight. This could increase the line-of-sight velocity dispersions, which are expected to scale with cloud depth in a turbulent medium (McKee & Ostriker 2007). There are no obvious differences in these kinematic properties between the high-contrast leaves and the rest of the leaves.

The dendrograms in Figures 7–9 show an apparently wide variety of hierarchical complexity between molecular tracers and subregions of L1451. In this section, we quantify the hierarchical nature of the gas with tree statistics that were introduced and explained in Paper I, so that the complexity between the molecules and subregions can be quantitatively compared. Specifically, we calculate the maximum branching level, mean path length, and mean branching ratio of the entire L1451 region for each molecule in a uniform way that accounts for the differences in the noise level of each data cube (see the Appendix). We then calculate those same statistics for individual features within each dendrogram.

To compare the tree statistics of the dendrograms from the different molecular tracers, we follow the Appendix and only consider leaves and branches above a 2.5σ_n intensity cut, where σ_n is the rms of the noisiest molecular data cube. The N₂H⁺ PPV cube has the highest noise level, at ∼0.14 Jy beam⁻¹, so the cut for each dendrogram is at ∼0.35 Jy beam⁻¹ and is represented as the horizontal dotted line in Figures 7–9. Only leaves that peak at least 2σ_n above the cut are still considered in the statistics—all other leaves are marked with an "x" in the dendrograms. A branch below the cut is discounted but if the leaf directly above it is more than 2σ_n above the cut, it is counted as a leaf with a branching level of zero (e.g., leaf 64 in Figure 7). This ensures that we can compare the hierarchical structure of emission from the different molecules independent of noise-level differences introduced from the observing setup. These considerations are not needed when comparing tree statistics from the subregions of a single dendrogram.

The path length statistic, defined only for leaves, is the number of branching levels it takes to go from a leaf to the tree base. The branching ratio statistic, defined only for branches, is the number of structures that a branch splits into immediately above it in the dendrogram. We will link these tree statistics to cloud fragmentation in the discussion to follow. A more evolved region with a lot of hierarchical structure will have a higher maximum branching level and larger mean path length than a young region that is starting to form overdensities and fragment. The branching ratio of a very hierarchical region will be smaller than a region fragmenting into many substructures in a single step. This is likely an overly simplistic view, because the molecular emission, and in turn the dendrograms and tree statistics, can also be affected by projection effects and line opacity. We briefly discussed these effects in Paper I; because they are extremely difficult to accurately model over such a large area, we use the simple view that hierarchy comes from fragmentation.²⁰

The tree statistics of each dendrogram from the different molecular tracers are reported in the first section of Table 7. The mean branching ratios are 3.8, 3.3, and 3.1 for the HCO⁺, HCN, and N₂H⁺ emission, respectively. We interpret this to mean that each molecule is tracing physical structures that are fragmenting in a hierarchically similar way (e.g., a structure is most likely to fragment into about three to four substructures). The mean path length of HCO⁺ is 1.0 level longer than that of HCN and N₂H⁺, and the maximum branching level of HCO⁺ is two more than that of HCN and three more than that of N₂H⁺. This trend in fragmentation levels is likely due to the ability of each tracer to detect material at different spatial scales and physical densities. As the effective excitation density of the tracer goes down from N₂H⁺ to HCN to HCO⁺, our observations are more sensitive to more widespread emission, which means we are sensitive to more levels of fragmentation extending from the higher-density leaves (detectable with all tracers) to the lowest detectable branches (most detectable with HCO⁺). Therefore, even though the dendrograms in Figures 7–9 look very different, a comparison of their tree statistics using a uniform noise-level cut, along with an understanding of what each molecule is tracing, produces a consistent picture of the hierarchical structure of dense gas in L1451 from ∼0.5 pc to sub-0.1 pc scales.

The tree statistics of subregions from individual molecular tracers are reported in the second section of Table 7. We compare the subregion statistics of Features A, B, C, and all other features with hierarchical complexity (e.g., Feature H in the HCO⁺ and HCN dendrograms). The sections of the dendrograms corresponding to the individual features are marked in Figures 7–9 with letter identifiers, and we only consider structures above those identifiers in this comparison. For example, in the HCO⁺ dendrogram, Features A and C merge together at the branch labeled "A+C," but we consider the statistics of the individual features only above labels "A" and "C." Features A, B, and C are the only features with any hierarchical complexity in the N₂H⁺ dendrogram, and are the ones with the most branching steps in the HCN and HCO⁺ dendrograms. They are also the most likely sites for current and future star formation, because they account for emission surrounding the only continuum source detections in the field: Feature B surrounds L1451-mm and Features A and C surround Bolocam 1 mm sources.

There is a trend of decreasing maximum branching level and mean path length from Feature B to A to C, and then to the other remaining features. Features A and B are more similar to each other than either is to the remaining other features, indicating that the gas in both features has fragmented a similar amount relative to the rest of the complex structure in the L1451 field. The mean path length and maximum branching level of Feature C bridge the gap between the maximum fragmentation amount seen in Feature A and B and the minimum fragmentation amount seen in the remaining other features.

We interpret the similarity in hierarchical branching levels between Features A and B to mean that these subregions have progressed to a similar stage along the evolutionary track of cloud fragmentation. We know that a young star or first core (L1451-mm) is forming within Feature B at or near the location of the maximum gas emission intensity (leaf 66, 30, and 10 for HCO⁺, HCN, and N₂H⁺, respectively). For Feature A, PerBolo 6 is at or near the location of maximum gas emission intensity in Feature A (leaf 15, 7, and 6 for HCO⁺, HCN, and N₂H⁺, respectively). We argue that with Feature A and B showing very similar tree statistics, Feature B having a confirmed compact continuum detection at its hierarchical peak, and Feature A having a single-dish continuum detection at its hierarchical peak, that a star is likely to form within Feature A. This argument can be extended to Feature C, making it the next most likely place for current or future star formation, followed by the even less fragmented features. Follow-up observations of the single-dish cores and other column density enhancements in these features will be useful for testing this expectation. The mean branching ratios between all features are similar for all molecules, indicating that all structures are fragmenting into a similar number of substructures at each branching step, regardless of how far a feature is along its evolution toward forming stars.

We used Herschel  160, 250, 350, and 500 μm observations of L1451 to derive the column density and temperature of the dust. There is no detectable emission at 70 μm toward this region. The Herschel images were corrected for the zero-level offset based on a comparison with Planck and DIRBE/IRAS data (Meisner & Finkbeiner 2015). In detail, we added offsets of 37.1, 42.6, 23.3, and 10.2 MJy sr⁻¹ to the Herschel 160, 250, 350, and 500 μm maps, respectively. The images were convolved to the angular resolution of the Herschel 500 μm band (∼36'') using the convolution kernels from Aniano et al. (2011), and were regridded to a common pixel size of 10''. The fitting was performed on a pixel-by-pixel basis with a modified blackbody spectrum of ${I}_{\nu }={\kappa }_{\nu }B(\nu ,T){\rm{\Sigma }}$, where κ_ν is the dust-mass opacity coefficient at frequency ν, B(ν, T) is the Planck function at temperature T, and ${\rm{\Sigma }}=\mu {m}_{p}N({{\rm{H}}}_{2})$ is the gas-mass column density for a mean molecular weight of μ = 2.8 (e.g., Kauffmann et al. 2008) assuming a gas-to-dust ratio of 100:1. We assume ${\kappa }_{\nu }=0.1\times {(\nu /{10}^{12}\mathrm{Hz})}^{\beta }$ cm² g⁻¹ (Beckwith & Sargent 1991) and β = 1.7.

The resulting column density and temperature maps are shown in Figures 13 and 14, respectively. We validated our spectral energy distribution (SED) fitting procedure by comparing our derived optical depth to that of the Planck Collaboration et al. (2014) model. Specifically, we re-ran our SED fits after smoothing the Herschel input maps to match the Planck Collaboration et al. (2014) resolution of 5', and found that our derived 300 μm optical depth agrees with that of the Planck-based model to within 5%, on average. The mean column density of L1451 that is enclosed within the 2.0 × 10²¹ cm⁻² contour in Figure 13 (the white contour that encircles all of the high column density regions) is 3.7 × 10²¹ cm⁻², with a standard deviation of 1.7 × 10²¹ cm⁻². The peak column density of 1.2 × 10²² cm⁻² occurs at the location of the Bolocam source, PerBolo 4. The mean temperature within the 15.0 K contour of Figure 14 is 14.0 K, with a standard deviation of 0.7 K. A minimum temperature of 11.9 K occurs at the location of L1451-mm.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Gas kinetic temperature measurements exist toward the four Bolocam sources in our field (Rosolowsky et al. 2008a). Those gas temperatures are ∼2–3.5 K lower than what we find by fitting the Herschel SEDs. Our results, and those from Planck Collaboration et al. (2014) that we compared to, assume emission from a single cloud layer. However, variations in temperature along the line of sight, typically caused by a warmer layer of foreground and background material surrounding a dense, cold star-forming region, are often present. This warmer cloud component can drive the fitted temperature up and the fitted column density down when only doing a single component fit.

To estimate the systematic overestimation of temperatures and underestimation of column densities toward the densest regions of clouds, we created a simple radiative transfer model with a cold layer at 9 K (representing L1451) between two warmer layers at 17 K (representing the foreground and background cloud material). If the warm layers have an A_V of 0.4, then the temperature in the cold layer regions with A_V ≳ 2 is overestimated by ∼2.5–5.5 K when using a single component fit (e.g., measure 12 K instead of the true 9 K). This matches the differences between our temperatures and the temperatures from the ammonia data. Furthermore, the measured column density values of a cold, dense region are predicted to be a factor of two lower than the true values in regions where the warm layers have A_V of 0.4 and the cold layer has A_V  ≳ 2. If the lower-resolution Herschel beam is not filled in regions of cold, dense gas, then this could further bias temperatures to be too warm and column densities to be too low. In Section 8, we discuss how these systematic uncertainties can affect our virial analysis results. An upcoming paper will present a detailed comparison of single-layer and multi-layer SED fitting across Perseus, quantifying the improvement that can be achieved for cold, dense regions when considering the hot, diffuse component along the line of sight.

The column density results in the previous section are angular resolution limited compared to our CLASSy maps, so it is difficult to estimate the mass within the smallest molecular structures we identified using the dendrogram analysis in Section 6. Therefore, we take the approach of first defining structures based on the dust data here, and then using the kinematic information within those structures to explore their energy balance in the next section.

We converted the ${N}_{{{\rm{H}}}_{2}}$ column density map to an extinction map using a conversion factor of ${N}_{{{\rm{H}}}_{2}}$/A_V = (1/2) × 1.87 × 10²¹ cm⁻² mag⁻¹ (Draine 2003). We then ran our non-binary dendrogram algorithm on the extinction map to define dust structures in the field. We used an rms and branching step of 0.2 A_V, and required local maxima to peak at least 0.4 A_V above the merge level to be considered a real leaf. The algorithm identified eight leaves and six branches in the region where we have molecular line data.

The extinction map, with dendrogram-identified structures, is shown in Figure 15. Properties of the dendrogram structures are listed in Table 8, including their coordinate, size, weighted mean temperature and column density, and total mass. The mean temperature and column density and total mass of each structure considers all of the emission interior to the structure (e.g., branch 9 includes emission from branch 9 and leaves 0 and 2). To calculate the total mass within each structure, we first converted the column density at each pixel to a solar mass unit:

$$\begin{eqnarray}&&M(\alpha ,\delta )={\mu }_{{{\rm{H}}}_{2}}{m}_{{\rm{H}}}{N}_{{{\rm{H}}}_{2}}(\alpha ,\delta )A,\end{eqnarray} \tag{ 1 }$$

where ${\mu }_{{{\rm{H}}}_{2}}$ is the molecular weight per hydrogen molecule (2.8), m_H is the mass of a hydrogen atom, ${N}_{{{\rm{H}}}_{2}}(\alpha ,\delta )$ is the column density at a pixel location, and A is the pixel area. As before, the assumed distance is 235 pc. We then totaled the mass enclosed within each structure. The mass of leaves ranges from 0.2 to 5.1 M_☉, and the mass enclosed within branches ranges from 3.0 to 15.9 M_☉. The mass interior to the yellow and purple contours in Figure 15 is 14.9 M_☉ and 12.4 M_☉, respectively.

Zoom In Zoom Out Reset image size

Table 8.  Dust Structure Properties

No.	R.A.	Decl.	Size	$\langle T\rangle$	$\langle {N}_{{{\rm{H}}}_{2}}\rangle$	Mtot	${\sigma }_{\mathrm{tot},{{\rm{N}}}_{2}{{\rm{H}}}^{+}}$	${\sigma }_{\mathrm{tot},\mathrm{HCN}}$	${\sigma }_{\mathrm{tot},{\mathrm{HCO}}^{+}}$
	(h:m:s)	(°:':'')	(pc)	(K)	(1021 cm−2)	(solar)	(km s−1)	(km s−1)	(km s−1)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)
Leaves
0	03:25:13.3	+30:19:05.8	0.23	13.6	5.1	5.1	0.34	0.39	0.49
1	03:25:00.2	+30:21:22.1	0.08	14.4	2.8	0.2	...	0.41	0.44
2	03:25:27.2	+30:21:43.6	0.17	13.2	5.0	3.1	0.27	0.29	0.40
3	03:24:35.5	+30:22:09.5	0.07	13.6	4.4	0.4	0.37	0.34	0.49
4	03:24:47.0	+30:23:11.2	0.14	13.9	4.0	1.4	...	0.54	0.64
5	03:24:26.8	+30:22:47.7	0.17	13.9	4.6	2.0	0.28	0.44	0.53
6a	03:25:08.7	+30:24:09.5	0.07	12.2	7.2	0.5	0.25	0.25	0.33
7	03:25:01.5	+30:24:27.9	0.07	12.7	6.4	0.5	0.26	0.28	0.40
Branches
8	03:25:01.5	+30:24:12.9	0.21	13.7	4.1	3.0	0.33	0.33	0.42
9	03:25:17.8	+30:20:07.0	0.41	14.2	3.1	13.1	0.44	0.47	0.54
10	03:24:29.2	+30:22:26.8	0.27	14.2	3.6	4.5	0.35	0.45	0.57
11	03:24:54.2	+30:23:38.9	0.33	14.1	3.3	6.2	0.34	0.39	0.51
12	03:24:43.1	+30:23:04.4	0.46	14.3	3.1	12.4	0.39	0.43	0.54
13	03:25:17.1	+30:19:52.8	0.51	14.3	3.0	14.9	0.47	0.50	0.56

Note. (1) Dust structure number as labeled in Figure 15. (2)–(3) Centroid position of dust structure determined from regionprops in MATLAB. (4) Geometric mean of the major and minor axis fit to the structure (used as diameter of structure). (5) Weighted mean temperature within the structure. (6) Weighted mean column density of H₂ within the structure. (7) Total mass within the structure. (8) N₂H⁺ velocity dispersion calculated from the integrated spectrum across the structure; see Section 8.2 for a discussion. (9) HCN velocity dispersion calculated from the integrated spectrum across the structure; see Section 8.2 for a discussion. (10) HCO⁺ velocity dispersion calculated from the integrated spectrum across the structure; see Section 8.2 for a discussion.

^aL1451-mm is located within this structure.

Download table as:  ASCIITypeset image

L1451 is the one CLASSy region with strong, widespread HCO⁺ and HCN that is not affected by outflows and significant self-absorption. This enabled a dendrogram analysis of all three molecules in Section 6, instead of just N₂H⁺, and lets us compare the identified molecular structures to the column density structure of L1451 presented in Section 7. Figures 16 and 17 have the dendrogram footprints of the lowest-level branches and all the leaves, respectively, overplotted on the column density map derived from Herschel data. We combined the footprints of all the contoured emission in Figures 16 and 17 to create a mask for the column density map to determine how well the molecular emission captures material at different column densities.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 18 shows the cumulative distribution functions for column density in regions with and without line emission (regions with line emission are defined as having at least a single molecular detection, while regions without line emission are defined as having no molecular detections). A threshold for star formation above A_V ∼ 8, or an H₂ column density of 7.5 × 10²¹ cm⁻², has been postulated based on the distribution of prestellar cores and protostars within the densest regions of molecular clouds (Johnstone et al. 2004; Hatchell et al. 2005; André et al. 2014 and references therein). Only 0.002% of dust regions without a molecular gas detection are above the threshold column density if we take our derived column densities as correct. If our measured column densities are uniformly underestimated by a factor of two from the true values, as was discussed in Section 7.1, still only 1% of dust regions without a molecular gas detection are above the threshold column density. Of the regions where we detect molecules, 90% are at column densities above 1.9 × 10²¹ cm⁻², with a maximum of 1.3 × 10²² cm⁻² and minimum of 8.9 × 10²⁰ cm⁻². Of the regions where we do not detect molecules, 90% are below 2.4 × 10²¹ cm⁻², with a maximum of 7.8 × 10²¹ cm⁻² and minimum of 3.6 × 10²⁰ cm⁻². This shows that spectral line observations using our suite of dense-gas tracer molecules are a great probe of the star-forming material in young regions above the threshold for star formation and down to column densities of a few ×10²¹ cm⁻².

Zoom In Zoom Out Reset image size

This result of molecular emission capturing most of the cloud material near and above the threshold of star formation, combined with the result that the branches in Figure 16 are fragmenting to form the leaves in Figure 17 in a similar hierarchical fashion for each molecule (see Section 6.3), suggests that the dendrogram-identified molecular structure is able to trace physical structure formation despite some biases due to chemistry and extinction. This provides observational evidence that structure formation precedes star formation in molecular clouds. Complex morphological structure in a turbulent cloud can be formed by turbulence-driven cascades from large-scale flows before the onset of star formation (Klessen & Hennebelle 2010, and reference therein), and it can be produced from internally driven turbulence from protostellar feedback (Carroll et al. 2009; Federrath et al. 2014; Nakamura & Li 2014). Our result helps to disentangle externally and internally driven structure which is important for demonstrating an observational case of complex, hierarchical dense-gas structure existing in a turbulent cloud at an epoch before internal feedback can affect the natal cloud environment.

We next use a virial analysis to assess the stability of L1451 structures against collapse. We know one star is currently forming in L1451, based on the detection of compact 3 mm emission at L1451-mm. However, it is unknown whether the majority of L1451 structures are dominated by gravitational potential energy and therefore unstable against collapse, or if they are dominated by internal pressure and therefore stable against collapse and possibly confined by large external pressure. For the virial analysis, we use the dust results in combination with the molecular data. In Section 7, we derived column densities and temperatures across L1451 at the angular resolution of the longest wavelength Herschel band (36''). We also found dendrogram-identified dust structures in the extinction map. In this section, we use that information for dust leaves along with CLASSy kinematic data from Section 5 to assess the energy balance of several structures in L1451.

The virial theorem is useful for describing the stability of structures. There are several approaches for applying the virial theorem to molecular cloud observations in the literature (e.g., Larson 1981; Bertoldi & McKee 1992; McKee & Zweibel 1992; Ballesteros-Paredes 2006; Kauffmann et al. 2013). We follow the formalism of Bertoldi & McKee (1992) and Kauffmann et al. (2013), who define the observed virial parameter of individual structures as

$$\begin{eqnarray}&&{\alpha }_{\mathrm{obs}}\equiv \displaystyle \frac{5{\sigma }_{\mathrm{tot}}^{2}R}{{GM}}\end{eqnarray} \tag{ 2 }$$

where σ_tot, R, and M are the one-dimensional velocity dispersion (including thermal and nonthermal gas motions), radius, and mass of the structure, respectively. Kauffmann et al. (2013) summarize that the critical value of the virial parameter for a spherical structure not supported by magnetic pressure is α ≈ 2. Structures with α ≪ 2 are supercritical and unstable against collapse unless they are supported by magnetic fields. Structures with α ≫ 2 are subcritical and stable against collapse—they should dissipate unless they are confined by high external pressure. As discussed in Kauffmann et al. (2013), this critical value includes the effects of density gradients within and outside the structure, as well as surface pressure from the immediate surrounding medium.

For the virial parameter definition in Equation (2), M is the mass of each source derived from our Herschel analysis and is reported in Column 7 of Table 8. A major source of systematic uncertainty in the calculated mass is the underestimated column density value derived from the Herschel SED fitting in Section 7.1. In Section 7.1, we discussed how using a single-temperature model (ignoring the warmer foreground and background components surrounding the colder L1451 region) for fitting the Herschel SEDs leads to estimated column densities that are about a factor of two lower than the true column densities. If the measured column densities are a factor of two lower than the true column densities, this directly leads to a mass underestimate by a factor of two and a virial parameter overestimate by a factor of two. Work has been done on the galactic scale to estimate dust-mass underestimation due to poor spatial resolution and temperature mixing in Herschel data. For example, Galliano et al. (2011) found that Herschel-derived dust masses can be ∼50% underestimated in the Large Magellanic Cloud (LMC). Even though our observations have much higher spatial resolution than observations of the LMC, this is a useful comparison because the warm-cold temperature contrast in local clouds like L1451 may be more extreme than on galactic scales. For core scales, Shetty et al. (2009) demonstrated that temperature variations along the line of sight will lead to inaccurate measurements of dust properties from SED fitting. We consider the likely systematic underestimation of mass from our SED fitting in the discussion of the results below.

The radius of each structure, R, in the virial parameter definition is half of the size reported in Column 4 of Table 8. The one-dimensional velocity dispersion of each structure, σ_tot, in the virial parameter definition includes thermal and nonthermal support against collapse. The thermal component of the velocity dispersion was calculated as:

$$\begin{eqnarray}&&{\sigma }_{\mathrm{th}}=\sqrt{\displaystyle \frac{{kT}}{\mu {m}_{{\rm{H}}}}},\end{eqnarray} \tag{ 3 }$$

where T is the mean temperature of the dust within the structure, μ is the mean molecular weight (2.33), and m_H is the hydrogen atomic mass. The nonthermal component of the velocity dispersion was calculated for each molecular tracer separately, as:

$$\begin{eqnarray}&&{\sigma }_{\mathrm{nth},\mathrm{tracer}}=\sqrt{{\sigma }_{\mathrm{obs},\mathrm{tracer}}^{2}-{\sigma }_{\mathrm{th},\mathrm{tracer}}^{2}}\end{eqnarray} \tag{ 4 }$$

where ${\sigma }_{\mathrm{obs},\mathrm{tracer}}^{2}$ is the observed velocity dispersion of the molecular emission within the structure, and ${\sigma }_{\mathrm{th},\mathrm{tracer}}$ is the thermal velocity dispersion of the molecule at the mean temperature found within the structure:

$$\begin{eqnarray}&&{\sigma }_{\mathrm{th},\mathrm{tracer}}=\sqrt{\displaystyle \frac{{kT}}{{\mu }_{\mathrm{tracer}}{m}_{{\rm{H}}}}}.\end{eqnarray} \tag{ 5 }$$

We determined ${\sigma }_{\mathrm{obs},\mathrm{tracer}}^{2}$ by masking the CLASSy molecular data cubes using the boundaries of each dust structure, generating an integrated spectrum for each molecule, and fitting for the velocity dispersion of that integrated spectrum. The final ${\sigma }_{\mathrm{total},\mathrm{tracer}}$ value was then calculated as:

$$\begin{eqnarray}&&{\sigma }_{\mathrm{tot},\mathrm{tracer}}=\sqrt{{\sigma }_{\mathrm{th}}^{2}+{\sigma }_{\mathrm{nth},\mathrm{tracer}}^{2}}\end{eqnarray} \tag{ 6 }$$

and is reported in Columns 8, 9, and 10 of Table 8 for N₂H⁺, HCN, and HCO⁺, respectively. The temperature overestimate from the single-temperature Herschel SED fitting discussed in Section 7.1 systematically increases the thermal component of the internal velocity dispersion of each structure, thereby increasing the dispersive term used to calculate α. The temperature error is less significant than the factor of two column density error discussed above, but still contributes a ∼7% systematic increase in the virial parameter if the measured temperature is overestimated by 2.5 K from the true value.

Virial parameters for each dust leaf structure are reported in Table 9, with unique values for each molecular tracer resulting from different velocity dispersions between the tracers. For each tracer we report the measured virial parameter using the masses derived from the Herschel SED fitting, and then report a value reduced by a factor of two to account for the likely systematic underestimation of the fitted column density values discussed above.

Two things are clear from the α values in Table 9. First, the α values typically vary by a factor of two between tracers, with N₂H⁺ and HCO⁺ data giving the lowest and highest virial parameters, respectively. This trend with the molecular tracer is not surprising, because it connects to the differences in spatial extent and kinematics between molecular emission discussed in Sections 4 and 5. HCO⁺ traces the most spatially extended emission on the plane of the sky, implying that it is tracing the most extended emission along the line of sight if the two spatial directions are correlated. This implication is strongly supported by HCO⁺ having the largest line-of-sight velocity dispersions. Because the virial parameter scales as the square of the velocity dispersion of each molecule, the larger HCO⁺ velocity dispersion leads to a larger virial parameter. This shows that the choice of the kinematic tracer has large impact on the interpretation of a given structure.

Second, all virial parameter values are greater than two if we use the masses derived from the single-temperature model of Herschel SED fitting. If we assume that the mass is underestimated by a factor of two and the virial parameter is overestimated by a factor of two (ignoring the systematic temperature uncertainty also increasing the virial parameter), then some of the structures have virial parameters near or less than the critical value of 2 for instability against collapse. Structures 1, 3, 4, and 7 have an α above two regardless of the molecular tracer used or assumption of mass. These structures have the weakest dense-gas emission in the field and/or appear the least centrally condensed in the column density maps, so it is reasonable that they appear to be the least unstable against collapse. Structures 0, 2, 5, and 6 have an α near or below two, appear centrally condensed in the column density map, and have strong molecular emission in the line maps. Structure 0 contains Per Bolo 4, structure 2 contains Per Bolo 6, structure 5 contains L1451-west, and structure 6 contains L1451-mm. These results support these structures (or at least the central parts of these structures) being unstable against collapse, and the L1451 region as a whole being at the onset of forming a few protostars.

We will extend this type of virial analysis in an upcoming work to compare virial parameters across CLASSy regions. This way, even if there are systematic and statistical uncertainties, all structures will have been observed and analyzed in the same way, making a relative comparison of virial parameters across cloud regions at different stages of evolution possible.

In this section we explore the properties of the two centrally condensed, roughly spherical cores in L1451. One is L1451-mm, which we know is a compact object containing a YSO or a FHSC (Pineda et al. 2011). The other is L1451-west, which was discovered using these observations. We summarize the morphology and kinematics of L1451-mm below, and then describe the properties of L1451-west and how they compare with L1451-mm.

The top row of Figure 19 shows the molecular and continuum features of L1451-mm, which is the only confirmed compact continuum core in L1451; all molecules trace strong emission near and surrounding its location. All molecules peak in integrated emission at a location slightly offset from the L1451-mm compact continuum core. The molecular emission that surrounds the core is more concentrated for N₂H⁺, and more widespread for HCN and HCO⁺. The N₂H⁺ centroid velocity field shows a gradient across the core, which Pineda et al. (2011) modeled as a rotating and infalling envelope using slightly higher angular resolution (5'') data; we measure a gradient of about 6 km s⁻¹ pc⁻¹ through the peak of integrated emission at a position angle of 120° east of north (see the white line in the top row panels of Figure 19). The velocity dispersion field is narrowest around the edges of the core with minima ∼0.07 km s⁻¹, and it shows an increase toward the core location with a peak ∼0.2 km s⁻¹. See Pineda et al. (2011) for detailed discussion and modeling of the possibility that L1451-mm is a dense core with a central YSO and disk, or a dense core with a central FHSC.

Zoom In Zoom Out Reset image size

The bottom row of Figure 19 shows molecular and continuum features of L1451-west. The N₂H⁺ integrated intensity map does not have a single peak of integrated emission at the center of the structure. Instead, there are two peaks in the southern half of the source ∼2.5 Jy beam⁻¹ km s⁻¹ and one peak in the northern half of the source ∼2.4 Jy beam⁻¹ km s⁻¹ (the red contours in the bottom-left panel of Figure 19 represent 1.1–2.5 Jy beam⁻¹ km s⁻¹ in 0.2 Jy beam⁻¹ km s⁻¹ increments). The source is considerably weaker in HCN, and not detected in HCO⁺. This apparent lower-abundance of HCO⁺ may indicate that this is a very cold, dense region of L1451, where CO is depleted. If CO is frozen out onto dust grains, then it is not able to destroy N₂H⁺, and N₂H⁺ remains abundant. Likewise, if CO is frozen out, it is not available to create HCO⁺ through collisions with ${{\rm{H}}}_{3}^{+}$ (Prasad & Huntress 1980). Our temperature map in Figure 14 does not suggest that dust around L1451-west is significantly colder than L1451-mm, but higher-resolution continuum observations with modeling to remove any warm component could reveal differences in the two sources that are not resolvable with the current data. It is also possible that the HCO⁺ ($J=1\to 0$) photons are being absorbed by HCO⁺ molecules in the lower-density foreground gas; observations of higher-J transitions are needed to test this.

The N₂H⁺ centroid velocity field of L1451-west shows a relatively complex structure compared with L1451-mm, with higher-velocity emission in the southern and northeastern sections of the core, and lower-velocity emission toward the center and west of the core. The velocity dispersion field of L1451-west is narrowest around the edges of the core with minima ∼0.07 km s⁻¹, and it peaks near ∼0.2 km s⁻¹ at the location of the lowest-velocity emission in the central part of the source. It is possible that infall is producing the increased velocity dispersion, in addition to any undetected outflows and/or unresolved disk rotation from source center. Infall would broaden the molecular emission toward the core center, with the blueshifted emission being brighter than the redshifted emission if the N₂H⁺ emission is optically thick. This could explain why the peak of velocity dispersion is at the location with the bluest centroid velocity. The signal-to-noise of the detectable HCN spectra are too low to look for evidence of infall motions, so this needs to be followed up with deeper observations of optically thick and thin lines.

There is no compact continuum emission detected toward L1451-west. The 3σ flux density limit for a point source (≲3'' = 700 au) in our observations is 3.9 mJy, corresponding to a mass limit of 0.08 M_☉ for the conversion from mJy to M_☉ discussed in Paper I. We use the location of peak velocity dispersion as a proxy for the location of a compact continuum core, if it exists; this peak is offset from the peak in the Herschel-derived column density map by about 16''.

We estimated the size and physical density of L1451-mm and L1451-west for a more detailed comparison. To determine the size, we used the MIRIAD imfit routine to fit a two-dimensional Gaussian to the N₂H⁺ integrated intensity map. The geometric means of the major and minor axes are 37'' and 32'' for L1451-mm and L1451-west, respectively. To determine the maximum physical density, we took the peak column density of each source within the fitted Gaussian, and assumed their depth was the same extent as their geometric mean across the plane of the sky. For both sources, the physical density is 7–8 × 10⁴ cm⁻³ in the 8'' beam. This density is sensible, considering that we observe N₂H⁺ and HCN ($J=1\to 0$) in L1451, which are molecular transitions with effective excitation densities on the order of 10⁵ cm⁻³.

We also compared the mass within 4200 au (18'' at d = 235 pc) to the intrinsic radius at 70% of peak intensity, following the analysis of starless cores presented in Kauffmann et al. (2008). We measured the mass using the Herschel data, and got 0.22 and 0.23 M_☉ for L1451-mm and L1451-west, respectively. We measured the radius using the N₂H⁺ integrated intensity, because the dust data is not as well resolved to derive accurate radius measurements. L1451-mm measured 8'' and L1451-west measured 12'' at 70% of the peak intensity. Compared to the sample of starless cores in Kauffmann et al. (2008), L1451-mm is more centrally dense than most starless objects, while L1451-west is near the upper limit of central density seen in starless objects.

All of these observational results point to L1451-west being similar to, but slightly less evolved than L1451-mm, yet more evolved than a prestellar core. Although deep continuum and spectral line observations will be needed to determine the true nature of these sources, they are some of the best evidence that L1451 is a region that is just starting to form its first stars.

Comparing the projected size of cloud structure to its depth along the line of sight gives a better understanding of the three-dimensional geometry of a region where stars are forming (e.g., a region that is primarily planar/sheet-like versus one that is more spherical). We described a statistical method to estimate the typical line-of-sight depth of cloud structures in Paper I. The method uses the spatial and kinematic properties of dendrogram objects presented in Tables 4–6. It assumes that the V_lsr variation (ΔV_lsr; Column 11) of a structure scales with its projected size (Column 9) in a turbulent medium. It also assumes that the mean nonthermal velocity dispersion ($\langle \sigma {\rangle }_{\mathrm{nt}}$) of a dendrogram structure scales with its depth along the line-of-sight in a turbulent medium. We calculate $\langle \sigma {\rangle }_{\mathrm{nt}}$ for each structure in all three molecular tracers by subtracting the thermal velocity dispersion of 10 K gas of the given tracer from the value reported in Column 12 of Tables 4–6: $\langle \sigma {\rangle }_{\mathrm{nt}}=\sqrt{\langle \sigma {\rangle }^{2}-{\sigma }_{\mathrm{th}}^{2}}$. With those assumptions, we create two size-linewidth relations using all the dendrogram structures (one with projected size versus $\langle \sigma {\rangle }_{\mathrm{nt}}$, and the other with projected size versus ΔV_lsr), and take the spatial scale where they cross as the typical depth of the region. See Paper I for the theoretical framework and numerical results that justify this method.

We have used the method in previous papers to argue that the typical depth of the N₂H⁺ emission in the CLASSy Barnard 1 and Serpens Main regions was 0.1–0.2 pc (see Paper I and Lee et al. 2014). Figure 20 shows the size-linewidth relations for each molecule. As in the analysis of Barnard 1 in Paper I and Serpens Main in Lee et al. (2014), the projected size-ΔV_lsr for each molecule has a positive slope (which is expected for a turbulent medium), while the projected size-$\langle \sigma {\rangle }_{\mathrm{nt}}$ relationship has a flatter slope (which is expected if all structures have a similar depth along the line of sight, independent of their projected size). The best-fit lines to the relations cross at a size-scale of 0.11 and 0.10 pc for N₂H⁺ and HCN, respectively, indicating that those molecules are tracing structures ∼0.1 pc in depth. The HCO⁺ fits cross near 0.40 pc, indicating that it is tracing structures with larger line-of-sight depths compared to the other molecules. This is consistent with HCO⁺ detecting larger-scale, lower density gas than the other molecules.

Zoom In Zoom Out Reset image size

These results show that the N₂H⁺ emission in L1451 has a similar depth as in Barnard 1 and Serpens Main. However, there are differences between L1451 and the other regions that lead us to different conclusions about the large-scale structure of L1451. In Paper I, we said that the depth of Barnard 1 (0.1–0.2 pc) is comparable to the largest size of individual N₂H⁺ dendrogram structures identified from the isolated N₂H⁺ hyperfine component (0.2–0.3 pc). But we followed up by discussing how the Barnard 1 has contiguous N₂H⁺ structure at parsec scales when considering the full N₂H⁺ emission instead of just the isolated hyperfine component. We used the estimated 0.1–0.2 pc depth and observed parsec projected size of Barnard 1 to conclude that the region is flattened at the largest detectable scales.

In projection, the full N₂H⁺ emission, submillimeter continuum emission, and other molecular emission from the L1451 region does not appear to have a contiguous structure across parsec scales like all other CLASSy regions. Instead, the emission is concentrated in a few major features, subparsec in size, that were identified in Section 4. Thus, we argue that L1451 is not a flattened, sheet-like region of dense gas and dust at parsec scales with connected substructure. In N₂H⁺, it appears more like a loose collection of dense concentrations that are ∼0.2 pc projected on the sky and ∼0.1 pc deep.

The physical density of the regions of L1451 with N₂H⁺ emission can be estimated from the cloud depth and column density. For a 0.1 pc depth of N₂H⁺ emission and a mean N(H₂) of 6 × 10²¹ cm⁻² measured from the Herschel data, the derived physical density is 2 × 10⁴ cm⁻³ in the regions of L1451 with N₂H⁺ emission. It is possible to excite N₂H⁺ at a range of physical densities, where the range depends on the molecular column density (and thus relative abundance of N₂H⁺ to H₂), the gas kinetic temperature, and the linewidth of emission. Using the RADEX program (van der Tak et al. 2007), we find that N₂H⁺ ($J=1\to 0$) can reach 1 K brightness temperatures with a physical density of 2.7 × 10⁴ cm⁻³ if the gas kinetic temperature is 12 K, the molecular column density of N₂H⁺ is 6 × 10¹² cm⁻², and the linewidth of emission is 0.3 km s⁻¹. This physical density is essentially the same as derived from our depth and column density measurements, providing a consistency check to our measurements.

We can compare the typical depth of the molecular emission to the projected size of the largest-scale structures of each molecule in Figure 16 to infer the three-dimensional morphology of the individual molecular structures. The projected size of N₂H⁺ structure 21 (the lowest level branch of Features A) is ∼0.17 pc, structure 18 (the lowest level branch of Features C) is ∼0.09 pc, and structure 17 (the lowest level branch of Features B) is ∼0.17 pc. For a typical N₂H⁺ depth of 0.11 pc, the structures have axis ratios of 1.5:1 and 0.8:1, with a mean of 1.3:1. The projected size of HCN structure 45 (the lowest level branch of Feature A) is ∼0.24 pc, structure 43 (the lowest level branch of Feature C) is ∼0.17 pc, and structure 40 (the lowest level branch of Feature B) is ∼0.19 pc. For a typical HCN depth of 0.10 pc, these structures have axis ratios between 2.4:1 and 1.7:1, with a mean axis ratio of 2.0:1. The projected size of HCO⁺ structure 111 (the lowest level branch connecting Features A and C) is ∼0.43 pc, and structure 112 (the lowest level branch connecting Features B, H, F, G, E, and I) is ∼0.46 pc. For our derived typical HCO⁺ depth of 0.40 pc, these structures have axis ratios of 1.1:1 and 1.2:1.

With axis ratios between 0.8:1 and 2.4:1 for individual structures that do not connect at larger scales, these results support the L1451 region being composed of discrete, high-density structures that are approximated as ellipsoids. In addition to not having a contiguous, flattened structure at parsec scales like the other CLASSy regions studied to date (Serpens Main and Barnard 1), none of the large-scale structures of L1451 appear to have clear filamentary substructures like all other CLASSy regions. These results differentiate L1451 from the other regions beyond the lack of protostars and outflows.

If L1451 is not a flattened, sheet-like region at parsec scales, if it does not have prominent filamentary structure, and if it has less active star formation than the other CLASSy regions, then it is natural to speculate what caused these differences. A cloud complex, like Perseus, that spans tens of parsecs will have different turbulent energies in different parts of the cloud. Simulations of turbulence-driven star formation with supersonic turbulence (e.g., Federrath & Klessen 2012) show that star formation within several-parsec scale clouds is clustered in regions where material has been compressed to high densities, leaving voids of star formation in other parts of the cloud. We speculate that the L1451 region of Perseus may not have been compressed as strongly by supersonic turbulence to form an overdense sheet-like structure at parsec scales like may have happened several parsecs away near the cloud regions that became Barnard 1 and NGC 1333. Without a comparable push to high-density that other more active regions of Perseus got, the L1451 region may have been predisposed to forming fewer stars.