Free Access
Volume 592, August 2016
Article Number A56
Number of page(s) 23
Section Interstellar and circumstellar matter
Published online 25 July 2016

© ESO, 2016

1. Introduction

Multiplicity is common in stars: 46% of the solar-type field stars (Raghavan et al. 2010) and more than 82% of the O- and B-type stars (Chini et al. 2012) are multiple stars. Multiple stars are responsible for some of the more interesting phenomena in evolved stars, for example in the dust and gas shells of evolved stars (Maercker et al. 2012; Decin et al. 2015), phenomena such as type Ia supernovae (SNe), blue stragglers and cataclysmic variables that are generated through mass transfer in close binaries. Multiples are also laboratories in which to test models of stellar physics and the products of star formation (Duchêne & Kraus 2013).

Chen et al. (2013) and Tobin et al. (2016) found that the frequency of multiplicity is highest for deeply embedded protostars and decreases to pre-main sequence and field stars in the separation range of 15 to 10 000 AU. These authors used Submillimeter Array (SMA) 1.3 mm and 850 μm archival data and Very Large Array (VLA) 8 mm and 1 cm observations, respectively. However, these surveys are incomplete toward small separations (<15 AU for the VLA and <600 AU for the SMA), and the derived frequency should be considered a lower limit. This clearly shows that stars are frequently born as multiple stellar systems.

While it is considered that fragmentation within the parent cloud is the mechanism through which multiples form, it is uncertain at which point in time and on what scale the fragmentation occurs. Models suggest one of two paths: turbulent fragmentation of the core (1600 AU scales, e.g., Offner et al. 2010), or gravitational instability of the protostellar disk (<500 AU scales, e.g., Stamatellos & Whitworth 2009a; Kratter et al. 2010). While some mechanisms are thought to produce coeval systems, turbulence can cause density enhancements that can lead to non-coevality in multiple protostellar systems. Dynamical ejections of close binaires can, on the other hand, yield apparently non-coeval systems.

Early studies at disk-scale (~100 AU) separations found that 15% out of 10 to 20 T Tauri and pre-main sequence binaries are formed of classical and weak-lined T Tauri stars, similar to mixed pairs in young binaries (Duchêne et al. 1999; Hartigan & Kenyon 2003). Classical T Tauri stars are generally considered to be younger and more actively accreting than weak-lined T Tauri stars (Duchêne et al. 1999; Kenyon & Hartmann 1995). Comparison of these binaries with isochrones showed that secondaries tend to be younger than primaries (Hartigan & Kenyon 2003), but it was suggested that this age difference would disappear with flatter isochrones. A larger study of 65 T Tauri stars in Ophiuchus, Taurus and Corona Australis also found classical and weak-lined T Tauri binaries, in agreement with earlier studies, as well as Class I and II binaries (McCabe et al. 2006) through comparison of color in K, L, [N] and 18 μm observations. This study noted that mixed pairs had a tendency of showing disks with low to no accretion, indicating different ages among the components, and supporting inside-out disk evolution.

Kraus & Hillenbrand (2009) studied the 36 known binaries in the Taurus-Auriga region (d ~ 145 pc) with separations >200 AU, known spectral types and flux ratios with the aim to probe the coevality of pre-main sequence binaries. Coevality of the sample of binaries was determined through comparison with a hybrid of two theoretical isochrones to estimate the ages of each component. Kraus & Hillenbrand (2009) found that two-thirds of the pre-main sequence binaries are coeval with a dispersion lower than 1.4 Myr (0.16 dex), with no trend between age and mass or separation, suggesting that coevality is a product of formation. It should be highlighted that only binaries were probed in Kraus & Hillenbrand (2009), which raises the question of whether the coevality frequency is different when higher order multiples are considered. This is related to the dynamic evolution of multiple systems, since higher order multiples tend to disintegrate more readily (Reipurth 2000), and fewer of them survive to main sequence stages (11% higher order multiples in solar-type stars, Raghavan et al. 2010).

While isochrones are considered the best technique to determine ages, age determination is plagued by large uncertainties, bias and the assumptions made to estimate the age, namely the definition of τ = 0 (Soderblom et al. 2014). For embedded systems, determining the age is even more difficult due to the lack of information on the spectral type and stellar luminosity. Using color, mass accretion rates and inner disk holes to determine evolutionary classfication and ages, while useful for T Tauri stars and even for a few Class I protostars (Duchêne et al. 1999; Hartigan & Kenyon 2003; McCabe et al. 2006), becomes difficult for deeply embedded sources, where near-infrared detections are often lacking and accretion can be more variable (Audard et al. 2014). A more viable focus therefore is to probe the relative evolutionary stages of the components of multiple systems. While the age coevality will not be probed, the evolutionary coevality, which sets the conditions for the system’s life, will be probed and can provide insight into the question.

The evolutionary stage of protostars is usually defined by the spectral energy distribution (SED) shape, infrared spectral index αIR, bolometric temperature Tbol and the ratio of submillimeter luminosity Lsubmm to bolometric luminosity Lbol, which reflects the ratio of stellar mass M to envelope mass Menv (Froebrich 2005). The SED peak will tend to move toward the shorter wavelengths as the protostar evolves and disperses its envelope, changing the shape of the SED. It is expected that the parameters derived from the SED also reflect this, for example, Tbol will increase as the protostar sheds its natal cocoon and αIR will decrease as the protostar moves from the embedded phases to Class II. As the envelope is dispersed, the submillimeter luminosity will decrease and therefore Lsubmm/Lbol will also decrease. In-depth studies of some individual embedded multiple protostellar systems suggest non-coevality in embedded systems (such as L1448N A and B: Ciardi et al. 2003; NGC 1333 SVS13: Chen et al. 2009; L1448C: Hirano et al. 2010; VLA1623: Murillo & Lai 2013) based on these criteria. However, the geometry, or in other words, the inclination and outflow cavity of the observed protostellar system, affects the shorter wavelength (70 μm) part of the SED. This in turn affects the derivation of parameters from the SED, some more strongly than others (Whitney et al. 2003; Robitaille et al. 2006; Crapsi et al. 2008), and the evolutionary stage classification (Enoch et al. 2009; Dunham et al. 2014). Studies of modeled protostellar SEDs demonstrate that accurately constraining the inclination of the source provides more accurate estimates of the derived parameters and thus of the evolutionary stage classification (Offner et al. 2012).

The inclination of the protostellar system with respect to our line of sight can be estimated from outflow observations and is derived with more precision from rotationally supported disk structures, if present. Protostellar systems alter their environment as they evolve, clearing out envelope material through widening of the outflow cavity (Arce & Sargent 2006), accretion and concentration of material onto the protoplanetary disk. Hence, the envelope and outflow can further constrain the evolutionary stage of the source through the chemical and physical structure of the envelope and core. As a consequence, to establish the evolutionary stage of a source and eventually the coevality of a multiple protostellar system, the SED, derived properties, inclination and environment must be accounted for.

In this work we present the construction and analysis of the SEDs of all identified protostellar systems in the Perseus molecular cloud, the largest sample of Class 0, I and II. Perseus is the main target of this work because it is a well-studied region whose multiplicity and environment are relatively well known. This provides data over a wide range of wavelengths and resolutions and both continuum and line emission towards most, if not all, of the region.

For this purpose, literature and archival data were used to construct the short (<70 μm) and long (>160 μm) wavelength regimes of the SEDs. Herschel Space Observatory Photodetector Array Camera and Spectrometer (PACS, Poglitsch et al. 2010) photometric maps were used to cover the peak of the SEDs (70, 100 and 160 μm), without which large uncertainties arise in the parameters derived from the SED. The Herschel PACS beamsize limits the range of component separations that can be probed to 7, which at the distance of Perseus (d ~ 235 ± 18 pc, Hirota et al. 2008, 2011) becomes ~1600 AU. Below this angular resolution, the fluxes of multiple systems are difficult to disentangle. Therefore the coevality frequency, system alignment and properties derived in this work, as well as the environment, will provide constraints for fragmentation models only at core scales (1600 AU).

In this paper we present the constructed SEDs of all identified and known systems in Perseus, with special focus on multiple protostellar systems. The SEDs are constructed from literature and Herschel PACS data, compared with canonical SEDs for different stages from Enoch et al. (2009). With this work, we aim to determine the frequency of coevality in multiple protostellar systems to provide constraints for multiple protostar formation scenarios. Section 2 defines the concepts used in this work. The data and sample studied in this paper, as well as the construction of SEDs and derivation of derived properties, are described in Sect. 3. The results, including an analysis of unresolved SEDs, are given in Sect. 4. Finally, the discussion and conclusions are given in Sects. 5 and 6, respectively.

thumbnail Fig. 1

Cartoon of the definitions used in this work. More evolved sources are represented by larger disks, wider outflow cavities and less envelope material.

2. Definitions

For consistency and clarity of the terms used throughout this work, definitions of certain terms are provided in this subsection. An illustration of these definitions is shown in Fig. 1.

Protostellar system is defined as a source and its surrounding environment composed of a disk, an envelope and a bipolar outflow.

Multiplicity or multiple is used to refer to a system consisting of two or more components or sources, regardless of whether they are stars or protostars. The terms binary, triples and higher order multiples are thus implicitly merged into this term.

Multiple protostellar system or simply multiple system here refers to two or more protostellar sources composing one system. Multiple systems are generally observed to share a common envelope and, in some cases, a common disk. We assume that a group of protostars is gravitationally bound, unless there is evidence to the contrary, if they were observed to have a common envelope in single-dish observations. An observed group of protostellar systems is considered a multiple when several observations and studies confirm its multiplicity through both continuum and molecular line emission.

Coevality is taken to mean the relative evolutionary stages of the components that make up a multiple protostellar system, accounting for the SED, derived properties, inclination and environment. Environment is taken here to mean the outflows, surrounding envelope, and disk(s), if any. A multiple system whose components show similar evolutionary stages is considered coeval, while a multiple system with different evolutionary stages is referred to as non-coeval.

Resolved multiple system is a system with confirmed multiplicity with separations 7 that can be resolved in the Herschel PACS maps.

Unresolved multiple system is a system with confirmed multiplicity down to 0.08 and separations <7 that cannot be resolved in the Herschel PACS maps.

Table 1

Sample of resolved multiple protostellar systems (separation 7).

Table 2

Sample of unresolved multiple protostellar systems (separation <7).

thumbnail Fig. 2

Herschel PACS stamps of resolved multiple protostellar systems in Perseus, except for NGC 1333 IRAS 4 which also shows IRAS 4A, an unresolved protobinary. Each stamp is 80× 80. 70, 100 and 160 μm are shown in blue, green and red, respectively. Blue symbols represent the components of a system, with circles denoting those with additional unresolved multiplicity and diamonds indicating those without (known) additional multiplicity. NGC 1333 IRAS 4A is marked in red.

3. Sample and data

3.1. Source sample

To study the coevality of multiple systems, the component protostellar systems must be identified. Perseus was chosen because of the large number of embedded young stellar objects in a single cloud at d< 300 pc. Our source sample list and coordinates are obtained from Tobin et al. (2016), who identified multiple systems in Perseus down to 15 AU separations in the VLA Disk and Multiplicity survey of Perseus protostars (VANDAM) survey. At the same time, the source sample was divided, based on the findings of Tobin et al. (2016), into three categories: resolved multiple systems and unresolved multiple systems with our 7 separation, and single protostars. The source sample is listed in Tables 1 and 2. Multiple systems and their components are referred to by their most common name. Components with designations PerXX up to 66 are shorthand for Per-emb XX from Enoch et al. (2009). Sources with designations EDJ2009-XXX refer to sources from Evans et al. (2009). Duplicated systems in Tables 1 and 2 arise because some multiple systems have components that have been observed to be close binaries with separations <7 (e.g., NGC 1333 IRAS 7 Per18, Tobin et al. 2016). Half of the wide resolved multiple systems in our sample have a close unresolved companion, that is, 8 out of 16 systems. Of these 8 systems, 3 (L1448N, NGC 1333 IRAS 2 and NGC 1333 IRAS 7) have two resolved components with an unresolved companion each. Confirmed single protostellar systems are listed and discussed in Appendix D.

The sources in our sample have all been confirmed to be protostars through studies at multiple wavelengths, ruling out background sources, AGB stars, or galaxies.

3.2. Literature data

Most star forming regions have been observed at infrared and (sub)millimeter wavelengths at different epochs and with varying resolutions. The first step to constructing SEDs is therefore a search of the available data in the literature. It needs to be noted, however, that even though there is much information in the literature, not all protostellar systems have been homogeneously observed or photometry reported, making it impossible to have all SEDs sampled at the same wavelengths.

The near- to mid-infrared regime of protostellar SEDs is well characterized from 2MASS and Spitzer Space Telescope observations with fluxes shortwards of 70 μm. The c2d catalog (Dunham et al. 2015) provides fluxes from 1.25 μm to 24 μm. The Spitzer 70 μm fluxes are not considered here given the large beam and saturation of the MIPS instrument for the 70 μm detector and the superior quality of the Herschel data. Sensitivity limits at each wavelength from the respective instruments are taken as upper limits for sources that lacked an entry in the c2d catalog. For NGC 1333, integrated fluxes at wavelengths <70 μm were obtained from the compiled catalog of Rebull (2015) after conversion from magnitude to mJy units.

Submillimeter and millimeter integrated fluxes were collected from diverse interferometric continuum surveys (e.g., Looney et al. 2000; Jørgensen et al. 2007; Chen et al. 2013; Yen et al. 2015) as well as works reporting fluxes for individual protostellar systems (e.g., Chen et al. 2009; Hirano et al. 2010; Palau et al. 2014). Careful selection of the fluxes from literature was made to ensure that as much emission could be recovered from the observations as possible, while at the same time the individual sources could be clearly and easily separated. The VANDAM survey (Tobin et al. 2016) provides fluxes from 8 mm to 1 cm for all sources in the Perseus star forming region. Interferometric observations are preferred over single-dish observations because of the resolution needed to separate the flux contribution from each component in a multiple system. The typical fraction of recovered flux varies by telescope configuration, sensitivity and structure being probed. Tobin et al. (2015) provided a comparison that gives an idea of the recovered flux in interferometric observations.

Although data from the literature can cover the near- to mid-infrared and (sub)millimeter regimes of the SED, the peak of the SED is not well sampled typically at 70 to 160 μm. The lack of a well sampled SED peak can seriously underestimate the derived parameters and evolutionary classification of a protostar, and in turn the coevality determination of a system. Herschel PACS data are therefore crucial to this work.

Table 3

StarFinder photometry parameters.

3.3. Herschel PACS photometric maps

Archival photometric maps from Herschel PACS from the Gould Belt Survey (André et al. 2010; Pezzuto et al. 2012) were obtained from the Herschel Science Archive for the entire Perseus region. The maps made with JScanmap were selected for performing the photometry (see Appendix A). From these data we can extract 70 μm Sν 160 μm integrated fluxes. Due to the resolution of Herschel PACS observations, fluxes from each component in a multiple protostellar system can be extracted only for systems whose projected separations are 7.

Star-forming regions tend to be clustered, hence, crowded-field photometry techniques are employed to best exploit the Herschel PACS maps. While aperture photometry is a simple and straightforward method, it is not a viable solution for crowded protostellar fields. Point spread function (PSF) photometry presents a better solution to the problem at hand. The IDL-based program StarFinder (Diolaiti et al. 2000) was employed to perform photometry on the Herschel maps. The PSF was extracted from the maps themselves to account for the specific observation mode, which cannot be achieved as easily with modeled ideal PSFs. Single isolated sources were used to extract the PSF, with moderate brightness, thus avoiding spikes and negative spots, and little to no surrounding nebulosity. The extracted PSFs provide beam sizes of 9.6, 7.2 and 12.8 for 70, 100 and 160 μm, respectively. StarFinder allows deblending of sources and setting a lower limit for the FWHM for source separation, which proves to be very useful in separating multiple systems from PACS maps.

Postcard maps of each sub-region of Perseus, measuring 44′ × 44′, were extracted for ease of photometry. For postcard maps from the same larger map, the same PSF was used, which means that we required the PSF to be extracted only once per map using the best single-source targets. To avoid an overestimation of the measured fluxes and facilitate source deblending, the extracted PSF was then masked by introducing an aperture factor. The StarFinder parameters for the photometry used in this work for each subregion and the typical flux uncertainty per wavelength are listed in Table 3. After PSF photometry was performed, PSF aperture and background corrections were applied to the raw fluxes. The values used for aperture correction are tabulated in Balog et al. (2014).

3.4. SED construction

Given the multiple names and identifiers each source has accumulated through surveys and literature, fluxes at different wavelengths were matched by means of the coordinates with a search radius of 4.5. The search radius was selected to be below the resolution limit of Herschel and similar to the FWHM lower limit for StarFinder, avoiding any confusion in source identification. Coordinates were obtained from (sub)millimeter interferometric observations given the higher angular resolution and because the source positions at these wavelengths are less likely to be contaminated by foreground stars (e.g., NGC 1333 IRAS 2B at λ ≲ 8 μm, Rodríguez et al. 1999) or scattered light.

Table 4

Statistics from the constructed SEDs.

Care was taken that fluxes at all available wavelengths for each SED were separated from the other protostars in their system. At 160 μm this criteria breaks down for systems with separations smaller than 9 to 10. In these cases, the flux is flagged as combined and is noted in the plotted SEDs. Upper limits are also flagged and noted with a different symbol in the plots.

L1448 IRS1, an unresolved multiple systems, has fewer than three points in the SED. The same situation occurs for 7 single protostellar sources, listed in Table D.1. Hence, these systems are not shown in Figs. 4 and D.1.

thumbnail Fig. 3

Constructed SEDs for resolved multiple protostellar systems. Filled circles denote the fluxes without contamination from nearby sources. Triangles indicate upper limits. Squares show combined fluxes. Each SED is overlaid with the template SEDs from Enoch et al. (2009) for comparison (dashed lines).

thumbnail Fig. 4

Constructed SEDs for unresolved multiple protostellar systems. Other details are the same as in Fig. 3.

Table 5

SED derived properties.

3.5. Source properties

Source properties derived from the constructed SED are expected to aid in the evolutionary stage classification. Constraining the peak of an SED improves the calculation of the protostellar system’s derived properties, which makes the Herschel PACS observations crucial for this task. Five parameters were derived for each constructed SED: infrared spectral index, bolometric temperature and luminosity and two luminosity ratios.

For the infrared spectral index αIR, the slope between 3 μm and 24 μm is given by where Fλ is the flux at a given wavelength λ. If the flux at 24 μm is absent, αIR is not reported. When one or more of the fluxes in this range is an upper limit, αIR is a lower limit.

The bolometric temperature Tbol is expressed as where Sν is the flux at a given frequency ν.

The bolometric luminosity Lbol is derived using where D is the distance. Submillimeter luminosity (λ ≥ 350μm) Lsubmm and far-infrared luminosity (λ ≤ 70 μm) Lfir were derived from the same equation using the corresponding wavelength ranges. Both Tbol and Lbol were derived from the SEDs using trapezoidal integration.

In addition, two luminosity ratios were taken: submillimeter to bolometric Lsubmm/Lbol and far-infrared to bolometric Lfir/Lbol. Both ratios were used since interferometric continuum observations resolve out much of the extended flux pertaining to the envelope, while the far-infrared fluxes from Herschel are expected to capture most of the envelope emission. These ratios are meant to reflect the envelope to central star mass ratio (André et al. 1993; Froebrich 2005), which is used to define the separate physical stages of protostars (Robitaille et al. 2006; Enoch et al. 2009). Deeply embedded sources are expected to have luminosity ratios higher than 0.005, while less embedded protostars tend to show ratios lower than 0.005.

3.6. Caveats

The results in this work are limited by the resolution of the Herschel PACS maps. Multiple systems with separations <7 cannot be resolved, making the frequency of non-coevality found in this work applicable to wider systems. Furthermore, the results obtained here can provide constraints for multiple protostellar system formation scenarios at the core scale (1600 AU).

Some multiple systems lack reported resolved submillimeter fluxes, which means that the derived properties are under- or overestimated. This affects the evolutionary classification derived from these parameters. Care must then be taken to consider this aspect when classifying the systems, and the relations between components of a system are more relevant than the actual quantities themselves.

4. Results and analysis

The constructed SEDs are presented in Fig. 3 for resolved systems and in Fig. 4 for unresolved systems. Flux uncertainties are in general comparatively small, hence when plotted, the errors are not much larger than the symbols used for plotting.

The constructed SEDs are analyzed with the aim to study the coevality of multiple systems. All the parameters typically used to identify a protostellar system’s evolutionary stage together with additional diagnostics are used in the classification. This is to ensure that there is as little bias as possible due to inclination, which tends to affect the derived SED parameters. In this section each method and the corresponding results are presented.

4.1. SED shapes

As the protostellar system evolves and clears out the envelope, the peak of the energy distribution shifts to shorter wavelengths. The wavelength at which the SED peaks can therefore be used as an indicator of the evolutionary stage. Average SEDs for the progressive evolutionary stages are shown in Fig. B.1. These SEDs were derived from the Spitzer c2d observations of a large sample of protostars in Perseus and Serpens by Enoch et al. (2009) and divided into classes based on Tbol. Figures 3 and 4 show the constructed SEDs for multiples compared with these average SEDs.

A quick look at the constructed SEDs makes it clear that several multiple systems have components with different SED shapes (e.g., IC 348 Per8+Per55, IC 348 Per32+EDJ2009-366), while others have components with similar SED shapes (e.g., NGC 1333 IRAS 5, IC 348 MMS2) or a combination (e.g., NGC 1333 IRAS 7, L1448 N, NGC 1333 SVS13, B1-b). The similar SEDs hint at coeval components, whereas non-coevality is suggested by the differing SEDs.

To obtain some simple statistics, we counted the systems and identified stages by eye in comparison to each other and to Fig. B.1. The frequency of non-coevality found in this way is listed in Table 4. Higher order multiples were counted twice, once for the first pair and then a second time for the pair compared to the third component. For example, NGC 1333 IRAS 7 was counted once as coeval and once as non-coeval, since Per18 and Per21 appear to have the same evolutionary stage, but are non-coeval relative to Per49. This generates a total of 21 systems where coevality is probed, in contrast to the 16 systems in our sample. We found that 7 of 21 systems (33 ± 10%) show non-ceovality: L1448 N, NGC 1333 SVS13, NGC 1333 IRAS 7, IC 348 Per8+Per55, B1-b, IC 348 Per32+EDJ2009-366 and NGC 1333 Per37. We did not set a maximum separation limit for a multiple system, but it is interesting to see the change in non-coevality frequency in our sample as a limit is set. Assuming the characteristic size of protostellar cores (30), we found that 6 of 15 systems (40 ± 13%) with separations 30 are non-coeval. An arbitrary separation limit of 20 shows 4 out of 14 systems (33 ± 14%) to be non-coeval. This means that the rate of non-coevality does not change significantly by limiting the separation of multiple systems.

For NGC 1333 Per37, the EDJ2009-235 component is not detected in the Herschel PACS maps, but is detected in the c2d and VANDAM surveys. Tobin et al. (2016) classified it as a Class II source, but based on the c2d fluxes, EDJ2009-235 appears to be an embedded source, closer in agreement with the classification from Young et al. (2015). A possible explanation for the discrepancy and its lack of Herschel PACS detection could be a highly extincted disk that might make a Class II source look much younger. The IC 348 systems Per8+Per55 and Per32+EDJ2009-366, which appear as proto-binaries at scales 7, have unresolved components (Table 2), which means they are higher order multiples. This would suggest that higher order multiples tend toward non-coevality.

thumbnail Fig. 5

Summed SEDs from resolved multiple systems, showing how unresolved multiple systems might be composed. The dashed lines are the same as in Fig. B.1.

4.2. Classification from derived properties

All the calculated parameters described in Sect. 3.5 and their errors are listed for each source in Table 5 for the resolved multiple systems. Comparing Tbol within the multiple systems indicates that the rate of non-coevality is much lower than found from the visual comparison of the SEDs. Based on Tbol, 6 multiple systems are found to be non-coeval (NGC 1333 IRAS 7, IC 348 Per8+Per55, B1-b, Per32+EDJ2009-366 and NGC 1333 Per37 twice). Marginal non-coevality, that is, one source being slightly younger than the other (e.g., early Class 0 and late Class 0), can be seen toward 4 systems (L1448 C, NGC 1333 SVS13, NGC 1333 Per6+Per10 and NGC 1333 Per58+Per65), while the remaining are quite coeval. The non-coevality frequency found based on Tbol is between 29 ± 10% and 48 ± 11%, with the latter value considering the marginally non-coeval systems in addition to the non-coeval systems.

Luminosity ratios, such as Lfir/Lbol and Lsubmm/Lbol, are expected to be indicators of the evolutionary stage, with values above 0.005 indicating a Class 0 source. However, for these ratios to provide reliable information, a well-sampled SED at λ ≥ 70 μm is required. Not all the systems we studied here are well sampled in this regime. For most systems the Lsubmm/Lbol ratio is underestimated or cannot be calculated at all, while the Lfir/Lbol can be calculated in all cases but is also underestimated due to the lack of submillimeter flux. Systems relatively well sampled in the submillimeter regime in our sample are NGC 1333 SVS13 and B1-bN & S. The Lsubmm/Lbol ratio for NGC 1333 SVS13 supports the non-coevality of this system, with the ratio for SVS13A being lower than 0.5%, the threshold for embedded sources, while SVS13B and C are (marginally) above the value, indicating they are embedded. In contrast, the Lfir/Lbol ratio shows values well above 0.5% for all three sources, with SVS13A and B having equal values and SVS13C showing a higher value. This would seem to suggest that SVS13C is more embedded than its northern companions, and SVS13A and B are less embedded. B1-bN & S have a Lsubmm/Lbol well above 0.5% and a high value for the Lfir/Lbol ratio, consistent with their deeply embedded condition. A non-coevality frequency cannot be derived from the luminosity ratios in this work because they are not well determined for almost all sources.

The infrared spectral index αIR, which is the slope of the SED between 3 and 24 μm, is assumed to be a good indicator of evolutionary stage, even when geometric effects affect the SED shape (Crapsi et al. 2008). Positive values indicate an embedded source (Class 0 and I), while negative values higher and lower than − 1.5 point toward Class II and III, respectively. NGC 1333 SVS13A and SVS13B lack a point at 24 μm and therefore the value was not calculated. SVS13C has a negative value, suggesting that it is a Class II source. For B1-bN & S, the northern source presents a negative value higher than − 1.5 while the southern source has a positive value, which would indicate that the northern component is Class II and the southern component is an embedded source. The reason for these discrepancies relative to the other parameters is the sensitivity of the spectral index to the fluxes between 3 and 24 μm. For example, NGC 1333 SVS13C’s flux at 24 μm drops below the fluxes at shorter wavelengths, causing the slope and consequently the spectral index to be negative. The same occurs for B1-bN, whose flux at 24 μm is lower than for B1-bS. Interestingly, this parameter causes systems such as NGC 1333 IRAS 7 and IC 348 Per8+Per55 to seem coeval. A coevality frequency cannot be derived from the values of αIR for the sample in this work, since it is not well determined for all sources.

The derived properties suggest the non-coevality frequency to be lower than that obtained from a visual examination of the SEDs, when considering only non-coeval systems and not marginally non-coeval ones. Still, the derived properties are sensitive to inclination and SED sampling, which biases their classification of evolutionary stage. Hence, these parameters should not be taken at face value, but instead be considered together with the SED shape and other properties of the multiple systems.

Table 6

Derived parameters of combined SEDs.

4.3. Resolved versus unresolved SEDs of multiples

The SED of an unresolved multiple system is composed of the sum of SEDs of its individual components. This led us to analyze to which extent an unresolved SED reflects the parameters and evolutionary stage classification of its compoennts. To do so, the SEDs of resolved systems were summed and the resulting shape and derived properties compared to those of the individual components. This analysis does not generate a method to separate unresolved SEDs, but will provide insight into the coevality of close multiple systems. Figure 5 and Table 6 show the results of combining the SEDs of four systems: NGC 1333 IRAS 7 Per18 and Per21, IC 348 Per55 and Per8, B1-bW and B1-bN, and NGC 1333 SVS13 A and B.

From this simple analysis we find that there are mainly three cases. The first is that if the two components have almost identical SEDs, then the combined SED will simply be doubled. This case is shown by NGC 1333 IRAS 7 Per18 and Per21. The second case is when the two components are non-coeval, then the SED will not follow a specific SED shape but will appear odd shaped with two peaks. This is illustrated by the combined SEDs of IC 348 Per8 and Per55, and B1-bW and B1-bN. The final case occurs when one component is noticeably dimmer and younger than the other, then the brightest component dominates the combined SED. NGC 1333 SVS13 A and B illustrate this scenario. Thus, different components can dominate different regions of the SED.

Figure 4 presents the SEDs of unresolved multiple systems. Examples of the second case from Fig. 5, such as L1455 FIR2, B1-a and Per62, are shown in Fig. 4. The first and third cases are next to impossible to identify without additional constraints, and systems such as IRAS 03292, IC 348 SMM22, EDJ2009-156, and HH211 could be examples of either case.

In all scenarios we find that Lbol for the unresolved SED is equal to the sum of both components. In contrast, parameters such as Tbol and Lfir/Lbol are an arithmetic average of the corresponding parameters of the two components. This, of course, affects the evolutionary classification of unresolved multiple systems. While taking the derived values and assuming each component contributes equally may be a good assumption in coeval cases, this could be an over- or underestimation of the true parameters in non-coeval systems.

5. Discussion

From the SED shapes of resolved multiple systems alone, we find a non-coevality frequency of 33 ± 10%. Higher order multiple systems contribute the most to this fraction, with all of five resolved systems (L1448 N, NGC 1333 SVS13, NGC 1333 IRAS 7, B1-b and NGC 1333 Per37) showing indications of non-coevality. The other two systems, IC 348 Per8+Per55 and Per32+EDJ2009-366, that appear as binaries at separations 7, show additional fragmentation at scales <7 in one of the components (Table 2), which also makes them higher order multiples. Binaries, on the other hand, tend toward coevality in our sample. The non-coevality frequency found here is similar to that found by Kraus & Hillenbrand (2009) of one-third in pre-main sequence stars. However, their frequency was found from binaries alone, whereas higher order multiple systems are responsible for this frequency in our study. Kraus & Hillenbrand (2009) also probed down to small separations (>200 AU) an order of magnitude lower than the separations probed here (1600 AU). The unresolved systems studied here that are suspected of non-coevality could account for the rest of that fraction. The question then arises whether the non-coevality frequency obtained here is real or a product of misalignment, that is, the difference in inclination (w.r.t. the line of sight) of each component in a multiple system.

Different SEDs do not necessarily indicate a non-coeval multiple system, but could also be due to geometrical effects, especially if the line of sight lies through the outflow cone. Inclination can alter the shape of the SED and derived parameters (Robitaille et al. 2006; Crapsi et al. 2008; Enoch et al. 2009) and is crucial to attain an accurate determination of the evolutionary stage (Offner et al. 2012). If the protostellar system is seen edge-on, the protostar is obscured and will seem younger. On the other hand, if the protostellar system is observed face-on, the protostar will be unobscured and appear more evolved. Hence, the alignment of multiple systems affects whether the differing SED shapes are a product of real non-coevality of inclination effects. Multiple system formation mechanisms suggest that ordered rotational fragmentation would produce systems with aligned inclination (Burkert & Bodenheimer 1993), whereas turbulent fragmentation is expected to produce random alignment (Offner et al. 2010). Work on pre-main sequence multiple systems shows that binaries have a tendency to be aligned, while higher order multiple systems are less likely to be aligned (Jensen et al. 2004; Monin et al. 2006). However, the numbers seen toward pre-main sequence multiple systems might not reflect the actual alignment at the time of formation because of the dynamical evolution (Reipurth 2000; Jensen et al. 2004). Although a handful of observations exist that show both aligned and misaligned multiple systems at every stage of evolution, there are no statistical numbers on the distribution of inclination and multiple system alignment at the time of formation.

The best method to obtain an accurate inclination estimate is through disk observations and modeling, but not all protostellar systems have confirmed and reported disks. Hence, the inclination of protostellar systems must be constrained through another technique.

5.1. Outflows

Outflows present a viable option, but they can provide only a broad inclination range and may not always be accurate, even more so in multiple systems where precession occurs due to companion perturbations (Fendt & Zinnecker 1998). The evolutionary stage of a protostar is closely linked to its outflow and circumstellar envelope. The envelope is dispersed as the protostar evolves and accretes part of its material (André et al. 1993), while the molecular outflow tends to become weaker with time and the outflow cavity broadens (Velusamy & Langer 1998; Arce & Sargent 2006).

A point to highlight here is that misalignment of outflow axes on the plane of the sky is not the same as misalingment of inclination with respect to the line of sight. The former does not affect the observed SED, while the latter has a strong effect on a protostar’s SED and derived parameters. Thus, we refer to the first case as the outflow position angle (PA) and the second case as inclination misalignment.

Lee et al. (2016) studied the outflows of 9 multiple systems in Perseus and found in all cases that the outflows of wide multiple systems have different PA, that is, that they are misaligned on the plane of the sky. However, strong indications of inclination misalignment were not found.

To compare with the other parameters used to determine the evolutionary stage of components in a system, we explicitly examine here the outflows of a few systems, focusing on signs that are expected to indicate evolutionary stage and inclination misalignment. Specific comments on each resolved multiple system treated in this work are given in Appendix C.

NGC 1333 SVS13A shows molecular outflow lobes that are wide and shell-like, while SVS13C exhibits a collimated outflow with indications of being in the plane of the sky, meaning that the disk is seen edge-on (Plunkett et al. 2013; Lee et al. 2016). This possible inclination misalignment might be the reason that this system appears to be non-coeval based on the SEDs. However, SVS13A has been suggested to be a transition class0/I object based on its association with a Haro-Herbig object and the shell-like morphology of its outflow, making it somewhat older than SVS13B and C.

B1-bN and S have outflows that appear parallel to each other, but the blueshifted lobes are in opposite directions, a tell-tale sign of inclination misalignment. Even though their outflows suggest inclination misaglinment (Gerin et al. 2015), the SEDs appear to be similar, indicating coevality, which is consistent with the results obtained from the analysis of their environment (Hirano & Liu 2014). The inclination misalignment may therefore be small. B1-bW is expected to be older based on the SED, which is why it is not detected in the submillimeter because it has too little envelope material that might be the product of stripping from a neighboring outflow or jet (Hirano & Liu 2014).

This shows that inclination misalignment of the systems does not always generate an apparent non-coevality. To do this, one of the components of a system would need to be significantly inclined, tending toward the line of sight along the outflow cavity so that the source would appear much older.

5.2. Alignment and coevality

To assess the frequency of apparent non-coevality due to misalignment versus real non-coevality, we made a simple statistical test with the Robitaille et al. (2006) SED model grid1. This grid of models offers 10 inclination angles, evenly sampled in cos(i), ranging from 87.13° to 18.19° with 0° being the disk seen face-on (i.e., looking down the outflow cone), which is ideally suited for our task. We constrained the number of models by choosing those for Class 0 sources, filtering by 2.0 MMenvelope ≤ 10.0 M and Menvelope > Mstar. The resulting set of model SEDs, a total of 1037 times 10 inclination angles, have stellar masses of up to 2.0 M and outflow cavity angles ranging from 15 to 30°. The fluxes are obtained for similar wavelengths as our observations, including the Herschel PACS fluxes, and an aperture of 7000 AU (30 at d ~ 235 pc).

Pairs of models with their respective inclinations are randomly drawn from the list of 10 370 Class 0 synthetic SEDs, resulting in 5185 pairs. This is done to simulate, in a simplified manner, multiple systems and compare their SEDs and derived parameters. The random paired models are separated into four groups determined by the difference in inclination angles, that is, their degree of alignment. The four groups are perfectly aligned (Δi = 0°), small misalignment (0° < Δi < 34°), large misalignment (34° < Δi < 69°) and perpendicular (Δi = 69°). Perpendicular alignment is not equal to 90° because of the available inclination angles of the models, but is instead the difference between the edge-on (87.13°) and face-on cases (18.19°). To obtain the frequency of pairs in each group, the times each case occurs were counted based on Δi and are listed in Table 7. Perpendicular alignment is the least likely case (2 ± 0.2%), with small alignment the most common (63 ± 0.6%). Examples of SEDs from each case are shown in Fig. 6.

Table 7

Random pair statistics based on synthetic SEDs.

thumbnail Fig. 6

Randomly paired model SEDs in the four groups based on alignment. The numbers indicate the inclination of the respective component. Note that the SED pairs that appear non-coeval peak at the same wavelength for both components, which is not expected from components at different evolutionary stages.

Apparent non-coevality was checked by first filtering with Tbol, assuming that for apparently non-coeval pairs the Tbol difference is larger than a factor of 3. The reason for using Tbol to filter the models is based on three points: i) Tbol tends to be sensitive to inclination; ii) the thresholds for evolutionary stage classification (late Class 0 to early Class I: 100 K; late Class I to Class II: 650 K); and iii) the non-coeval resolved multiple systems in our sample identified from Tbol have a ratio of about 6 or higher. A factor of 3 was chosen to ensure that Class 0 and I pairs are also included, since for example a component with a Tbol = 50 K may appear non-coeval with a companion having a Tbol > 150 K. The paired SEDs filtered this way were then inspected by eye to subtract the SED pairs that did not appear non-coeval. The frequency of apparent non-coevality due to misalignment is found to be 17% ± 0.5%. Examining the frequency of apparent non-coevality in each of the four groups, we found that large misalignment and perpendicular have the most common occurrences of apparent non-coevality. This is mainly due to one component being face-on or close to face-on, in combination with the outflow cavity opening angle, causing one component to appear older. This was also suggested by the results of the outflows of multiple systems.

A characteristic of SEDs for different evolutionary stages is that the peak of the SED shifts to shorter wavelengths as the envelope is dispersed (see Fig. B.1). While the SEDs can appear non-coeval as a result of inclination effects, the inspection by eye of these SED pairs revealed that the peak around λ ~ 100 μm, characteristic of Class 0 sources, does not significantly shift to shorter wavelengths as a result of inclination, even in the face-on case. Examples of this are shown in Fig. 6. In other words, even though one of the SEDs in the pair appeared more evolved than the companion, this apparently older source retained its peak around λ ~ 100 μm. This is in contrast to the SEDs of non-coeval resolved multiple systems discussed in this work. For example, in NGC 1333 IRAS 7, the peak of Per49 is located around 5.8 μm, whereas the peaks of Per18 and Per21 are around 100 μm. The same is true for B1-b, where B1-bW has a peak at 8 μm, while B1-bN & S peak at around 100 μm. On the other hand, NGC 1333 SVS13 and L1448 N might appear non-coeval as a result of misalignment, given that the SED peaks of all three components are at about the same wavelength. For NGC 1333 SVS13, the outflow and continuum detections of this object suggest that it might be transitioning to the Class I stage and therefore be slightly non-coeval with its companions.

5.3. On coevality and non-coevality

The frequency of non-coevality found in the sample of resolved multiple systems studied in this work can be safely assumed to be due to real non-coevality and not solely to misalignment, since most sources are expected to present small misalignment rather than one component close to face-on. Non-coevality in our resolved sample is exhibited by higher order multiples, except for IC 348 MMS2, a triple with a component at a separation of ~3 toward the western source (Table 2) that appears coeval. Proto-binaries, on the other hand, tend toward coevality. Hence, protostellar siblings most of the time form and probably evolve simultaneously. This presents some interesting constraints to multiple system formation mechanisms and also raises questions.

For a multiple system to be non-coeval, the companion must either be formed after the first source or binary. In other words, fragmentation in the core must occur after the initial collapse and formation of the first protostar or protobinary. Hydrodynamical simulations predict that heated gas reduces the chance of further fragmentation, with only the cold gas tending to fragment (Stamatellos & Whitworth 2009b; Offner et al. 2010; Bate 2012). A possible explanation is that gas heating occurs along the outflow cavity while the dense envelope reduces heating of the surrounding gas, allowing further fragmentation to occur. Turbulence, on the other hand, is thought to be able to produce non-coevality through random density enhancements in the core.

For the coeval systems found in our resolved sample, fragmentation of the core would have occurred during the initial collapse and the system remained stable enough to hinder any further fragmentation and formation of younger companions, either through heated gas or lack of density enhancements and strong enough turbulence. Observational evidence for gas heating along the outflow cavity walls has been provided by van Kempen et al. (2009) and Yıldız et al. (2015), including some of the multiple systems studied here.

Dynamical evolution, that is, the interaction that occurs in multiple systems, can cause these systems to evolve non-coevally, for example by expelling one of the companions, considerably reducing its envelope and thus truncating its accretion of material (Reipurth et al. 2010; Reipurth & Mikkola 2012). However, for embedded systems, not enough time has elapsed for dynamical evolution to play a major role in the appearance of the system. External factors, such as neighbouring outflows and jets, can affect a protostar in a system, for example by stripping material or possibly triggering further fragmentation. While these mechanisms, dynamical evolution and external factors, are not formation mechanisms, they can alter the conditions of a multiple system and cause it to evolve non-coevally.

Given that only about a third of the multiple systems present non-coevality, the question then arises which factor or factors contribute to making some regions fragment and collapse even more while others do not. Probing the distribution of heated gas around multiple and single protostellar systems, in the case of multiples for both coeval and non-coeval, could provide insight.

6. Conclusions

This work presents the constructed SEDs of known resolved multiple protostellar systems in Perseus with separations 7. The SEDs were constructed from Herschel PACS photometric maps and Spitzer c2d catalogs, together with fluxes from the literature for the longer wavelengths. The properties were then derived from the observed SEDs. The SED shape and derived parameters were taken together with literature work on the envelope and outflow of these systems to determine the coevality of multiple systems. The literature, both from observations and models, lacks statistics on the frequency of alignment of multiple systems, but work is ongoing. A simple test whether the different SEDs might be due to misalignment was carried out in this work by randomly pairing model Class 0 sources from model SED grids with different inclinations with respect to the line of sight (Robitaille et al. 2006) and then counting the frequency of apparently non-coeval systems. The results of this work can be summarized in the following points.

  • 1.

    From our sample of resolved multiple protostellar systems,which have separations 7, a coevality frequency of 66 ± 10% is found, suggesting that most wide multiples are born together.

  • 2.

    From the observed SED shapes alone, a non-coevality frequency of 33 ± 10% is found, with higher order multiples being responsible for this percentage. Random pairing of model SEDs indicates that the frequency of apparent non-coevality that is due to misalignment of the components’ inclinations is 17 ± 0.5%, with most occurrences in systems with large misalignments and perpendicular orientations. But most pairs tend toward small misalignment (63 ± 0.6%). This indicates that the observed non-coevality toward multiple systems in our sample is not due to misalignment, but is instead real.

  • 3.

    Derived properties, such as Tbol and αIR, suggest that the non-coevality frequency may be lower (21 ± 9%). However, the parameters derived from the observed SEDs produce contradicting results. Physical parameters derived from the synthetic SEDs demonstrate that these parameters are very sensitive to inclination or do not provide clear-cut evolutionary stage separations. As previously found in Offner et al. (2012), the parameters are therefore unreliable for evolutionary stage classification, unless the inclination is well constrained.

  • 4.

    Unresolved SEDs from multiple systems are not always dominated by the primary or brightest component, but can present an odd double-peaked shape that is due to non-coeval components. This can alter the fraction of non-coevality, but we did not take this into account because of the high uncertainty.

Owing to the limit on multiple system separations in this work, our results can only place constraints on formation mechanisms at core scales (1600 AU). Higher order multiples show a stronger tendency to be non-coeval. This suggests that fragmentation at core scales can occur at different times, thus generating these types of systems.

The main conclusion of this work is then that real non-coeval multiple protostellar systems exist in the early stages of protostellar formation at core scales, which is most of the time due to formation, and can be enhanced by dynamical evolution (e.g., component ejection or envelope stripping by external influence, such as an outflow). Several questions then arise. What causes some cores to fragment and collapse at different times? Which role do the already formed protostars play on the formation of their companions? Does temperature play an important role, and if so, to which extent? Future work, both from observations and models, is needed to address these questions.


Retrieved Oct. 2015 from


Astrochemistry in Leiden is supported by the Netherlands Research School for Astronomy (NOVA), by a Royal Netherlands Academy of Arts and Sciences (KNAW) professor prize, and by the European Union A-ERC grant 291141 CHEMPLAN. This work is supported by grant 639.041.335 from the Netherlands Organisation for Scientific Research (NWO) and by the Netherlands Research School for Astronomy (NOVA) and by the Space Research Organization Netherlands (SRON). J.J.T. is currently supported by grant 639.041.439 from The Netherlands Organisation for Scientific Research (NWO). D.F. acknowledges support from the Italian Ministry of Science and Education (MIUR), project SIR (RBSI14ZRHR).


  1. André, P., Ward-Thompson, D., & Barsony, M. 1993, ApJ, 406, 122 [NASA ADS] [CrossRef] [Google Scholar]
  2. André, P., Men’shchikov, A., Bontemps, S., et al. 2010, A&A, 518, L102 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Arce, H. G., & Sargent, A. I. 2006, ApJ, 646, 1070 [NASA ADS] [CrossRef] [Google Scholar]
  4. Audard, M., Ábrahám, P., Dunham, M. M., et al. 2014, Protostars and Planets VI, eds. H., Beuther, R. S., Klessen, C. P., Dullemond, & T. Henning (Tucson: Univ. Arizona Press), 387 [Google Scholar]
  5. Balog, Z., Müller, T., Nielbock, M., et al. 2014, Exp. Astron., 37, 129 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Barsony, M., Ward-Thompson, D., André, P., & O’Linger, J. 1998, ApJ, 509, 733 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bate, M. R. 2012, MNRAS, 419, 3115 [NASA ADS] [CrossRef] [Google Scholar]
  8. Burkert, A., & Bodenheimer, P. 1993, MNRAS, 264, 798 [NASA ADS] [CrossRef] [Google Scholar]
  9. Chen, X., Launhardt, R., & Henning, T. 2009, ApJ, 691, 1729 [NASA ADS] [CrossRef] [Google Scholar]
  10. Chen, X., Arce, H. G., Zhang, Q., et al. 2010, ApJ, 715, 1344 [Google Scholar]
  11. Chen, X., Arce, H. G., Zhang, Q., et al. 2013, ApJ, 768, 110 [NASA ADS] [CrossRef] [Google Scholar]
  12. Chini, R., Hoffmeister, V. H., Nasseri, A., Stahl, O., & Zinnecker, H. 2012, MNRAS, 424, 1925 [NASA ADS] [CrossRef] [Google Scholar]
  13. Ciardi, D. R., Telesco, C. M., Williams, J. P., et al. 2003, ApJ, 585, 392 [NASA ADS] [CrossRef] [Google Scholar]
  14. Crapsi, A., van Dishoeck, E. F., Hogerheijde, M. R., Pontoppidan, K. M., & Dullemond, C. P. 2008, A&A, 486, 245 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Curtis, E. I., Richer, J. S., Swift, J. J., & Williams, J. P. 2010, MNRAS, 408, 1516 [NASA ADS] [CrossRef] [Google Scholar]
  16. Decin, L., Richards, A. M. S., Neufeld, D., et al. 2015, A&A, 574, A5 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Diolaiti, E., Bendinelli, O., Bonaccini, D., et al. 2000, in Adaptive Optical Systems Technology, ed. P. L. Wizinowich, SPIE Conf. Ser., 4007, 879 [Google Scholar]
  18. Duchêne, G., & Kraus, A. 2013, ARA&A, 51, 269 [Google Scholar]
  19. Duchêne, G., Monin, J.-L., Bouvier, J., & Ménard, F. 1999, A&A, 351, 954 [NASA ADS] [Google Scholar]
  20. Dunham, M. M., Crapsi, A., Evans, II, N. J., et al. 2008, ApJS, 179, 249 [NASA ADS] [CrossRef] [Google Scholar]
  21. Dunham, M. M., Stutz, A. M., Allen, L. E., et al. 2014, Protostars and Planets VI, eds. H., Beuther, R. S., Klessen, C. P., Dullemond, & T. Henning (Tucson: Univ. Arizona Press), 195 [Google Scholar]
  22. Dunham, M. M., Allen, L. E., Evans, II, N. J., et al. 2015, ApJS, 220, 11 [NASA ADS] [CrossRef] [Google Scholar]
  23. Eislöffel, J., Froebrich, D., Stanke, T., & McCaughrean, M. J. 2003, ApJ, 595, 259 [NASA ADS] [CrossRef] [Google Scholar]
  24. Enoch, M. L., Evans, II, N. J., Sargent, A. I., & Glenn, J. 2009, ApJ, 692, 973 [NASA ADS] [CrossRef] [Google Scholar]
  25. Evans, II, N. J., Dunham, M. M., Jørgensen, J. K., et al. 2009, ApJS, 181, 321 [NASA ADS] [CrossRef] [Google Scholar]
  26. Fendt, C., & Zinnecker, H. 1998, A&A, 334, 750 [NASA ADS] [Google Scholar]
  27. Froebrich, D. 2005, ApJS, 156, 169 [NASA ADS] [CrossRef] [Google Scholar]
  28. Gerin, M., Pety, J., Fuente, A., et al. 2015, A&A, 577, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Hartigan, P., & Kenyon, S. J. 2003, ApJ, 583, 334 [NASA ADS] [CrossRef] [Google Scholar]
  30. Hiramatsu, M., Hirano, N., & Takakuwa, S. 2010, ApJ, 712, 778 [NASA ADS] [CrossRef] [Google Scholar]
  31. Hirano, N., & Liu, F.-C. 2014, ApJ, 789, 50 [NASA ADS] [CrossRef] [Google Scholar]
  32. Hirano, N., Ho, P. P. T., Liu, S.-Y., et al. 2010, ApJ, 717, 58 [NASA ADS] [CrossRef] [Google Scholar]
  33. Hirota, T., Bushimata, T., Choi, Y. K., et al. 2008, PASJ, 60, 37 [NASA ADS] [CrossRef] [Google Scholar]
  34. Hirota, T., Honma, M., Imai, H., et al. 2011, PASJ, 63, 1 [NASA ADS] [CrossRef] [Google Scholar]
  35. Hsieh, T.-H., Lai, S.-P., Belloche, A., Wyrowski, F., & Hung, C.-L. 2015, ApJ, 802, 126 [NASA ADS] [CrossRef] [Google Scholar]
  36. Hull, C. L. H., Plambeck, R. L., Kwon, W., et al. 2014, ApJS, 213, 13 [NASA ADS] [CrossRef] [Google Scholar]
  37. Jensen, E. L. N., Mathieu, R. D., Donar, A. X., & Dullighan, A. 2004, ApJ, 600, 789 [NASA ADS] [CrossRef] [Google Scholar]
  38. Jørgensen, J. K., Bourke, T. L., Myers, P. C., et al. 2007, ApJ, 659, 479 [NASA ADS] [CrossRef] [Google Scholar]
  39. Jørgensen, J. K., van Dishoeck, E. F., Visser, R., et al. 2009, A&A, 507, 861 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Kenyon, S. J., & Hartmann, L. 1995, ApJS, 101, 117 [NASA ADS] [CrossRef] [Google Scholar]
  41. Kratter, K. M., Matzner, C. D., Krumholz, M. R., & Klein, R. I. 2010, ApJ, 708, 1585 [NASA ADS] [CrossRef] [Google Scholar]
  42. Kraus, A. L., & Hillenbrand, L. A. 2009, ApJ, 704, 531 [Google Scholar]
  43. Launhardt, R., Stutz, A. M., Schmiedeke, A., et al. 2013, A&A, 551, A98 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Lee, C.-F., Hirano, N., Palau, A., et al. 2009, ApJ, 699, 1584 [NASA ADS] [CrossRef] [Google Scholar]
  45. Lee, C.-F., Hasegawa, T. I., Hirano, N., et al. 2010, ApJ, 713, 731 [NASA ADS] [CrossRef] [Google Scholar]
  46. Lee, K. I., Dunham, M. M., Myers, P. C., et al. 2015, ApJ, 814, 114 [NASA ADS] [CrossRef] [Google Scholar]
  47. Lee, K. I., Dunham, M. M., Myers, P. C., et al. 2016, ApJ, 820, L2 [NASA ADS] [CrossRef] [Google Scholar]
  48. Looney, L. W., Mundy, L. G., & Welch, W. J. 2000, ApJ, 529, 477 [NASA ADS] [CrossRef] [Google Scholar]
  49. Maercker, M., Mohamed, S., Vlemmings, W. H. T., et al. 2012, Nature, 490, 232 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  50. McCabe, C., Ghez, A. M., Prato, L., et al. 2006, ApJ, 636, 932 [NASA ADS] [CrossRef] [Google Scholar]
  51. Monin, J.-L., Ménard, F., & Peretto, N. 2006, A&A, 446, 201 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Murillo, N. M., & Lai, S.-P. 2013, ApJ, 764, L15 [Google Scholar]
  53. Offner, S. S. R., Kratter, K. M., Matzner, C. D., Krumholz, M. R., & Klein, R. I. 2010, ApJ, 725, 1485 [NASA ADS] [CrossRef] [Google Scholar]
  54. Offner, S. S. R., Robitaille, T. P., Hansen, C. E., McKee, C. F., & Klein, R. I. 2012, ApJ, 753, 98 [NASA ADS] [CrossRef] [Google Scholar]
  55. Palau, A., Zapata, L. A., Rodríguez, L. F., et al. 2014, MNRAS, 444, 833 [NASA ADS] [CrossRef] [Google Scholar]
  56. Pezzuto, S., Elia, D., Schisano, E., et al. 2012, A&A, 547, A54 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Plunkett, A. L., Arce, H. G., Corder, S. A., et al. 2013, ApJ, 774, 22 [NASA ADS] [CrossRef] [Google Scholar]
  58. Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Raghavan, D., McAlister, H. A., Henry, T. J., et al. 2010, ApJS, 190, 1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Rebull, L. M. 2015, AJ, 150, 17 [NASA ADS] [CrossRef] [Google Scholar]
  61. Reipurth, B. 2000, AJ, 120, 3177 [CrossRef] [Google Scholar]
  62. Reipurth, B., & Mikkola, S. 2012, Nature, 492, 221 [Google Scholar]
  63. Reipurth, B., Mikkola, S., Connelley, M., & Valtonen, M. 2010, ApJ, 725, L56 [NASA ADS] [CrossRef] [Google Scholar]
  64. Robitaille, T. P., Whitney, B. A., Indebetouw, R., Wood, K., & Denzmore, P. 2006, ApJS, 167, 256 [NASA ADS] [CrossRef] [Google Scholar]
  65. Rodríguez, L. F., Anglada, G., & Curiel, S. 1997, ApJ, 480, L125 [NASA ADS] [CrossRef] [Google Scholar]
  66. Rodríguez, L. F., Anglada, G., & Curiel, S. 1999, ApJS, 125, 427 [NASA ADS] [CrossRef] [Google Scholar]
  67. Sandell, G., & Knee, L. B. G. 2001, ApJ, 546, L49 [NASA ADS] [CrossRef] [Google Scholar]
  68. Soderblom, D. R., Hillenbrand, L. A., Jeffries, R. D., Mamajek, E. E., & Naylor, T. 2014, Protostars and Planets VI, 219 [Google Scholar]
  69. Stamatellos, D., & Whitworth, A. P. 2009a, MNRAS, 392, 413 [Google Scholar]
  70. Stamatellos, D., & Whitworth, A. P. 2009b, MNRAS, 400, 1563 [NASA ADS] [CrossRef] [Google Scholar]
  71. Tafalla, M., Kumar, M. S. N., & Bachiller, R. 2006, A&A, 456, 179 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Tobin, J. J., Dunham, M. M., Looney, L. W., et al. 2015, ApJ, 798, 61 [NASA ADS] [CrossRef] [Google Scholar]
  73. Tobin, J. J., Looney, L. W., Li, Z.-Y., et al. 2016, ApJ, 818, 73 [NASA ADS] [CrossRef] [Google Scholar]
  74. van Kempen, T. A., van Dishoeck, E. F., Güsten, R., et al. 2009, A&A, 507, 1425 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Velusamy, T., & Langer, W. D. 1998, Nature, 392, 685 [NASA ADS] [CrossRef] [Google Scholar]
  76. Walawender, J., Bally, J., Kirk, H., & Johnstone, D. 2005, AJ, 130, 1795 [NASA ADS] [CrossRef] [Google Scholar]
  77. Walawender, J., Bally, J., Kirk, H., et al. 2006, AJ, 132, 467 [NASA ADS] [CrossRef] [Google Scholar]
  78. Whitney, B. A., Wood, K., Bjorkman, J. E., & Wolff, M. J. 2003, ApJ, 591, 1049 [Google Scholar]
  79. Yen, H.-W., Koch, P. M., Takakuwa, S., et al. 2015, ApJ, 799, 193 [NASA ADS] [CrossRef] [Google Scholar]
  80. Yıldız, U. A., Kristensen, L. E., van Dishoeck, E. F., et al. 2015, A&A, 576, A109 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Young, K. E., Young, C. H., Lai, S.-P., Dunham, M. M., & Evans, II, N. J. 2015, AJ, 150, 40 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: HIPE map makers and photometry

Table A.1

Aperture photometry results for 70 μm.

The three HIPE map-makers (High Pass Filter; Jscan map; MADmap) were tested to determine the best map for performing photometry. The test was made only on the 70 μm maps for Perseus. The method used involved performing aperture photometry on the source and the surrounding background at four positions. The aperture (12′′) used was the same for all source and background measurements. Ten sources were selected from different regions of Perseus, ranging from isolated to clustered sources.

The flux was calculated in the following manner: where is the raw flux, Bi is the background flux, Acorr is the aperture correction factor, and Fsource is the background corrected flux. Aperture correction values were taken from Balog et al. (2014). For an aperture of 12′′, the correction factor is of 0.802.

Table A.1 lists the background and aperture corrected results. This shows that the difference between maps is not significant. We have adopted JScanmap for photometry.

Appendix B: Evolutionary stage classification

The SEDs of the resolved systems are compared to average SEDs obtained by (Enoch et al. 2009; Fig. B.1) to determine by eye whether the multiple systems are coeval or not.

Physical parameters derived from the SED are known to be sensitive to inclination. This is confirmed by the properties derived from the model SEDs. The derived properties were calculated for the model SEDs in the same way as for the observed SEDs. Comparing the derived parameters, we find that Tbol varies widely with inclination, while Lfir/Lbol are independent of inclination. These results confirm previous work on the subject (Jørgensen et al. 2009; Launhardt et al. 2013). Figure B.2 shows the parameters for all the Class 0 models versus inclination.

thumbnail Fig. B.1

Average SEDs derived by Enoch et al. (2009). The classifications are defined based on the bolometric temperature. These average SEDs are used to compare to the constructed SEDs in this work.

thumbnail Fig. B.2

Physical parameters derived from the model SEDs from the grid of Robitaille et al. (2006) compared to the inclination w.r.t. line of sight. Edge-on is defined as 90°. The dashed red lines indicate the evolutionary stage boundaries. Top: bolometric temperature Tbol varies largely as it inclines more toward face-on (0°). The blue lines are models with a small difference in Tbol with inclination; these have dense envelopes. Green and magenta show the models with an increase in Tbol when the inclination is close to face-on. Middle: the infrared spectral index from 2 to 24 μm oscillates without relation to inclination. This is most likely an effect of the ice feature at 8 μm. Bottom: the luminosity ratio Lfir/Lbol is independent of inclination, but there is no clear boundary for systems at different evolutionary stages.

Appendix C: Resolved multiple systems

Appendix C.1: NGC 1333

NGC 1333 SVS13: located at the heart of NGC 1333, SVS13 is a quintuple system with components SVS13A1 and 2, VLA3, SVS13B and C (Tobin et al. 2016). Based on the velocity field, VLA3 and SVS13B are suggested to be a binary even though VLA3 is closer to SVS13A (Chen et al. 2009). One of the components in SVS13A is expected to be the driving source of HH 7-11 (Rodríguez et al. 1997; Looney et al. 2000). SVS13A is observed to have a prominent molecular outflow in the SE-NW direction with a moderate inclination angle and wide opening angles (Plunkett et al. 2013). These outflow characteristics together with the presence of a centimeter source and being the exciting source for HH 7-11 (Rodríguez et al. 1997) suggest that SVS13A is a Class 0/I transition object. SVS13C may be driving a N–S outflow possibly directed along the plane of the sky, but its outflow emission may be confused with other outflows (Plunkett et al. 2013). SVS13B does not show a clear outflow, which could be due to confusion with the outflow of SVS13A (Plunkett et al. 2013). SVS13 has been suggested many times to be non-coeval because components SVS13B and C are more embedded than SVS13A, which has an optical counterpart (Looney et al. 2000; Chen et al. 2009). SVS13B is not detected in the Herschel PACS 70μm maps, but is detected at λ ≥ 100μm, suggesting it is deeply embedded.

NGC 1333 IRAS 2: located in the west of the NGC 1333 region, the IRAS 2 system is composed of sources AC (Sandell & Knee 2001), although source C is expected to be a starless core. IRAS 2B is known to be confused with a field star at λ ≤ 8 μm (Rodríguez et al. 1999). IRAS 2A is typically classified as Class 0, while IRAS 2B is considered Class I, but here we find them to be coeval. IRAS 2A is well known because of its spectacular quadrupole outflow. The N–S outflow has a shell-like structure, while the E–W outflow is more collimated (Plunkett et al. 2013). Tobin et al. (2015) resolved components with a separation of 0.6 toward IRAS 2A and suggested that the southern component drives the E–W outflow, while the northern component drives the N–S outflow. The outflow of IRAS 2B runs parallel to the N–S outflow of IRAS 2A (Plunkett et al. 2013).

NGC 1333 IRAS 7: three systems in a common core (CLASSy N2H+ observations) of about 30 diameter. Per49 is located in the SE edge of the dust and gas core, while Per18 is located at the peak with Per21 13 to the SW. Per18 and Per49 were found to be close binaries with separations 0.3 (Tobin et al. 2016). This system of sources is associated with the Haro-Herbig object HH6. The outflow of this system has been less frequently observed, with candidate outflow lobes proposed to extend to around the NE of the SVS13A outflow (Plunkett et al. 2013). 12CO observations with CARMA show outflow lobes associated with Per21, but no clear outflow signatures toward the other sources.

NGC 1333 IRAS 4: the region contains several systems and IRAS 4B can be resolved at the resolution of Herschel, but IRAS 4A cannot. IRAS 4B (also referred to as IRAS 4C, Looney et al. 2000, and IRAS 4B2, Hull et al. 2014) is not detected in the Herschel observations, which may mean that it is still deeply embedded. Hull et al. (2014) detected outflows toward both sources, with IRAS 4B showing an N–S outflow and IRAS 4B driving a weak E–W outflow.

NGC 1333 IRAS 5: located at the western edge of NGC 1333, this protobinary composed of Per52 and Per63 does not show a prominent molecular outflow (Curtis et al. 2010).

NGC 1333 Per58+Per65: located to the north of the NGC 1333 core, there are indications of outflow from these sources, although their orientations are unclear because of confusion (Curtis et al. 2010).

NGC 1333 Per37: located along a filament in the northernmost region of NGC 1333, this triple protostellar system, identified by (Tobin et al. 2016), does not appear to have a strong outflow in the observations of Curtis et al. (2010).

Appendix C.2: L1448 and L1455

L1448 C: is a protobinary in the south of L1448, with the southern component having a projected position in the outflow of the northern source. L1448 C-N shows a more prominent high-velocity collimated outflow traced in 12CO and SiO, with a low velocity conic cavity observed in 12CO, while L1448 C-S shows a much weaker low velocity 12CO collimated outflow (Hirano et al. 2010). Both outflows are not aligned, but do not show signs of significant misalignment along the line of sight either.

L1448 N: also commonly known as L1448 IRS3, this system is a sextuple, with the B and C sources containing three and two components, respectively (Lee et al. 2015). All three sources drive molecular outflows, with components B and C almost parallel to each other and source A perpendicular to B and C (Lee et al. 2015). The observed outflows suggest no significant misalignment along the line of sight. The outflow from L1448 C was suggested to induce fragmentation in this core (Barsony et al. 1998).

L1448 IRS2: Tobin et al. (2016) found L1448 IRS2 to be a close binary (separation ~0.7). SCUBA 850 μm observations showed a continuum peak to the east of L1448 IRS2, which Chen et al. (2010) referred to as IRS2E, and together with SMA and Spitzer upper limits proposed to be a first core candidate. Although SCUBA 850 μm observations suggest a shared envelope, the separation of these two sources is 46 and therefore cannot compose a multiple system. The outflow of IRS2 is observed to be conical in the SE-NW direction, while the suggested outflow of IRS2E is composed of only the red-shifted lobe directed toward the SW (Hull et al. 2014; Chen et al. 2010).

Appendix C.3: IC 348

IC 348 Per8+Per55: this triple protostellar system shows a large-scale jet directed in the north-south direction and is associated with two Haro-Herbig objects, HH841 and HH842 (Walawender et al. 2006). Based on the large-scale jet, there is no indication of multiple components.

IC 348 MMS: a protobinary located in the southwest of IC 348 and associated with a strong north-south outflow where the overlapping of redshifted emission at the tip of the blueshifted lobe is suggested to be due to a change in environment and not a product of inclination; see Tafalla et al. (2006). The outflow is driven by the Class 0 western source, MMS2, and is also associated with HH797. The eastern source, MMS2E, was suggested to be a Class 0 proto-brown dwarf driving a weak outflow in the NE–SW direction by Palau et al. (2014). Both sources are found to be coeval, as suggested by Palau et al. (2014).

IC 348 SMM2: also referred to as Per16+Per28, this protobinary shows jet emission in the east-west direction for both components and a large-scale S-shape bend in the flow (Walawender et al. 2006). The outflows appear to be short and clumpy (Eislöffel et al. 2003; Walawender et al. 2006).

IC 348 Per32+EDJ2009-366: the H2 emission around these sources appears to be clumpy, with bow shocks pointing east and possibly not all clumps belonging to the respective sources (Eislöffel et al. 2003). Per32 is a close protobinary (~6) (Tobin et al. 2016) and was also found to be a low luminosity object (Dunham et al. 2008; Hsieh et al. 2015).

Appendix C.4: B-1

B1 Per6+Per10: the outflows of this protobinary are not well determined owing to their location with several neighbouring

systems whose outflows become entangled (Walawender et al. 2005; Hiramatsu et al. 2010). Per10 might drive a red-shifted outflow lobe directed NW, but it is uncertain whether the lobe belongs to Per10 or to a source about ~70 to the north (Hiramatsu et al. 2010). Per6 is also known as SMM3 and Per10 as SSTc2d J033314.4+310711.

B1-b: a triple system composed of Per41 to the west and B1-bN and B1-bS to the north and south, respectively. B1-bN and S are found to be very young based on 7 to 1.1 mm continuum data and lack of Spitzer detections at λ ≤ 24μm (Hirano & Liu 2014). The outflows for both sources are directed E–W, but the location of the blueshifted lobes suggests that the two sources are misaligned, since the north source has the outflow directed to the west, while the south source has it directed to the east (Gerin et al. 2015). Per41 appears to be older, with no detection in the millimeter regime. It is suggested that there is no emission in the millimeter because the envelope of this system is being striped off by neighboring outflows, which makes it appear to be more evolved (Hirano & Liu 2014).

Appendix D: Single protostellar systems

Single protostellar sources indentified in Perseus (Tobin et al. 2016) are listed in Table D.1 and the constructed SEDs are shown in Fig. D.1. Not all sources have an SED, either due to lack of fluxes in the literature, non-detection in the Herschel PACS maps, or both. This is denoted in the last column of Table D.1.

Table D.1

Sample of single protostellar systems.

thumbnail Fig. D.1

Constructed SEDs for single protostellar systems. Other details are the same as in Fig. 3.

Appendix E: Herschel catalog

This appendix contains the Herschel PACS flux catalog for the Perseus star forming region obtained in this work through PSF photometry with StarFinder (Diolaiti et al. 2000). The fluxes have been background and aperture corrected, with the aperture correction values tabulated in Balog et al. (2014).

Table E.1

Herschel PACS protostellar fluxes for Perseus.

All Tables

Table 1

Sample of resolved multiple protostellar systems (separation 7).

Table 2

Sample of unresolved multiple protostellar systems (separation <7).

Table 3

StarFinder photometry parameters.

Table 4

Statistics from the constructed SEDs.

Table 5

SED derived properties.

Table 6

Derived parameters of combined SEDs.

Table 7

Random pair statistics based on synthetic SEDs.

Table A.1

Aperture photometry results for 70 μm.

Table D.1

Sample of single protostellar systems.

Table E.1

Herschel PACS protostellar fluxes for Perseus.

All Figures

thumbnail Fig. 1

Cartoon of the definitions used in this work. More evolved sources are represented by larger disks, wider outflow cavities and less envelope material.

In the text
thumbnail Fig. 2

Herschel PACS stamps of resolved multiple protostellar systems in Perseus, except for NGC 1333 IRAS 4 which also shows IRAS 4A, an unresolved protobinary. Each stamp is 80× 80. 70, 100 and 160 μm are shown in blue, green and red, respectively. Blue symbols represent the components of a system, with circles denoting those with additional unresolved multiplicity and diamonds indicating those without (known) additional multiplicity. NGC 1333 IRAS 4A is marked in red.

In the text
thumbnail Fig. 3

Constructed SEDs for resolved multiple protostellar systems. Filled circles denote the fluxes without contamination from nearby sources. Triangles indicate upper limits. Squares show combined fluxes. Each SED is overlaid with the template SEDs from Enoch et al. (2009) for comparison (dashed lines).

In the text
thumbnail Fig. 4

Constructed SEDs for unresolved multiple protostellar systems. Other details are the same as in Fig. 3.

In the text
thumbnail Fig. 5

Summed SEDs from resolved multiple systems, showing how unresolved multiple systems might be composed. The dashed lines are the same as in Fig. B.1.

In the text
thumbnail Fig. 6

Randomly paired model SEDs in the four groups based on alignment. The numbers indicate the inclination of the respective component. Note that the SED pairs that appear non-coeval peak at the same wavelength for both components, which is not expected from components at different evolutionary stages.

In the text
thumbnail Fig. B.1

Average SEDs derived by Enoch et al. (2009). The classifications are defined based on the bolometric temperature. These average SEDs are used to compare to the constructed SEDs in this work.

In the text
thumbnail Fig. B.2

Physical parameters derived from the model SEDs from the grid of Robitaille et al. (2006) compared to the inclination w.r.t. line of sight. Edge-on is defined as 90°. The dashed red lines indicate the evolutionary stage boundaries. Top: bolometric temperature Tbol varies largely as it inclines more toward face-on (0°). The blue lines are models with a small difference in Tbol with inclination; these have dense envelopes. Green and magenta show the models with an increase in Tbol when the inclination is close to face-on. Middle: the infrared spectral index from 2 to 24 μm oscillates without relation to inclination. This is most likely an effect of the ice feature at 8 μm. Bottom: the luminosity ratio Lfir/Lbol is independent of inclination, but there is no clear boundary for systems at different evolutionary stages.

In the text
thumbnail Fig. D.1

Constructed SEDs for single protostellar systems. Other details are the same as in Fig. 3.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.