{'fraction_non_alphanumeric': 0.03802976912560414, 'fraction_numerical': 0.024647141509252304, 'mean_word_length': 4.438684679023981, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Gravitational fragmentation caught in the act: the filamentary Musca molecular cloud

Filamentary structures are common in molecular clouds. Explaining how they fragment to dense cores is a missing step in understanding their role in star formation. We perform a case study of whether low-mass filaments are close-to hydrostatic prior to their fragmentation, and whether their fragmentation agrees with gravitational fragmentation models. For this, we study the 6.5 pc long Musca molecular cloud that is an ideal candidate for a filament at an early stage of fragmentation. We employ dust extinction mapping in conjunction with near-infrared data from the NEWFIRM instrument, and 870 um dust continuum emission data from the LABOCA instrument, to estimate column densities. We use the data to identify fragments from the cloud and to determine the radial density distribution of its filamentary part. We compare the cloud's morphology with 13CO and C18O line emission observed with the APEX/SHeFI instrument. The Musca cloud is pronouncedly fragmented at its ends, but harbours a remarkably well-defined, 1.6 pc long filament in its Center region. The line mass of the filament is 21-31 Ms pc^-1 and FWHM 0.07 pc. Its radial profile can be fitted with a Plummer profile that has the power-index of 2.6 \pm 11%, flatter than that of an infinite hydrostatic filament. The profile can also be fitted with a hydrostatic cylinder truncated by external pressure. These models imply a central density of 5-10 x 10^4 cm^-3. The fragments in the cloud have a mean separation of 0.4 pc, in agreement with gravitational fragmentation. These properties, together with the subsonic and velocity-coherent nature of the cloud, suggest a scenario in which an initially hydrostatic cloud is currently gravitationally fragmenting. The fragmentation has started a few tenths of a Myr ago from the cloud ends, leaving its center yet relatively non-fragmented, possibly because of gravitational focusing in a finite geometry.


Introduction
Filamentary structures have been known to be common morphological features of the interstellar medium (ISM) for decades (e.g., Barnard 1927;Schneider & Elmegreen 1979;Dobashi et al. 2005;André et al. 2014). Observational works examining the physical properties of filaments have built a picture in which filamentary clouds fragment into denser clumps and/or dense cores that can further collapse to form new stars (e.g., Bally et al. 1987;Mizuno et al. 1995;Myers 2009;André et al. 2010André et al. , 2014Jackson et al. 2010;Schneider et al. 2010;Kainulainen et al. 2013). This phenomenology is apparent across a wide range of size-and mass-scales, spanning ∼0.1-100 pc in size and ∼1-10 5 M ⊙ in mass (e.g., Pineda et al. 2010;Hacar & Tafalla 2011;Schneider et al. 2010;Hill et al. 2011;Hacar et al. 2013;Goodman et al. 2014;Ragan et al. 2014). The new, sensitive dust emission observations by the Herschel Space Observatory have also manifested the universality of filamentary structures in molecular clouds, reinforcing a hypothesis that they may provide a dominant pathway to star formation in molecular clouds (for a review, see André et al. 2014). Explaining how filamentary structures evolve and form stars is then a necessary step in understanding the exact role of this pathway.
Analytical studies have established a viable framework for the existence of filamentary structures, regardless of their origin, as hydrostatic, self-gravitating, thermally pressurised cylinders (Ostriker 1964;Fischera & Martin 2012a). These cylinders can fragment in time-scales of ∼Myr when perturbed (Inutsuka & Miyama 1992, 1997Fischera & Martin 2012a), and even small density perturbations in the filament can collapse faster than the global collapse takes place in it (Pon et al. 2011). This gravitational fragmentation provides one viable model to explain the "dense cores" observed in filaments. The details of the fragmentation process in this model depend most crucially on the line mass of the filament at its initial condition, i.e., at the equilibrium configuration. The line mass at equilibrium, and hence the fragmentation process, can be affected by the possible pressure confinement by the surrounding environment (Fischera & Martin 2012a), support provided by magnetic fields (Fiege & Pudritz 2000a,b), accretion of gas onto the filament (Heitsch 2013a,b), temperature profile (Recchi et al. 2013), and rotation (Recchi et al. 2014). For a more comprehensive review, we refer to André et al. (2014).
One key open question related to filaments is whether i) the filaments are close to hydrostatic objects prior to their fragmentation, in which case they first have to detach from the turbulent regime of the cloud and then fragment, or ii) the evolutionary picture of filaments is dynamic and hydrostatic conditions are never prevalent in the cloud. This question can be addressed by first carefully identifying filaments that are at the earliest stages of their evolution, possibly not yet fragmented, and by comparing their structure with the predictions of the hydrostatic equilibrium models. However in practice, doing this is severely hampered by the fact that non-fragmented filaments that carry the imprints of the initial condition are rare. This alone indicates that the phase in which a filament is non-fragmented and yet CO-emitting (so that it would be identified as a molecular cloud) is relatively short. While several recent works have made efforts to estimate the structure of filaments in nearby clouds (e.g., Arzoumanian et al. 2011;Fischera & Martin 2012b;Juvela et al. 2012b;Palmeirim et al. 2013), none of these works to our knowledge has identified and analysed specifically non-fragmented filaments that would allow an assessment of the filament structure close to its initial condition.
Another key question is how the filaments fragment into dense cores that form stars. Analytic models predict the initially hydrostatic, infinitely long filament to fragment through gravitational instability at a wavelength of roughly four times the filament scale-height (Inutsuka & Miyama 1992), in a timescale comparable to the sound-crossing time (order of a Myr, e.g., Fischera & Martin 2012a). For finite-sized filaments, it has been hypothezised that the local collapse time-scale varies as a function of the position along the filament, and that the local collapse proceeds most rapidly at the ends of the filament (Bastien 1983;Pon et al. 2011). While testing these predictions is an active goal of observational studies (e.g., Jackson et al. 2010;Schmalzl et al. 2010;Busquet et al. 2013;Kainulainen et al. 2013;Takahashi et al. 2013), it is of great interest to find and characterise the fragmentation in as simple as possible filaments to isolate different processes playing a role.
The main technical difficulty in addressing the questions above is the general difficulty of tracing the distribution of molecular gas in the interstellar medium: all commonly-used tracers (dust emission, extinction, molecular line emission) have their specific problems and none of them can alone trace the wide range of column densities present in filaments. An additional difficulty is introduced by the small spatial scale of the filaments and their fragments. Filaments typically show widths that are on the order of 0.1 pc (e.g., Barnard 1927;Arzoumanian et al. 2011); fragments are similar in size. Consequently, single-dish observations can only study the structure of filaments in adequate detail in the most nearby molecular clouds (distances 200 pc). As a result, it remains unknown what the physical structure of filaments is like prior to their fragmentation and how exactly their fragmentation proceeds from that initial state to star formation.
To make progress in characterising the early evolution and fragmentation of filaments, we target in this study the remarkably long (∼6.5 pc), highly rectilinear Musca molecular cloud. The cloud is located in the Musca-Chamaeleon cloud complex (see Fig. 1) at the distance of ∼150 pc (Knude & Hog 1998). The cloud has been identified to harbor only a very small number of cold cores given its spatial extent (Vilas-Boas et al. 1994;Juvela et al. 2010Juvela et al. , 2011Juvela et al. , 2012b, and it contains only one T Tauri star candidate, attesting to its low current star formation activity (Vilas-Boas et al. 1994). We have recently shown using CO line observations that Musca is entirely subsonic and velocitycoherent (Hacar et al. submitted). Contrasting this apparent quiescence, we have also shown that the probability distribution of column densities in the cloud shows a power-law-like excess of high-column densities, a feature common for star-forming molecular clouds (Kainulainen et al. 2009(Kainulainen et al. , 2011. Similarly, the cloud's volume density distribution is more top-heavy than that of other nearby starless clouds, suggesting that Musca may be close to active star formation (Kainulainen et al. 2014), or possibly in transition between non-star-forming and star-forming states.
The above properties make Musca a perfect candidate to be in the earliest phases of fragmentation. In this paper, we present a detailed study of the column density structure of Musca. We use archival near-infrared (NIR) observations and 870 µm dust emission observations to map column densities through the cloud and to characterise its internal structure. We then place the results in the context of hydrostatic equilibrium models, establish their suitability to describe the cloud, study its fragmentation, and finally propose an evolutionary scenario for the cloud.

NIR data and dust extinction mapping
The Musca cloud was observed with the NIR wide-field camera NEWFIRM (Autry et al. 2003) mounted on the Blanco 4-m telescope at the Cerro Tololo Inter-American Observatory between 2011 April 11 and 13. We retrieved the raw data and associated calibration files from the NOAO public archive. Six pointings (5 covering the entire cloud and 1 off the cloud) were obtained in the J, H and K s broad-band filters. Each pointing was observed in dither mode. Table 1 gives a summary of the observation strategy. The instrument footprints are shown in Fig. 1, right panel; note that only four frames fall along the main part of Musca shown in the figure.
The individual raw images were processed using an updated version of Alambic (Vandame 2002), a software suite developped and optimized for the processing of large multi-CCD imagers. Alambic includes standard processing procedures such as overscan and bias subtraction, flat-field correction, bad pixel masking, CCD-to-CCD gain harmonization, sky-subtraction, and de-stripping. Aperture and PSF photometry were extracted from the individual images using SExtractor (Bertin & Arnouts 1996) and PSFEx (Bertin 2011). The individual catalogs were then astrometrically registered and aligned on the same photometric scale using Scamp (Bertin 2006), and the final mosaic produced using SWarp (Bertin et al. 2002). The photometric zero-points of the merged catalog were derived using the 2MASS catalog. The final absolute astrometric accuracy is expected to be better than 0. ′′ 1 and the absolute photometric accuracy (rms of the individual images zeropoints) of the order of 5%. The NEWFIRM data were complemented with J, H, and K s band photometric data from the 2MASS archive (Skrutskie et al. 2006) to cover a wider field towards Musca. The area covered by the NEWFIRM frames cover the area used in all analyses of this paper; the 2MASS data is merely used to complete the outskirts of a wide-field extinction map (Fig. 1, right panel).
We used the JHK S band data in conjunction with the nicest color-excess mapping technique (Lombardi 2009) to estimate the dust extinction through the Musca cloud. For the details of the technique we refer to Lombardi (2009). The technique is based on comparing NIR colors (H − K and J − H) and surface density of the stars that shine through the cloud to those  from a field that is assumed to be free from extinction. We note that since the photometric depths of the NEWFIRM and 2MASS data are different, these parameters were determined separately for the two data sets. In implementing nicest, a relative extinction law between JHK s bands has to be assumed. We chose the empirical extinction law of Cardelli et al. (1989) We emphasize that the technique estimates the dust extinction in NIR. We convert the estimates of NIR extinction to visual extinction by adopting the ratio between optical depths in NIR and visual wavelengths (Cardelli et al. 1989) The visual extinctions are finally used to estimate the total hydrogen (atomic + molecular) column densities using the ratio measured by Savage et al. (1977) and Bohlin et al. (1978) (see also Vuong et al. 2003;Güver &Özel 2009) The final extinction map of Musca is presented in Fig. 2, with higher-detail blow-ups shown in Fig. 3. The map has the resolution of 45 arcsec and uncertainty at low column densities of σ(A V )=0.2 mag. We analyse the column and volume density structure of the cloud using this map in Sections 3.1 and 3.2.

LABOCA dust emission data
In addition to extinction maps, we carried out in May 2010 870 µm continuum observations of the Musca cloud using the LABOCA bolometric camera (Siringo et al. 2009) at the APEX telescope. The observations consisted of fourteen individual raster-spiral maps combined to obtain a large-scale mosaic of approximately 2.2 by 0.4 deg 2 covering the total extension of this cloud. Guided by our previous extinction maps of Musca (Kainulainen et al. 2009), the location, coverage, and overlap of these submaps were selected to optimise the sensitivity of the final mosaic along the main axis of the cloud. Skydips plus flux calibrations and pointing corrections in nearby sources were carried out every 1.5-2 hours. Using the BoA software (Schuller 2012), the data reduction process followed the standard method described by Belloche et al. (2011) in order to recover the extended emission in the cloud. Convolved to a final resolution of 30 arcsec, the resulting map presents a typical rms of 17 mJy beam −1 and is shown in Fig.2.

SHeFI 13 CO and C 18 O line emission data
We employ in the discussion of this paper 13 CO (J = 2 − 1) and C 18 O (J = 2 − 1) line emission data of the Musca cloud from Hacar et al. (submitted), to which we refer for the description of the observations and data reduction. The data were gathered with the SHeFI instrument at the ESO/APEX telescope in Chile. The spatial resolution of the data is 28. ′′ 5 and the velocity resolution 0.10 km s −1 or 0.17 km s −1 (the receiver backend was changed during the observations). The spectra typically have a main beam temperature rms of 0.14 K. Line emission was observed in 305 positions along a curve that approximately traces the crest of the Musca cloud (see Fig. 1). The data were used to determine the central velocities, V LSR , and velocity dispersions, σ v , along the filament. They were also used to estimate the non-thermal contributions to the observed line-widths, σ NT and σ NT,corr (opacitycorrected), assuming a temperature of 10 K. For the details of the derivation of these values we refer to Hacar et al. (submitted).

Results
3.1. Column density structure of the Musca filament Figure 2 shows the column density map of the Musca cloud derived from the NIR dust extinction data. Figure 3 shows zoomins of the map. The highest column densities in the map are about 35 × 10 21 cm −2 . The mean column density above 2 × 10 21 cm −2 is N H = 4.5×10 21 cm −2 . Figure 2 also shows the contours of the LABOCA 870 µm emission data. The emission contours show good morphological correspondence with the NIR-derived column density map, evidencing that no column density peaks remain undetected from the NIR-derived data, e.g., because of the saturation of extinction or resolution effects. In general, the relationship between the extinction-and emission-derived column densities shows a high scatter, likely due to the spatial filtering of the dust emission data and their lower signal-to-noise compared to extinction data; the relationship is presented in Appendix A. We estimated the total mass of the cloud from the NIRderived column densities by summing up all column densities inside the contours 1 of N H = {1.5, 3, 4.5} × 10 21 cm −2 (see Fig. 3 for the contours). This resulted in the total masses of M = {281, 174, 132} M ⊙ . We used the total masses to calculate equivalent 2 line masses, M e l , of the entire Musca cloud. The length of Musca for all thresholds is about l = 6.5 pc, resulting in M e l = {43, 27, 20} M ⊙ pc −1 . The inclination angle, i, of the cloud is not known, and for reference, we adopt here i = 0. The inclination affects the line masses through a factor of cos i, and hence, the values we give represent upper limits. Note that in addition to the inclination uncertainty, the uncertainty of the line mass is significantly affected by distance uncertainty, which is at least some ∼10%.
The column density map ( Fig. 2) reveals numerous details within the large-scale Musca cloud. We divided the cloud into three regions, namely South, Center, and North, based on their morphology. The detailed column density maps of these regions are shown in Fig. 3.
The column density map (and emission data) show that Musca contains several small-scale condensations. To have a simple, quantitative census of these structures, we identify "fragments" in the observed region. We emphasize that our purpose here is merely to have a first look at the small-scale substructures in the cloud; a simple fragment definition will suffice for this purpose. With this in mind, we first smooth the resolution of the emission data down to the resolution of the extinction data (45 ′′ ) and identify fragment candidates as local emission maxima. In identifying the local maxima, we use contour steps of 0.03 mJy. If the fragment candidate is smaller than 45 ′′ in diameter, it is rejected. As a second step, the fragment candidates are inspected against the local maxima of the extinction map. Local maxima are determined from the extinction data using contour intervals of 1.5×10 21 cm −2 . If the lowest contour that outlines a fragment candidate does not overlap with the lowest contour outlining an extinction maximum, the candidate is rejected. The pixel with maximum emission value is recorded as the location of the fragment. As a result of these steps, our fragments are essentially local emission maxima that have counterparts in the extinction map. The procedure results in identifying 21 fragments from the observed area.
We estimated the masses associated with the fragments by summing up all column densities from the NIR-derived data within circular apertures. The radii of the apertures, listed in Table 2, were chosen by eye. The fragment masses span 0.3-7.1 M ⊙ ( Table 2). We note that all column densities within the apertures were included in the calculation, i.e., possible "background" was not considered; it is not our goal to analyse the internal structure of the fragments or their relation to their surroundings. Finally, we roughly estimated the mean densities of the fragments from their area and mass, assuming spherical symmetry (Table 2).
In the following, we describe the specific properties of the South, Center, and North regions separately. The South region is clearly fragmented into five core-like structures, the most massive of which contain several solar masses. We measured the separation between the fragments as distances between neighbouring fragments; in case of the South (and Center) this is trivial, because the fragments are arranged into a line along the filament crest. The mean separation in South is 0.78 ± 0.6 pc. Note that the mean is strongly affected by the long distance between fragments #3 and 4; without the separation between those two fragments the mean peak separation is 0.50 ± 0.2 pc. The total mass of the South region above N H = {1.5, 3, 4.5} × 10 21 cm −2 are M = {96, 88, 74} M ⊙ . The region is about 2.7 pc long, resulting in the equivalent line masses of M e l = {36, 33, 27} M ⊙ pc −1 . Finally, we note that fragments #4 and #5 practically compose the object studied as "a filament" by Juvela et al. (2012b). It appears to us that this region is clearly fragmented, and therefore, the axisymmetric cuts across the region measure in effect a combination of the profile of the fragments and that of the collapsing filament, not the profile of the filament itself.
The Center region of the Musca cloud consists of a 1.6 pc long filament that shows less sub-structure than the South region and has a relatively constant peak column density of N H ≈ 10 −15 ×10 21 cm −2 along its crest. The line mass of this filament using the column density thresholds of N H = {1.5, 3, 4.5} × 10 21 cm −2 results in M l = {31, 26, 21} M ⊙ pc −1 . Thus, its line mass is practically the same as that of the entire Musca cloud on average and that of the South and North (see below) regions. The FWHM of the filament is 0.07 pc, measured from its mean column density profile (see Section 3.2). There are five fragments along the main filament and one about 0.5 pc off the filament. The fragments along the filament have very low contrast to their surroundings compared to those in the South (see the quantification of the column density contrast later in this Section). The mean separation of the five fragments along the crest is 0.35 ± 0.13 pc (fragment #9 was not included in the analysis). We consider physical models for the Center filament in Section 3.2.

Fig. 2.
Column density map of Musca (colour scale), derived using near-infrared dust extinction mapping and data from NEWFIRM and 2MASS. The spatial resolution of the map is FWHM = 45 ′′ and the sensitivity σ(N H ) ≈ 0.4 × 10 21 cm −2 . The white boxes indicate three regions that are shown in blow-ups in Fig. 3. The contours show the 870 µm continuum emission contours, observed with the APEX/LABOCA instrument and smoothed down to the resolution of the extinction data. The contours are drawn at the intervals of 0.03 mJy/beam, which is the 3-σ rms noise of the smoothed data.
The North region of the Musca cloud is clearly fragmented, similarly to the South region. In the south-west corner of the region, there is a complex area containing four fragments. This region is characterised by two velocity components, likely resulting from two overlapping structures (Hacar et al. submitted). In the middle of the North region, there is a curious structure that consists of two thin (N H ≈ 3 × 10 21 cm −2 ) structures parallel to each other (a bubble-like structure), and of two peaks (fragments #19 in northwest and fragments #17 and 18 in southeast) that are elongated almost in perpendicular direction to those two structures. Unfortunately, we did not cover this region with molecular line data, and thus the nature of the structure remains un- Fig. 3. Blow-ups from the NIR-derived column density map of the Musca cloud (Fig. 2). The contours start at N H = 1.5 × 10 21 cm −2 and the step is 1.5 × 10 21 cm −2 . The lowest contour is drawn with a dotted line. The dashed lines in the Center region indicate the 1.6 pc velocity-coherent filament whose radial structure is modelled in Section 3.2. The numbered signs refer to the fragments identified within the cloud (see Section 3.1).
clear. There is a strong local maximum at the north-east corner of the region. This peak coincides with the location of the only candidate young star (T Tauri star candidate, Vilas-Boas et al. 1994) in the cloud. The North region harbours ten fragments ( Table 2). These fragments are not aligned along a line as clearly as fragments in Center and South. Acknowledging this caveat, we calculated the separations along two routes that connect the fragments (in the order of the fragment number): one including fragment #19 and excluding fragments #17 and 18, and the other excluding fragment #19 and including fragments #17 and 18. This resulted in mean separations of 0.24 ± 0.16 pc and 0.25 ±0.17 pc, respectively. The total masses of the region above We emphasize that the North and South regions contain highcontrast density enhancements compared to the smoother, filamentary Center region; even though we identify "fragments" from the Center region, their column density contrast against their surrounding, more elongated gas component is very low. We further quantify this morphological difference by measuring the standard deviation of column density (extinction) values along the approximate crest of the cloud (see Fig. 1 right panel). The histograms of the extinction values along the crest in different regions are shown in Fig. 4. The figure shows that the distributions of North and South regions are clearly different from that of the Center region: the North and South peak at around N H ≈ 6 × 10 21 cm −2 and have strongly skewed distributions, while the Center peaks at around N H ≈ 13 × 10 21 cm −2 and has a narrow distribution. The relative standard deviations of the extinction values along the cloud crest are: 55% for South, 49% for North, and 15% for the Center. This shows that the column density variation in the Center is clearly smaller, by a factor of ∼3, than in the North and South.

Physical models for the Center filament
The column density maps of the Center region (Figs. 2, 3) revealed an about 1.6 pc long, relatively non-fragmented filament. In this section, we examine the radial structure of this filament and consider physical models for it.
We estimated the mean column density profile of the Center filament from radial cuts across it. For this purpose, we first defined the filament crest by connecting the local maxima within the N H = 8 × 10 21 cm −2 contour with lines. The filament crest is shown in Fig. B and the points defining it are given in Table B.1. We then sampled the radial profile with cuts perpendicular to the crest. The positions of the cuts along the crest were separated by 4 pixels to ensure that the profiles at the crest position are not correlated. Consequently, the 1.6 pc long filament was sampled with 20 profiles. Figure 5 shows the resulting mean profiles, given separately for the west and east sides of the filament. Importantly, the profiles show that the filament is not symmetric, but the west side of it is steeper than the east side. The CO spectra show two velocity components on the east side of the filament crest (cf., Hacar et al. submitted), suggesting that another, faint structure overlaps with the main crest at the east side. This asymmetry is also seen at larger scales in the extinction map of the region (see Fig. 1). Therefore, we consider that the west side profile reflects the true filament profile more correctly. It is noteworthy that even if Musca seems to represent one of the simplest filamentary clouds we know, its radial profile is asymmetric with respect to the main axis when resolved in high detail. Similarly, as shown in Section 3.1, small scale variations are present even in the simple Musca filament. These phenomena are likely present in more complex filamentary clouds, such as L1495-B213 (e.g., Schmalzl et al. 2010) and IC5146 (e.g., Arzoumanian et al. 2011). This advices caution when analysing profiles of complex filaments by averaging individual profiles and need to be addressed for the correct Notes. (a) The peak value of dust emission from the fragment is more than 45 ′′ offset from the associated extinction peak. The structure is considered as a fragment, because the contours outlining the local maximum in emission overlap with the contours outlining a local maximum in extinction. The mass estimate is based on extinction data; the aperture for the mass calculation was centered on the nearest extinction maximum. (b) The contours defining the local maximum in emission map overlap with two local extinction maxima. The mass estimate is based on extinction data; the aperture for the mass calculation was centered on the nearest extinction maximum. (c) The fragment is offset by about 0.5 pc from the filament crest. Consequently, we do not include the fragment into the fragmentation analysis.
physical interpretation of the derived parameters. In the following, we quantify the radial structure of the Center filament in Musca by fitting different physically motivated models to its observed profile.

Isolated, infinite filament in hydrostatic equilibrium
We first modelled the 1.6 pc long filament in the Center region in the context of an infinite, isolated cylinder in equilibrium between thermal and gravitational pressures (Ostriker 1964). Such an equilibrium configuration has the line mass where c s is the sound speed. This indeed is close to the observed line mass of the Center filament (M l = 21-31 M ⊙ ). We note that it is possible that non-thermal motions provide additional support for a filament, making its critical line mass higher. An estimate of this "turbulent support" can be calculated from Eq. 4 by replacing the sound speed with an effective sound speed, c s,eff , that results from the total line width in the cloud. The non-thermal line-width of C 18 O) is about 0.7c s , yielding c s,eff = √ 1 + 0.7 2 c s = 1.22c s . The critical line mass is affected by the square of this, yielding 24 M ⊙ pc −1 .
The volume density profile of the infinite, hydrostatic filament in equilibrium is (Ostriker 1964) where r 0 is the scale radius as defined by Ostriker (1964), ρ c the central density, and the exponent p=4. Note that if the exponent is not 4, the radial profile does not describe an equilibrium solution. Using the generic, Plummer-like profile is however useful, because it allows for quantification of r 0 , ρ c , and p through a fit to an observed profile, and from therein, examination of the validity of the model. The scale radius is coupled to the central density and temperature, T , through The column density profile corresponding to the volume density profile is (e.g., Arzoumanian et al. 2011) In the equation, ) p/2 , with i as the inclination angle of the filament, which we assume to be zero. Figure 5 shows the best fit of Eq. 7 to the observed mean profile, achieved by treating ρ c , r 0 , and p as free parameters, fixing the fitting range to r = [0, 0.6] pc, and then minimising the chi-square between the data points and the model. The inverse variances of the 20 individual profile cuts were used as the weights of the data points in the fit. The best-fit model for the west-side profile has the exponent p = 2.6±11%, the scale radius r 0 = 0.016 ± 33% pc, and the central density n c = 54000 ± 50% cm −3 ( Table 3). The quoted errors are the one-sigma uncertainties of the fit. One should note that fitting Eq. 7 in this manner suffers from degeneracy between ρ c , r 0 , and p (Juvela et al. 2012a;Smith et al. 2014). This degeneracy likely causes the large uncertainties in r 0 and ρ c . The slope p suffers less from the degeneracy and in general can be constrained with lower uncertainty (Smith et al. 2014).

Isothermal filament in equilibrium, confined by external pressure
We also modelled the filament in the Center region as a pressureconfined (and truncated), isothermal filament in hydrostatic equilibrium (Fischera & Martin 2012a). The column density profile of this model is In the equation the dimensionless coordinate x is the radius, r, normalized by the outer radius of the filament, R f , i.e., x = r/R f , and f cyl = M l /M l,crit is the ratio of the filament's line mass to the critical line mass, (see Eq. 4 for M l,crit ). The relative shape of the column density profile depends only on the parameter f cyl and on the normalisation of the radial profile, i.e., on R f . The temperature and line mass affect the profile shape implicitly through f cyl . In addition, the normalisation of the profile is affected by the confining external pressure p ext .
The applicability of this model can be immediately examined through f cyl , which in the context of this model cannot be larger than unity. In Section 3.1, we estimated the line mass of the filament to be M l = {31, 26, 21} M ⊙ when using the column density thresholds of N H = {1.5, 3, 4.5} × 10 21 cm −2 to compute the line mass. The critical line mass without the contribution from nonthermal motions (16 M ⊙ pc −1 ) is slightly below this range, and the critical line mass with the non-thermal contribution (24 M ⊙ pc −1 ) overlaps with it. One should note that the observed line masses are upper limits because of the assumption of zero inclination. Distance uncertainty is also at least some 10%. We also do not know exactly the temperature structure of the filament, which also affects M l,crit . Finally, in the context of this model, the filament is pressure-confined by the surrounding gas; this surrounding gas possibly contributes to the observed M l value (we address this further below). We conclude that the observed line masses are such that considering the pressure-confined filament model is reasonable.
We fitted Eq. 8 to the observed data points at radii below 0.6 pc by letting f cyl and p ext change as free parameters, again using the inverse variances of the 20 profiles to weight the data points. This did not yield an acceptable fit. However, allowing also the fitting radius, i.e., the radius of the filament, to change as a free parameter and ignoring the profile outside it resulted in a reasonable fit (see Fig. 5). The best fit has the parameters f cyl = 0.66 ± 12%, and p ext = 12 ± 37% × 10 4 K cm −3 , and R f = 0.12 ± 7% pc where the uncertainties result from the chi-square fit ( Table 3). The external pressure is consistent with pressures needed to confine isolated globules modelled as Bonnor-Ebert spheres (∼2-20×10 4 K cm −3 , Kandori et al. 2005), and a factor of a couple higher than inferred for parsec-scale filaments (∼1-5×10 4 K cm −3 , Fischera & Martin 2012a,b).
One has to recognise that the above model does not include any contribution to the observed column densities from the surrounding gas that is required to provide the confining pressure (see the discussion in Recchi et al. 2014). Thus, the use of Eq. 8 by itself may not be physically entirely meaningful. In practice, the confining pressure could arise from the thermal pressure in a two-phase medium, in which the cooler, denser component is surrounded (and confined) by the less dense, warm component. Accretion of gas into the filament could also possibly provide confining pressure (Heitsch 2013a,b). Turbulent ram pressure could also provide confining pressure, although to provide confinement, it should be relatively isotropic, which may be unrealistic. Regardless of the nature of the pressure, a surrounding gas component is needed, and it can affect the observed column density structure.
To explore the effect of the surrounding gas to the parameters of this model, we made an experiment in which the pressureconfined filament is embedded inside a cylindric envelope (e.g., Hernandez et al. 2012;Recchi et al. 2014). This cylindric envelope has a power-law radial density distribution, ρ e (r) = Cr −2 , where C is a proportionality constant. The column density profile of the envelope results then from integrating ρ e (r) along the line of sight (see, e.g., Recchi et al. 2014), yielding where R e and R f are the radius of the envelope and filament, respectively. At the limit r → 0 the equation has the value C(−R −1 e + R −1 f ), which gives the column density of the envelope at r = 0, i.e., N 0 . The total column density of the model is the sum of the envelope and the embedded filament (Eq. 8). This model was fitted to the observed profile with f cyl , p ext , R f , and C as free parameters, assuming that the envelope goes to zero at 0.6 pc, i.e., R e = 0.6 pc (Fig. 6). We stress here that our goal is to probe how the filament parameters are affected by an envelope, not to find the best possible envelope model, nor even a physically meaningful one. It is evident that the envelope model we use is rather arbitrary, and other envelopes may well produce equally good fits.
We find that the observed profile is best fit when a filament is surrounded by an envelope that has the central column density N 0 = 2.4 × 10 21 cm −2 . The parameters of the embedded filament, i.e., f cyl and p ext (see Table 3) are not significantly different from those of the filament without the envelope. In conclusion, a pressure-confined hydrostatic filament, embedded in a low-density envelope, seems to provide a reasonable explanation for the entire observed column density profile.

Relationship between the column density and velocity structure in Musca
We briefly compare in this section the column density structure of the Musca cloud to the large-scale kinematic data from Hacar et al. (submitted). While Hacar et al. discusses the properties of the global velocity field in Musca, we put here those data in the context of the column density structure and fragments of the cloud. Figure 7 shows the NIR-derived column densities, 13 CO and C 18 O central velocities (V lsr ), non-thermal velocity dispersions (σ NT ), and non-thermal opacity-corrected velocity dispersions (σ NT,corr ) for the ∼300 line-of-sights observed along the Musca cloud (see Fig. 1). Note there are some double-peaked spectra at the edges of the cloud; Fig. 7 shows the V lsr and σ NT values for both these components, obtained via fits of separate Gaussian components (see Hacar et al. submitted, for details). The locations of the fragments that are closer than 100 arcsec of the cloud crest are also shown in the figure. As found out by Hacar et al. (submitted), Musca is velocity-coherent and subsonic at densi-ties probed by the C 18 O data (sub-/transsonic at densities probed by the 13 CO data) practically throughout its length. The velocity gradient along the South region is about dV lsr /dr ≈ 0.5 km s −1 pc −1 . The South region also shows local velocity gradients, up to ∼3 km s −1 (cf., Hacar et al. submitted), that clearly coincide with the column density fragments. The North and Center regions show a velocity gradient of dV lsr /dr ≈ 0.15 km s −1 pc −1 , with the only local gradients coinciding with the densest fragment in the North region (#21) and the complex clump of multiple maxima (close to fragment #13). In particular, the two most pronounced fragments in the Center region (#6 and 7 in Table  2)  Summarising, the above properties demonstrate the fact that even if Musca as a whole is subsonic and velocity-coherent, the opposite ends of it show slightly more active internal dynamics than the extremely quiescent Center filament. This activity is likely related to the fragments present in the ends of the cloud.

Musca: caught in the middle of gravitational fragmentation
What do our observations tell about the nature of Musca and its fragmentation? We have demonstrated that the physical structure of the extremely quiescent Center filament of Musca is in a good agreement with the predictions for the line mass, radial density distribution, and velocity structure of the pressureconfined hydrostatic equilibrium model. The cloud contains 21 fragments (according to our definition of fragments), with the higher-contrast ones located in the South and North regions. Specifically, the North and South have a factor of three higher column density variance compared to the Center region. Some of the fragments contain several solar masses, and thus have the potential to be pre-stellar cores. Importantly, the equivalent line masses of the North and South are equal with the line mass of the Center filament, giving rise to a hypothesis that they may have fragmented from an initial condition that resembles the relatively non-structured Center filament.
To further examine this hypothesis, we consider the prediction of the gravitational fragmentation model (the pressureconfined filament, Fischera & Martin 2012a) for the fragmentation length and time-scale of a cylinder that has the properties of the Center filament. In this model, filaments with line masses less than the critical line mass can fragment when perturbed. Filaments with line masses larger than critical (that are not equilibrium solutions to begin with) should collapse into a spindle rapidly. The wavelength of the fastest growing perturbations is between λ m = 5 − 6.3 × FWHM cyl (Fischera & Martin 2012a). The lower limit corresponds to the case f cyl → 1 and the upper limit to f cyl → 0. We measured the FWHM of the Center filament to be 0.07 pc, which results in the fragmentation scale of λ m ≈ 0.35 − 0.44 pc. The observed mean peak separation of the fragments in all regions is 0.41 ± 0.37 pc and the median is 0.34 pc (mean is 0.32 ± 0.18 pc and median 0.33 pc without the long separation between fragments #3 and 4), overlapping with the predicted interval. The results are practically the same regardless of how we compute the separations in the North region (see Section 3.1). We also note that fragment #9 is not included in this analysis because of its location off the filament axis. In conclusion, it appears that also from the fragmentation scale pointof-view, the North and South regions may have originated from an initial condition that resembles the Center filament.
The fragmentation time-scale of the infinite, pressure-confined cylinder depends on f cyl , T , and p ext (Fischera & Martin 2012a), and is between where τ Kp = kT/(4πµm H Gp ext ) ≈ 0.2 Myr. This gives the range τ m ≈ 0.2 − 0.8 Myr when using the parameters of the bestfitting pressure-confined model (Table 3). We note the caveat that the above predictions for the fragmentation scale and timescale are for infinite cylinders; it may be that the predictions differ in finite geometry. Why has the fragmentation proceeded strongly at the ends of the Musca cloud, but not at its center? One possible mechanism in the context of gravitational fragmentation is provided by finite-size effects in cylindric geometry (as opposed to the infinite models of Ostriker 1964 andMartin 2012a). In finite filaments, the acceleration of gas depends on the relative position along the filament, being largest at the ends of it (Bastien 1983;Pon et al. 2011Pon et al. , 2012Clarke & Whitworth 2015, see also Burkert & Hartmann 2004 for a study in a sheet geometry). These higher accelerations may translate into shorter local collapse time-scales, and therefore, into formation of fragments at the ends of the filament prior to its centre (Pon et al. 2011(Pon et al. , 2012. For example, Pon et al. (2011, see their Figs. 5 and 6) suggests that the local collapse might occur a factor of two-tothree faster at the ends of the filament than at its center. Similarly, Clarke & Whitworth (2015) found in their numerical simulations that the first fragments form at the ends of the filament, regardless of its aspect ratio. While there are some quantitative differences between the predictions of Clarke & Whitworth (2015) and Pon et al. (2011Pon et al. ( , 2012, both models predict that long (but finite-sized) filaments fragment starting from their ends. Consequently, they lay out the basis for explaining the fragmentation at the ends of Musca as a result of the end-dominated collapse caused by higher accelerations, and therefore, shorter local collapse time-scales.
In the models of Pon et al. (2011Pon et al. ( , 2012 and Clarke & Whitworth (2015), high-aspect ratio filaments collapse globally and the fragments formed at the ends of the cloud fall into its centre; this can be used to gain further insight into the evolution of Musca under this framework. The time-scales for the global collapse (time-scale for the end fragments to reach the cloud centre) are (following the notation of Clarke & Whitworth 2015) and t pon = 0.44 A (Gρ 0 ) −1/2 A 5 0.98 A 1/2 (Gρ 0 ) −1/2 A 5, where G is the gravitational constant and ρ 0 the central density.
The aspect ratio of Musca depends on how we define its width. If we adopt the ∼0.6 pc extent of the Plummer profile (see Fig. 5), the aspect ratio 6.5 pc / 0.6 pc ≈ 5 follows; if we adopt the 0.1 pc radius of the pressure-confined filament, the aspect ratio is 6.5 pc / 0.2 pc ≈ 32. This interval of A = 5 − 33 is very conservative, as it represents the extreme choices for the filament width. For the central density, we use 10 5 cm −3 . The global collapse time-scale for the chosen range of aspect ratios is then between ∼0.5-3 Myr, depending on the model. We showed earlier that the fragmentation time-scale of an infinite pressure-confined filament is between 0.2-0.8 Myr. Even if these estimates are highly uncertain, they suggest that the fragmentation time-scale may be significantly shorter than the global collapse time-scale. This is in agreement with the fact that we do not observe Musca to be globally collapsing, but rather forming multiple fragments along it. This, in turn, suggests that the age of Musca is rather closer to its fragmentation time-scale than its global collapse time-scale, i.e., only some tenths of a Myr. Alternatively, additional processes must be supporting Musca against global collapse. Instead of the gravitational focusing, another possibility to explain the differences in fragmentation along Musca may be differences in the initial physical properties along it. In principle, the fragmentation time-scale anti-correlates with the line mass (e.g., Fischera & Martin 2012a); however as shown earlier, the equivalent line masses of the North and South are the same as the line mass of the Center, stating that they should collapse in the same time-scale. It could also be hypothesised that difference in the dynamics of the regions plays a role: the North and South regions show significantly higher global velocity gradients than the Center region, which might affect fragmentation in them. Whether such an initial velocity gradient could promote fragmentation remains to be explored via theoretical works. Finally, the time-scale also correlates with the size of the initial density perturbations (Pon et al. 2011). It cannot be ruled out that the North and South regions would have had higher density variance along them compared to the Center, inducing faster fragmentation.
Ultimately, one would like to compare the time-scales derived above with other age estimators in the cloud. However, the only direct sign-post of age in Musca is the single T Tauri star candidate in the North region (fragment #21). Its presence indicates that fragmentation in the North region (and supposedly in South) must have taken place on the order of ∼0.5 Myr ago (Class 0+1 protostellar lifetime, Dunham et al. 2014), with an admittedly large uncertainty. We showed earlier that the timescale of gravitational fragmentation in Musca is 0.2-0.8 Myr, coinciding with this estimate. Further studies targeting, e.g., the chemical ages at the different locations along the filament would provide additional, independent measures of age.
Altogether, the morphology of Musca, its fragmentation, radial structure, and velocity characteristics suggest that gravitational fragmentation is taking place in situ in Musca, moulding it towards star formation. This is in a qualitative agreement with our earlier results regarding the cloud's column density structure: in Kainulainen et al. (2009) we analysed the probability density function of column densities (N-PDF) in Musca, showing that it exhibits (at least a start of) a power-law-like tail at high-column densities, similarly with actively star-forming clouds. Further, in Kainulainen et al. (2014) we showed with the help of volume density modelling that Musca currently contains only a very small amount of gas (∼6 M ⊙ ) at densities high-enough for star formation. From the point-of-view of its column and volume density structure, the cloud must become more condensed locally to form stars; this could happen via the on-going fragmentation and gravitational collapse. Musca represents a perfect target for studying this: the simplicity of the cloud is ideal for isolating the different processes driving fragmentation.

Conclusions
We analyzed in this work the column density structure and fragmentation of the 6.5 pc long, rectilinear Musca molecular cloud that shows very little signs of on-going star formation. We examined the structure of the cloud with the help of near-infrared dust extinction data and 870 µm dust emission data. We also examined the large-scale kinematic structure of Musca using 13 CO and C 18 O line emission data (from Hacar et al. submitted). Our conclusions are as follows.
1. The column density data show that the Musca cloud is fragmented into 21 fragments that have masses between ∼0.3-7 M ⊙ . The Center region of Musca harbours a remarkably well-defined, 1.6 pc long filament that contains only lowcontrast fragments. The filament shows only 15% relative column density variations along its crest, in stark contrast to ∼50% variations in the more pronouncedly fragmented ends of the cloud. We emphasize that the relative column density variations provide a physically motivated, simple measure of how fragmented a filament is; a low degree of fragmentation is necessary for characterising the initial condition of fragmenting filaments. 2. The line mass of the Center filament is 21-31 M ⊙ pc −1 , i.e., close to the critical value of a cylinder in hydrostatic equilibrium. The equivalent line masses (total mass divided by length) of the North and South regions are practically the same. 3. The radial density structure of the Center filament is well described with either i) a Plummer profile that has a power-law exponent of 2.6, scale radius of 0.016 pc (FWHM of 0.07 pc), and central density of 5 ×10 4 cm −3 , or ii) a hydrostatic filament that is confined by external pressure imposed by a surrounding envelope. In this case the ratio of the line mass to the critical line mass is 0.66, central density is 10 5 cm −3 , and external pressure is ∼10 5 K cm −3 . We note that even this simple and well-resolved filament shows an asymmetric column density profile with respect to the main filament axis. This advices caution when analysing profiles of filaments that are more complicated than Musca. 4. The fragmentation length of the cloud, ∼0.4 pc, agrees with the prediction of the gravitational fragmentation of an infinite, pressure-confined cylinder. The fragmentation timescale predicted by the model is 0.2-0.8 Myr and the global collapse time of the cloud 0.5-3 Myr. 5. The ends of Musca show somewhat more variation in their C 18 O velocity structure than the extremely quiescent Center region, both in the central line velocities and the scatter of non-thermal line widths. The central line velocities peak at the locations of the fragments, suggesting motions associated with the fragments. 6. The above characteristics give rise to a hypothesis for the evolutionary scenario of Musca: an initially 6.5 pc long, close-to hydrostatic filament is in the middle of gravitational fragmentation. Fragmentation started a few tenths of a Myr ago, and it has so far taken place only at the ends of the cloud, possibly because of the shorter local collapse timescales compared to the cloud centre ("gravitational focusing"). The relatively non-fragmented Center region may still represent the close-to hydrostatic initial condition of the filament, regardless of whether the fragmentation at the ends of the cloud is driven by gravitational focusing or some other process.

Fig. A.1.
Relationship between emission-derived (N(LABOCA)) and extinction-derived (N(NIR)) column densities. The red symbols show the data in the area of the most massive fragment in the cloud (#5). The dashed line shows the one-to-one relationship.