Open Access
Issue
A&A
Volume 628, August 2019
Article Number A33
Number of page(s) 16
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201935545
Published online 06 August 2019

© J.-F. Robitaille et al. 2019

Licence Creative Commons
Open Access article, [Article], under the terms of the Creative Commons Attribution License ([Article]), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

A good statistical characterisation and morphological analysis of the interstellar medium (ISM) is important for many astrophysical studies. Identifying the general gas density distribution of molecular clouds as a function of their hosted star formation activity allows us to recognise the dominating physical process of the region and thus to make a link between the ISM structure formation and the emergence of star formation activity. A detailed decomposition of the signal received at different wavelengths is also fundamental in order to characterise correctly the properties of the different foreground Galactic components, such as the temperature and column density of ISM gas, and the extragalactic components, such as the cosmic infrared background (CIB) and the cosmic microwave background.

For these reasons, a reliable morphological analysis of interstellar maps is needed. During the last decades, some statistical tools became the foundations of many theories of ISM structure formation and star formation, such as the Fourier power spectrum (Crovisier & Dickey 1983; Green 1993; Miville-Deschênes et al. 2003), the probability distribution function (PDF; Padoan et al. 1997b; Federrath et al. 2008; Schneider et al. 2013; Burkhart et al. 2017), and the Δ-variance (Stutzki et al. 1998; Ossenkopf et al. 2008a,b), a method that filters and averages the structures of different sizes l in a map to produce a spectrum showing the relative amount of structure as a function of structure size. Column density PDFs and the Δ-variance slope have proven to be strongly dependent on the type of forcing (compressive or non-compressive) present in turbulent medium, and are sensitive to turbulence intermittency (Federrath et al. 2009, 2010). From theory and molecular cloud simulations, it has been proposed that turbulent motions are the main cloud-shaping mechanism and that they produce a lognormal low (column) density PDF (e.g. Padoan & Nordlund 2002; Hennebelle & Chabrier 2008; Ward et al. 2014; Burkhart 2018) followed by a power-law tail due to self-gravitating gas (e.g. Kritsuk et al. 2011; Girichidis et al. 2014). This scenario is supported by observations using Herschel dust column density maps or extinction maps (e.g. Schneider et al. 2013, 2015; Kainulainen et al. 2013; Pokhrel et al. 2016; Alves et al. 2017), while a pressure governed power-law tail is proposed by Kainulainen et al. (2011) and Tremblin et al. (2014).

In addition, fluctuations of some physical properties in the ISM, such as the density, can be so large that the average value provides inadequate information. For instance, in their analysis of the statistical properties of the line centroid velocity in a turbulent, compressible, but gravitationless simulation, Lis et al. (1996) found that the global PDF of centroid velocity is close to Gaussian. However, PDFs of the centroid velocity increments, i.e. the two-point statistics for a centroid velocity map separated by a distance l, show non-Gaussian wings increasing towards small values of l. The same behaviour has also been measured on CO data for the Polaris and Taurus fields (Hily-Blant et al. 2008).

The dilution effect of averaging large density or velocity fluctuations over a field is also affecting the Fourier power spectrum and the Δ-variance analysis. As noted by Panopoulou et al. (2017), even if the distributions of structure widths in a field, in simulated data or in observations, demonstrate clearly the existence of a dominating characteristic scale, the spatial power spectrum analysis still show a unique power law attributed to a scale-free medium. The power spectrum analysis of the Polaris field observed by Herschel is a good example (Miville-Deschênes et al. 2010; Men’shchikov et al. 2010) where the field is dominated by highly contrasted structures and where the spatial power spectrum has a single power law without any break pointing at a characteristic scale. Similar results were obtained on Polaris using the Δ-variance method (Stutzki et al. 1998), which shows the same information content as the Fourier power spectrum (Bensch et al. 2001; Ossenkopf et al. 2008a). In addition, a new technique recently introduced by Ossenkopf-Okada & Stepanov (2019) that compares the power of isotropic and anisotropic structures shows that Polaris has an almost scale-free filamentary spectrum. Schneider et al. (2011) and Elia et al. (2014) concluded, by applying the Δ-variance on nearby molecular clouds and regions across the Galactic plane, that despite the presence of characteristic scales, the underlying cloud structure is self-similar. This discrepancy between the common scale-free medium measured in the ISM and the presence of highly contrasted filamentary structures remains a fundamental issue in our understanding of the density distribution of the ISM. From these results we can conclude that the typical power spectrum analysis used to measure the hierarchical nature of the ISM fails to identify the typical sharp transitions in density and filamentary structures in the ISM. These structures are nevertheless physically important since they are crucial to the mechanisms of star formation (e.g. Elmegreen et al. 2001).

Historically, Crovisier & Dickey (1983) and Green (1993) measured the first angular power spectrum of HI emission directly from interferometric data. Even if Green (1993) admitted that various structures in the ISM, such as sheets and filaments, dominate atmultiple spatial scales, the fact that there was no preferred angular scale measured in the HI emission power spectrum was interpreted as a sign that turbulence must play a significant role in the hierarchical structure of the ISM. In order to test the fractal nature of the ISM, Elmegreen et al. (2001) compared HI emission maps of the Large Magellanic Cloud (LMC) with fractal models made from the inverse Fourier transform of random complex-number noise in the u-v plane multiplied by a power law. These models are often called fractional Brownian motion (fBm). The resulting fractal models were also exponentiated in order to reproduce the log-normal PDF usually obtained in simulations of magnetohydrodynamical turbulence. However, even if the fractal model respects similar one-point and two-point statistics compared to the LMC (the PDF and the power spectrum analysis, respectively) the model lacks all the usual structures associated with the ISM, such as filaments, holes, and shells. Fractal models fail to reproduce the common sharp transition in intensity seen in the ISM. Recently, Elia et al. (2018) has shown that fBm models are not a good approximation of the ISM, but that multifractal analysis offers a more complete characterisation of molecular cloud structures.

In the light of these past studies, it is important to recall some fundamental properties of fully developed turbulence. As can be seen in hydrodynamical simulations, a turbulent field can be described as a superposition of some random distribution, as produced by fBm simulations, and a set of localised and coherent structures, which also demonstrate a spatial hierarchical scaling (Farge 1992). These coherent structures are sometimes identified as a manifestation of intermittency. Federrath et al. (2010) summarised the signature of intermittency in three manifestations: (1) non-Gaussian wings in density and/or velocity PDFs; (2) anomalous scaling of the higher-order (p ≥ 4) structure functions of the velocity field and velocity increments, implying that the statistics are increasingly non-Gaussian at small scales; and (3) coherent structures of intense vorticity and of strong shocks. In this paper intermittency is considered in a broad sense as irregularities and alternation in the spatial statistical distribution of ISM properties and more specifically for density fluctuations in the case of the present study. This definition corresponds closely to “the dual nature of molecular clouds” described in the review of Falgarone et al. (2004), where the diffuse component, traced by the 12CO (J =1–0) line emission, is fractal and highly dynamical and the coherent (as opposed to fractal) component, traced by mid-infrared absorption and submillimetric dust thermal emission, is accurately described by a network of filaments and dense cores.

As demonstrated by Elmegreen et al. (2001) for HI emission and Miville-Deschênes et al. (2007) for dust far-infrared emission, if exponentiated fBms are able to reproduce the non-Gaussian wings of log-normal PDFs, these mock fractal simulations fail to reproduce the typical coherent structures in the ISM. Furthermore, Robitaille et al. (2014) have shown, by applying for the first time the segmentation method described in the present paper, that exponentiated fBms also fail to reproduce the non-Gaussian wings of PDFs measured as a function of spatial scales. In the same way as for the PDF analysis of centroid velocity increments, Robitaille et al. (2014) showed that dust emission at 250 μm also has more important non-Gaussian wings towards small spatial scales.

In this paper we present a novel decomposition technique based on complex-wavelet power spectrum analysis that we have developed to perform an in-depth analysis of the ISM signal1. This new technique can be applied on any data set, for example column density and velocity centroid maps. Contrary to the Fourier power spectrum, this new technique is sensitive to the dense and coherent filamentary structures. By merging the multiscale PDF analysis with the power spectrum analysis, the technique exposes, in the case of density fluctuations, the dual or even plural nature of molecular clouds and how the diffuse medium is linked to the dense coherent structures. A new complementary method, also using complex-wavelet transforms, the Reduced Wavelet Scattering Transform, was also recently proposed by Allys et al. (2019). Compared to our technique, this method focuses on the detailed statistical description of non-Gaussian structures in the ISM rather than on the extraction and on the impact of such structures on the Fourier power spectrum.

Since this paper focuses on the transition between the two regimes of non-coherent and coherent structures in the ISM, our in-depth analysis is performed on the Polaris region, which represents the early stages of star formation activity in a molecular cloud. In future works, this analysis technique will be applied to centroid velocity maps and numerical simulations where the intermittentbehaviour of the density and velocity field can be compared.

The paper is organised as follows. An overview of the power spectrum analysis is presented in Sect. 2 and the wavelet power spectrum and the Multiscale non-Gaussian Segmentation technique (MnGSeg) are presented in Sect. 3. In Sect. 4 we applied MnGSeg to the Herschel column density image of Polaris located at 150 pc, identified the signature of the CIB, and revealed a characteristic scale. Finally, the results are discussed in Sect. 5, and a conclusion is presented in Sect. 6.

2 Assumptions behind the power spectrum analysis

The classical power spectrum is usually calculated in Fourier space. A schematic representation of 2D Fourier space, or the u-v plane, is shown in Fig. 1. The Fourier transform decomposes the signal f(x)2 into a linear combination of Fourier coefficients defined as (1)

where the wavenumber k describes the spatial frequency content of the signal or the image. Each Fourier coefficient is a complex number from which can be calculated the amplitude, , and the phase, ϕ = arctan(Im ∕ Re). The power spectrum analysis of a signal is a statistical measure of the amount of power |A|2 as a function of the spatial frequency k. In an ideal world, the experiment leading to the structure formation in the ISM would be reproduced several times under the same initial condition in order to average the different outcomes and to get an adequate statistical sample. This methodology is obviously impossible to achieve in our context. Consequently, we are forced to assume the ergodicity of the medium, so that the local intensity fluctuations averaged over many samples is equal to the spatial average of intensity fluctuations of one realisation. Usually, for the Fourier power spectrum, the information in 2D Fourier space is averaged over the azimuthal angles θ shown in Fig. 1, so that (2)

According to the Kolmogorov theory of turbulence (Kolmogorov 1941), the velocity power spectrum of an isothermal, subsonic, and non-compressive turbulent medium follows a power law over the spatial frequencies (3)

where γ is the power-law index. The power spectrum is equivalent to the Fourier transform of the second-order structure function (Boldyrev 2002). The pth-order structure function is defined as (4)

where f(x) referred generally to the fluid velocity v(x) has a power law of ζp = p ∕ 3. For the second-order longitudinal structure function a power law of ζ2 = 2 ∕ 3, which relates to the Fourier power spectrum index as γ = ζ2 + 1 = 5 ∕ 3. For a 3D incompressible and isotropic turbulent medium, the Fourier power spectrum then becomes P3D (k) ∝ k−2k−5 ∕ 3k−11 ∕ 3. Experimental evidence suggests that ζp is lower than p ∕ 3 for p ≥ 4 (She & Leveque 1994; Boldyrev 2002; Padoan et al. 2004). This implies that the velocity fluctuations are increasingly non-Gaussian on small scales, a phenomenon also referred to as inertial-range intermittency (Frisch 1995).

The scale-free property (i.e., a scaling function with a unique power law) of the turbulence assumes that the velocity components or the density distribution of a gas dominated by turbulent motions are random variables. Consequently, all the structural properties of a turbulent medium is contained in its power spectrum. The same arguments were used for the development of the Δ-variance analysis (Stutzki et al. 1998), which shows the same information contained in the power spectrum (Ossenkopf et al. 2008a). It can be shown that this method is equivalent to the Fourier power spectrum smoothed by the filter spectrum at each scale or size l (Farge 1992; Bensch et al. 2001).

Although the Fourier power spectrum is useful to describe the intensity distribution of a map, notably when intensity fluctuations as a function of spatial scales are largely random, it fails to provide an accurate description of the distribution of more complex medium. In the case of the ISM gas distribution, where compressive mechanisms and the intermittency of turbulence produced dense filamentary structures, intensity distributions as a function of spatial scales are no longer random nor isotropic, as is the case for instance in fractal simulations. These inhomogeneities in the medium also break the ergodic assumption. The spatial average of intensity fluctuations is no longer representative of intensity fluctuations occurring locally in the map.

The Fourier power spectrum loses all the local information associated with intensity fluctuations as a function of the spatial scales. Because the trigonometric functions associated with the Fourier coefficients ( in Eq. (1)) oscillate forever, all the information content of f(x) is completely delocalised (Farge 1992). As a solution to this limitation, Robitaille et al. (2014) shows that the analysis of the wavelet power spectrum of an image allows access not only to the spatial frequency content of the signal, but also to the information on the localised intensity fluctuation in an image as a function of the spatial scales. The gain in information is substantial and can be used to localise the transitions to high-intensity regions, perhaps associated with important changes in the main physical mechanisms at play.

In contrast to the Δ-variance spectrum, which averages all the wavelet coefficients as a function of spatial scale or size l and therefore loses all the local intensity fluctuation information, the MnGSeg technique, as primarily described by Robitaille et al. (2014), isolates the random and isotropic component of a map as a function of the spatial scales by analysing the PDF of complex-valued and directional wavelet coefficients before analysing the power spectrum. This method has the advantage of separating the map component that satisfies the ergodic assumption from the dense and anisotropic structures, such as the ubiquitous interstellar coherent filaments, which normally bias the Fourier power spectrum analysis. These two components refer to “the dual nature of molecular clouds” described by Falgarone et al. (2004). In addition to the non-biased power spectrum, which can measure the true scale-free nature of a map, this multiscale segmentation technique allows us to separate density structures contributing to the non-Gaussian part of the PDFs, i.e. structures that correspond to intermittency in the simple isotherm model, but in the context of the more complex ISM, may also correspond to self-gravitating structures, for instance.

thumbnail Fig. 1

Schematic representation of 2D Fourier space, where u and v are the two dimensions, k is the wavenumber, and θ is the azimuthal angle.

Open with DEXTER

3 The MnGSeg technique

This sectionreintroduces the MnGSeg procedure described by Robitaille et al. (2014) and improves on it.

3.1 Wavelet power spectrum

Wavelet transforms are designed to analyse local fluctuations in a signal. The wavelet transform is obtained with the convolution of a map f(x) with a family of translated and dilated wavelets generated from the mother wavelet function ψ(x): (5)

As a result, a low or a high wavelet coefficient means that at the position x and spatial scale l, the signal has a small or a large variation compared to the mean value of the signal. As was true for the Fourier transform, the wavelet transform respects the Plancherel relation (6)

where (7)

and (8)

Equation (7) represents the global wavelet power spectrum. This relation is true for all wavelet functions ψ(x). However, because some functions have a better resolution in frequency space, some wavelets estimate the power spectrum of a signal more accurately. Kirby (2005) showed that the Morlet wavelet is the best wavelet function for reproducing the Fourier power spectrum. The Morlet wavelet is complex-valued and anisotropic. It is defined in Fourier space as (9)

where the constant |k0| is set to to ensure that the admissibility condition3 is almost met (Kirby 2005). As defined in Eq. (9), in Fourier space the complex Morlet wavelet is equivalent to a Gaussian kernel that can easily sample spatial frequencies as a function of the azimuthal angle θ (see bottom left panel of Fig. 2). In contrast to a real isotropic wavelet, complex-valued wavelets with an azimuthal dependency allow us to estimate the true power, as defined in Eq. (2).

With this additional azimuthal dependency, it is also possible to estimate the power spectrum of an image by integrating over θ instead of over x as in Eq. (7). Because of the finite azimuthal resolution of the Morlet wavelet, the integration is exchanged for a discrete summation over a limited number of angles. Kirby (2005) showed that the optimal angle interval for an efficient and uniform sampling is (see top right panel of Fig. 2). Equation (7) thus becomes (10)

where are the Morlet wavelet coefficients for map f(x) and Nθ = Δθδθ is the number of directions θ needed to sample Fourier space over the range Δθ. Since, for a real image, quadrants 3 and 4 represented in Fig. 1 are redundant, they are the complex conjugates respectively of quadrants 1 and 2, only the angles in quadrant 1 and 2 need to be sampled, which leads to Δθ = π (see bottom right panel of Fig. 2). The convolution operation for the wavelet transform can be done directly in Fourier space, so that , where denotes the inverse Fourier transform and represents the Fourier transform of a “daughter” Morlet wavelet dilated at a given scale l and rotated at an azimuthal angle θ.

Compared to the Fourier power spectrum analysis, the complex Morlet wavelet power spectrum analysis, as defined in Eq. (10), is not only dependent on the spatial scale l, but also on the map position x. This property provides a far more complete description of intensity fluctuations as a function of spatial scale in a map. From PW (l, x) we can recover the global wavelet power spectrum by averaging the power over all positions x, (11)

where Nx = Nx × Ny corresponds to the number of pixels in the map. By converting the spatial scale l to the Fourier wavenumberk using k = |k0|∕l, we can compare the global wavelet power spectrum of Eq. (11) directly with the Fourier power spectrum defined in Eq. (2).

thumbnail Fig. 2

Top left panel: real part of the Morlet wavelet. Top right panel: Morlet wavelet rotated over an azimuthal angle of π, also called the Fan wavelet (Kirby 2005). Bottom left panel: Morlet wavelet in Fourier space. Bottom right panel: Morlet wavelet rotated over an azimuthal angle of π in Fourier space.

Open with DEXTER

3.2 Non-Gaussian segmentation

The complex-wavelet power spectrum allows us to analyse the local variation in intensity and global power as a function of spatial scales. This complete description of the intensity fluctuations in the map can be used to isolate the random component linked to the scale-free nature of the ISM. The residues of this segmentation procedure does not satisfy the randomness and the ergodicity assumptionsand are called coherent structures because of their fundamental properties of being spatially correlated across the scales.

As for the previous segmentation analysis by Robitaille et al. (2014), in order to separate these two components we use the coherent vorticity extraction algorithm (Nguyen van yen et al. 2012; Azzalini et al. 2005). This iterative algorithm was initially developed as a method for determining the optimal denoising threshold among wavelet coefficients. In many cases, denoising consists in deleting the wavelet coefficients of a noisy signal whose modulus is below a threshold, usually found at small scales,and reconstructing the denoised signal from the remaining coefficients (Azzalini et al. 2005). In our case, the algorithm is applied at every spatial scale and as a function of the azimuthal direction. For our analysis, the noisy random coefficients are considered to be the scale-free component of the map and the other set of coefficients is considered to be the coherent structures of the map. The algorithm is defined as follows: let Φ be the threshold splitting the non-Gaussian terms from the Gaussian terms in the wavelet coefficient distribution, and the function indicator. The threshold Φ is first estimated according to the variance (12)

where (13)

and (14)

The iterative calculation then converges to an optimal value of the threshold Φ which allows us to separate outliers from randomly distributed wavelet coefficients. The sequence is defined by (15)

where q is a dimensionless constant controlling how restrictive the definition of non-Gaussianities is. The first study using the MnGSeg technique chose the value of q following two criteria: the normal distribution of the power for the Gaussian features and the unique power law of its power spectrum. In the present paper, q is dynamical and thus dependent on the spatial scale. We considered that the amount of non-Gaussianities produced by turbulence intermittency and/or compressive physical processes (e.g. shocks) as a function of scale is unknown and that the non-Gaussianity contributions vary from one scale to another. Consequently, in order to adjust the parameter q to its optimal value at each spatial scale, when the algorithm converges to an optimal value for a threshold Φ, we calculate the skewness of the Gaussian wavelet coefficient distribution (i.e. the third moment of the distribution; see the skewness definition in Eq. (24)). If the skewness is greater than 0.4, the value of q is diminished by 0.1. This operation is repeated until convergence of the parameter q. The value 0.4 is justified by the fact that the distribution skewness is evaluated on the absolute value of complex-valued wavelet coefficients, and consequently that the distribution is not Gaussian, but Rician. Rician distributions are not completely symmetrical and have a non-zero skewness.

3.3 Power spectra analysis and reconstruction

After the convergence of the extraction algorithm, two sets of wavelet coefficients are obtained, the Gaussian set and the non-Gaussian set, also called the coherent set, . Then Eqs. (10) and (11) can be applied in order to calculate the power spectrum of both sets, PG (k) and PC (k). Since both power spectra are obtained from the squared amplitude of independent wavelet coefficients, the total power spectrum, equivalent to the Fourier power spectrum is simply the linear combination of the segmented components: (16)

As demonstrated by Robitaille et al. (2014), images corresponding to both set of coefficients can also be reconstructed. Originally, the reconstruction formula used the same synthesising wavelet as the analysing wavelet, the Morlet wavelet. However, thanks to the redundancy of continuous wavelets, many reconstruction formulas exit for a wavelet decomposed signal (Farge 1992). J. Morlet found empirically that even the delta function can be used to reconstruct the signal. In that case the reconstruction formula becomes (17)

where μ0 is the mean value of the original map and Cδ is a correction factor. The reconstruction was found to be optimal when the scale separation in logarithm Δ ln l is set to δθ, the same separation as for the azimuthal angle. This separation allows the construction a quasi-orthogonal set of wavelet coefficients representing the signal. The reconstruction is not perfect, but it is sufficient for most applications. The correction factor Cδ is set to σrσ0, where σr and σ0 are the standard deviations of the reconstructed map and the original map, respectively.

Because of the linearity of Eq. (17) and the linearity of wavelet transforms, (18)

where, (19)

and (20)

Here fG(x) and fC(x) are the reconstructed maps of the Gaussian and coherent wavelet coefficient sets, respectively.

4 Application on Polaris flare column density map

In this section, we apply the MnGSeg technique on the Polaris flare region. The data are briefly presented before doing the full analysis.

4.1 Data

We used the SPIRE 250 μm and the column density map of Polaris derived from Herschel imaging data taken as part of the HGBS and the “Evolution of interstellar dust” key programs (André et al. 2010; Abergel et al. 2010; Miville-Deschênes et al. 2010). The column density map was produced following a procedure described in most Herschel papers (see e.g. Palmeirim et al. 2013; Schneider et al. 2015; Könyves et al. 2015) and adopting a mean molecular weight per hydrogen molecule μH2 = 2.8 (Kauffmann et al. 2008). The column density map has an angular resolution of 36′′, corresponding to the half-power beam width resolution of Herschel/SPIRE 500 μm data and is estimated to be accurate to better than ~50% (e.g. Roy et al. 2014). The pixel size is 14′′. The 250 μm map has an angular resolution of 18′′ and a pixel size of 6′′. The Polaris column density map is shown in Fig. 3.

4.2 Fourier and wavelet power spectra

The Fourier and wavelet power spectra were calculated on the Polaris flare region located at high Galactic latitude (b ~ 25°). This region has no ongoing star formation activity. Only prestellar cores and unbound starless cores have been detected so far (André et al. 2010; Ward-Thompson et al. 2010).

The power spectrum analysis of this region was performed by Miville-Deschênes et al. (2010) on the three wavelengths observed by the SPIRE instrument, 250, 350, and 500 μm, and the IRAS 100 μm map. All power spectra (once corrected for the noise, the SPIRE beam, and the point sources contribution) show a straight power law with, within the uncertainties, a similar power law of − 2.7. This measurement also agrees with the previous power law measured by Stutzki et al. (1998) using the Δ-variance method on this low-density region properly traced by CO.

Figure 4 presents our power spectrum analysis of the Polaris flare region presented in Fig. 3. To reduce the map edge effects, the Fourier transform and wavelet transforms were done on a map ~ 1.25 times larger than the original one, where the frame pixel values are zero and an apodisation was applied over 3% of the original map edges. The mean pixel value of the map was subtracted prior to the apodisation to reduce the gap between the intensity of the signal and the zero value pixel frame. In order to produce the wavelet power spectrum, Eqs. (10) and (11) were calculated on the ~1.25 times larger map for every scale corresponding to the diamond symbol in Fig. 4. Miville-Deschênes et al. (2010) previously modelled the Polaris power spectrum as (21)

where (22)

The factor Γ(k) is the telescope transfer function, N(k) is the noise level, and P0 models the excess of power at small scales induced by point sources and the CIB associated with unresolved infrared galaxies at high redshift. In order to measure the power associated with Psky(k), the original power spectrum is first subtracted by the noise level and then divided by the telescope transfer function. For this analysis, the noise level is estimated by the last point of the wavelet power spectrum at k ≈ 1.75 arcmin−1. The fit values for the Fourier power spectrum are listed in Table 1. The fit is estimated between 0.05 and 0.8 arcmin−1. The measured power law is shallower than the previous measurements made on individual wavelength maps by Miville-Deschênes et al. (2010) (see Sect. 4.5 for a corresponding power spectrum analysis of the Herschel 250 μm map).

The wavelet power spectrum closely matches the Fourier power spectrum except for small deviations that are noticeable near the noise level.

thumbnail Fig. 3

Polaris column density map derived from Herschel imaging data taken as part of the HGBS and “Evolution of interstellar dust” key programs. The brightest structure with the highest column density value located on the south-western part of the field has been labelled as MCLD 123.5+24.9 and as the “saxophone” by Schneider et al. (2013).

Open with DEXTER

4.3 Intermittency

The good correspondence between the Fourier and wavelet power spectra validates Eq. (11) and suggests that the Fourier power spectrum is sensitive only to the mean variation of column density over the map as a function of spatial scales. However, as is seen in many numerical simulations and as is generally measured in column density PDFs of star formation regions, molecular cloud dense structures produce a tail at large densities on the column density distribution.These tails do not generally have a significant impact on the mean value of a statistical distribution. Large skewness has also been predicted as a function of spatial scales on centroid velocity increments by Federrath et al. (2010, their Fig. 9) for solenoidal and compressive forcing. PDFs should be close to Gaussian distributions (in semi-log plots) on large scales and present exponential tails on smaller scales. Concerning Polaris, they concluded that the kurtosis values measured on centroid velocity increments (CVIs) were compatible with intermittency of solenoidal (incompressible) forcing. This measure of intermittency was in good agreement with the CVIs of 12CO(1–0) IRAM map by Hily-Blant et al. (2008, their Fig. 4).

In this paper, we use the spatial scale filtering property of wavelet transforms as an alternative to the two-point statistics used to calculate CVIs. In addition, the spatial scale filtering was performed on a column density map rather than on a velocity centroid map. This choice allowed us to perform the multiscale analysis on a wider range of spatial scales and to investigate the intermittentbehaviour of the density field.

We recall that in the case of compressible and supersonic turbulence, the velocity field and density fluctuations become strongly coupled (Kritsuk et al. 2007). Large velocity-shears indeed produce intermittent structures in the velocity field that can follow boundaries of high-density structures traced in a column density map (Hily-Blant & Falgarone 2009). In addition, Federrath et al. (2010) showed that the PDFs of the logarithm of the density are usually roughly consistent with log-normal distributions for both solenoidal and compressive forcings. They attributed the non-Gaussian higher-order moment deviations of the distributions, such as the skewness and the kurtosis, to turbulence intermittency. They suggested that these deviations can be caused by head-on collisions of strong shocks or oscillations in very low-density rarefaction waves. In agreement with this interpretation, Hily-Blant et al. (2008) showed that for Polaris the centroid velocity structures associated with increments of CVI tails closely follow the boundaries of the optically thin 12CO emission traced by the broad 12CO line wings. Finally, Miville-Deschênes et al. (2003) showed that power spectra of integrated emission and centroid velocity fields of the high Galactic latitude HI cirrus Ursa Major have similar 3D spectral index and that their spatial fluctuations share similar structures, despite the moderate correlation. These results support the fact that for statistical measurements over a region, the amount of velocity fluctuations as a function of scale has an impact on the gas density fluctuations, even if velocity and density fluctuations are not expected to be perfectly correlated.

In Fig. 5, we display the normalised PDFs of the squared amplitude of complex-wavelet coefficients as a function of spatial scale. This plot corresponds to the intermittency measure as defined by Farge (1992): (23)

As described by Farge (1992), if I(l, x) = 1 for all x and for all l, then it means that there is no flow intermittency. In that case, each location x would have the same power spectrum, which corresponds to what we expect from a Fourier power spectrum. Figure 5 shows extreme departures from the mean power value at almost every spatial scales above 0.025 arcmin−1. The intermittency measure shows that many locations on the maps contribute more than 30 times the average over x to the Fourier power spectrum for a broad range of scales. The intermittency PDFs in Fig. 5 are calculated for a constant bin of 0.5; however, the number of statistically independent pixels varied as a function of spatial scales. The number of independent pixels is determined following the relation n =(Nx × Ny) ∕ l2. Figure 6 shows the intermittency measure in linear scale for three spatial frequencies smaller than 0.1 arcmin−1, where the number of independent pixels is respected.

The skewness of the intermittency PDFs was calculated as a statistical measurement of the increasing intermittency as a function of the spatial scale. The skewness is defined as (24)

where σl is the standard deviation of the intermittency measure I(l, x) for the given scale l. Figure 7shows the skewness value ς for each spatial scale converted in wavenumber k. It can be seen that ς(k) increases exponentially towards smaller scales until the scales become dominated by the noise level. The skewness value follows the fitted relation ς(k)= (103 ± 1) × k1.31 ± 0.08. A small deviation from the exponential fit is present at 0.02≲ k ≲ 0.08 arcmin−1, which corresponds to 0.1 ≲ l ≲ 0.3 pc. This deviation will be discussed further in Sect. 4.6.

thumbnail Fig. 4

Fourier (solid lines) and wavelet (diamonds) power spectra of the Polaris flare region presented in Fig. 3. The black line and black symbols represent the original power spectrum, P(k) in Eq. (21), while the blue line and blue symbols represent Psky (k) in Eq. (22).

Open with DEXTER
Table 1

Fit values for the column density map Fourier power spectrum.

thumbnail Fig. 5

Intermittency measure I(l, x) as a functionof the spatial frequency k. Frequencies from ~0.01 to ~ 0.75 arcmin−1 are plotted, and approximately correspond to scales where enough statistically independent points are available on large scales and where the signal is not dominated by noise on small scales.

Open with DEXTER
thumbnail Fig. 6

Intermittency measure in linear scale for three spatial frequencies. Compared to Fig. 5, the number of bins for each spatial scale respects the relation n = (Nx × Ny) ∕ l2, where k = |k0| ∕ l.

Open with DEXTER
thumbnail Fig. 7

Skewness, ς(l), for the intermittency measures defined in Eq. (23) and shown Figs. 5 and 6.

Open with DEXTER

4.4 Multiscale non-Gaussian segmentation

This sectionpresents the result of the MnGSeg technique applied on the Polaris flare region. Power spectra of both components are shown in Fig. 8, the Gaussian and the coherent parts, where the noise level is subtracted and the telescope transfer function is divided. As seen in the intermittency PDFs of wavelet coefficients in the previous section, intermittency in density starts to appear on large spatial scale leading to the segmentation of coherent features by the MnGSeg technique from k ≳ 0.02. It corresponds to a spatial scale of ~0.4 pc. The coherent features then dominate on all lower scales and follow closely the total power law measured in the Fourier power spectrum. The Gaussian part of the segmentation also has a power law followed by a complete flattening on small scales. According to the previous model of the Fourier power spectrum described in Eqs. (21) and (22), it is now possible to propose a more detailed model taking into account respectively the Gaussian and coherent power spectra, PG (k) and PC (k), (25)

and (26)

Since the noise component of the signal usually respects a distribution close to a normal distribution, it should now be present entirely in the Gaussian segmented part. This argument should also be generally true for the CIB component, except for the brightest and most nearby galaxies which are found in the small-scale non-Gaussian component. The CIB signature can be easily associated with the flattened part of the Gaussian power spectrum where coherent structures, on the same spatial scales, are still contributing to the coherent power law (i.e. for 0.3 ≤ k ≤ 1.0 arcmin−1). A comparison with the CIB component extracted from the Polaris region is done at 250 μm in the next section. The values for the power spectra fits, according to Eqs. (25) and (26) are summarised in Table 2. Both wavelet power spectra were fitted between 0.03 ≤ k ≤ 1.0 arcmin−1 after subtracting the noise level for the Gaussian part only, and divided by the telescope transfer function for both spectra. The coherent power-law fit corresponds, within the uncertainties, to the power law estimated for the total Fourier power spectrum. On the other hand, the Gaussian power-law fit is steeper and closer to the Kolmogorov power law of − 11∕3 for a 3D turbulent medium.

The Gaussian and coherent reconstructed maps following the relations 19 and 20 are shown in Fig. 9. The Gaussian map shows smooth features dominated by large-scale density fluctuations, as measured by the steeper power law of its power spectrum. The small-scale fluctuations are dominated by a granular component characterised by the flat part of the Gaussian power spectrum and dominated by the CIB signal. The PDFs associated to the Gaussian map are compared to a fBm in Fig. 10. Both maps are filtered for k ≳ 0.3 arcmin−1 in order to filter the CIB and the instrumental noise in the case of the Polaris map. The Gaussian PDFs show that the segmentation algorithm successfully removed most of the intermittent coefficients compared to the original intermittency measure in Fig. 5. After the segmentation, most of the intermittency values are in the range 0 < I(l, x) ≲ 2.5 and, in contrast with the original intermittency measure, all the distributions are now centred on 1.0 and are similar across the wavenumber k. Compared to the fBm simulation, the PDFs associated with the Gaussian density fluctuations of Polaris are broader.This may come from the fact that the Polaris region shows more areas of low density than the fBm. In fractal analysis, this aspect refers to the lacunarity of the structures (Mandelbrot 1982). Even if two fractal images share the same power law, their appearance can differ according to their different distribution I(l, x). In this case themain difference is that the fBm image, constructed in Fourier space, is made only of stationary fluctuations, which is not necessarily the case for the ISM density fluctuations.

The coherent map in Fig. 9 is dominated by the denser elongated structures. According to its power spectrum in Fig. 8, the non-Gaussian fluctuations are present on a broad range of spatial scales. A closer inspection of the coherent map indeed indicates that most of the filamentary structures are embedded in larger structures that have been identified as non-Gaussianities. No particular break is visible in the coherent power spectrum, which could be interpreted as no spatial scales being predominant for the non-Gaussianities. This result is in agreement with the recent analysis of Polaris filamentary structures done by Ossenkopf-Okada & Stepanov (2019), where the Polaris Flare shows an almost scale-free filamentary spectrum. However, as defined in Eqs. (10) and (11), the power spectrum is only sensitive to the mean value of the power distribution. As can be seen in Fig. 11, the intermittency measure for the non-Gaussian wavelet coefficients is largely skewed and not centred on I(k, x) = 1.0, which means that the coherent wavelet power spectrum, as is the case for the Fourier power spectrum, is not directly representative of the underlying power distribution. However, the unique power law associated with the non-Gaussianities is still a very interesting result and can be an indication that the non-Gaussianities as a function of scales are linked to the inertial range of the Gaussian part. This result might possibly be a direct measurement of the inertial-range intermittency detected on a density map of a star-forming region. This hypothesis is discussed further in Sect. 5.

thumbnail Fig. 8

Segmented power spectra for the Polaris flare column density map. The total Fourier power spectrum shown in Fig. 4 is represented by the solid black line. The red diamonds show the power spectrum for the Gaussian self-similar part of the map. The blue triangles show the power spectrum for the non-Gaussian coherent part of the map. The black dashed line represents the fits of the curves and the green dashed line represents the CIB and contribution level of the sources. The 3D Kolmogorov power law of − 11∕3 is plotted with the red dash-dotted line.

Open with DEXTER
Table 2

Fit values for the column density map Gaussian and the non-Gaussian power spectra.

thumbnail Fig. 9

Gaussian (left) and coherent (right) reconstructed column density maps following Eqs. (19) and (20); the correction factor Cδ has been multiplied to both maps, and the mean values of the initial map μ0 has been distributed to both maps.

Open with DEXTER
thumbnail Fig. 10

Top left panel: Gaussian reconstructed map filtered for scales k ≳ 0.3 arcmin−1. Top right panel: fractional Brownian motion simulation with the same power law as the Gaussian part of Polaris. The fBm is also filtered for k ≳ 0.3 arcmin−1. Bottom left panel: intermittency measure as defined in Eq. (23) for the Gaussian segmentation of wavelet coefficients. Bottom right panel: intermittency measure of the fBm simulation on the same spatial scales as Polaris.

Open with DEXTER
thumbnail Fig. 11

Intermittency measure as defined in Eq. (23) for the non-Gaussian wavelet coefficients.

Open with DEXTER

4.5 CIB measurement and extraction

The CIB corresponds to high-redshift starburst galaxies that are unresolved by far-infrared and submillimetre observations(Puget et al. 1996; Fixsen et al. 1998; Lagache et al. 1999). The energy peak of this signal is around 200 μm. Because the redshift distribution of these dusty star-forming galaxies is relatively broad, the signal also changes as a function of the observed wavelength (Béthermin et al. 2012). For this reason, using the column density map is not appropriate for analysing the CIB signature we detected in the Polaris flare region. A second segmentation analysis has therefore been done on the Herschel 250 μm map alone. Accurate CIB measurements using IRAS, Planck, and Herschel observations have already been done (Planck Collaboration XVIII 2011; Pénin et al. 2012; Viero et al. 2013a,b). Our analysis on the 250 μm map will be compared with the Viero et al. (2013a) results on the Multi-tiered Extragalactic Survey (HerMES; Oliver et al. 2012).

The segmented (Gaussian and coherent) power spectra for the Herschel 250 μm map are presented in Fig. 12. The spectra were calculated in the same way as for the column density map. The fit values are listed in Table 3. The general shapes of the spectra are similar to those calculated for the column density map in Fig. 8. The fit values for the total Fourier power spectrum corresponds within uncertainties with the values fitted by Miville-Deschênes et al. (2010) for the Herschel250 μm map of Polaris covering a larger and slightly different field of view. Contrary to the column density map, the power for the Gaussian power spectrum decreases quickly for k ≳ 1 arcmin−1. This difference can be attributed to the smaller beam pixel size at 250 μm.

Viero et al. (2013a) measured (from the combination of different Herschel fields and extended sources masked at k = 1.406 arcmin−1) a power of (8.54 ± 0.31) × 103 Jy2 sr−1. Their spectra are corrected for a cirrus power law fixed to γ = 3.0. Here, we considered the power law of the Gaussian part only as the cirrus signal. Following Eq. (25), we find a power of (7.6 ± 0.7) × 103 Jy2 sr−1. This value of the CIB power corresponds very well to that evaluated previously by Viero et al. (2013a). It is important to recall that in our case the CIB signal was dominated by the foreground emission of Polaris. As we can see in Fig. 12, the small-scale CIB power was dominated by the power of coherent structures and the MnGSeg method succeeded nonetheless to extract the mean power of this relatively faint signal over the map. Figure 13 shows on a subregion of the Polaris flare that the spatial extraction of the CIB signal is also possible using the MnGSeg decomposition. A more in-depth analysis of the CIB signal is beyond the scope of this paper, but the MnGSeg method is a good strategy for this analysis.

thumbnail Fig. 12

Segmented power spectra for the Polaris flare map at 250 μm. The total Fourier power spectrum is represented by the solid black line. The red diamonds show the power spectrum for the Gaussian self-similar part of the map. The blue triangles show the power spectrum for the non-Gaussian coherent part of the map. The black dashed line represents the fits of the curves and the green dashed line represents the CIB and sources contribution level. The fit values are listed in Table 3.

Open with DEXTER
Table 3

Fitvalues for the total 250 μm map, the Gaussian, and the non-Gaussian power spectra.

4.6 Ratio coherent/Gaussian

As shown in Fig. 7 and described in Sect. 4.3, the skewness of the intermittency PDFs increases exponentially towards smaller scales until the scales become dominated by the noise level. Another way to look at the intermittency measure is to compare the two components by calculating the ratio of the coherent power spectrum to the Gaussian power spectrum. The ratios for the column density and the 250 μm maps are plotted in Fig. 14. The ratios are corrected for the noise and CIB/sources contributions and, according to Eqs. (25) and (26), it is defined as (27)

For both maps, the power spectra ratio have a power law shape, as defined in Eq. (27), with a bump centred at k ≈ 0.05 arcmin−1, which corresponds to ~0.15 pc. This bump was also seen in Fig. 7 for the skewness as a function of scales. The dashed curve in Fig. 14 shows the exponential curve as defined in Eq. (27) using the fitted power-law values for the column density mapin Table 2. The dot-dashed curve shows the exponential for the power spectra ratio fitted for 0.077 ≲ k ≲ 0.24 arcmin−1 and has a power law of 1.2 ± 0.1. Assuming that the power-law fit shown in Fig. 8 were also affected by this bump around k ≈ 0.05 arcmin−1, the corrected Gaussian power law would become γr + γC = 3.6 ± 0.1, which corresponds to the Kolmogorov power law for a 3D incompressible turbulent medium.

This excess of power is present in both the column density map and the 250 μm map, which confirms that the excess is real and not an artefact from the column density calculation. Figure 15 shows the Gaussian and coherent reconstructed maps for 0.025 ≲ k ≲ 0.077 arcmin−1, which corresponds to 0.11 ≲ l ≲ 0.33 pc. As expected, the Gaussian map shows density fluctuations that are more uniform than the coherent map. For this range of spatial scales, the excess of power can be attributed to contrasted structures that can also easily be identified in Fig. 3, as the saxophone in the south-west region of the field. However, the excess of power in the power spectra ratio is also present when the saxophone, the brightest, high column density region is excluded. This means that more quiescent regions of the Polaris Flare, like the one defined by Schneider et al. (2013), also exhibits this characteristic scale.

thumbnail Fig. 13

Different segmentations applied on a subregion of the Polaris map at 250 μm. (a) Original map. (b) Gaussian map reconstructed from all spatial scales. (c) Coherent map reconstructed from all spatial frequencies. (d) Original map without the frequencies associated with the flattened and decreasing part of the Gaussian power spectrum related to the CIB (i.e. 0.9 ≲ k≲ 2.2 arcmin−1). (e) Gaussianmap without frequencies 0.9 ≲ k ≲ 2.2 arcmin−1. (f) Sum of the Gaussian frequencies 0.9 ≲ k ≲ 2.2 arcmin−1 dominated by the CIB signal.

Open with DEXTER
thumbnail Fig. 14

Ratio of the coherent to the Gaussian power spectra for the column density and the 250 μm maps. The dashed curve shows the exponential curve defined in Eq. 27 using the fitted power-law values for the column density map in Table 2. The dot-dashed curve shows the exponential for the power spectra ratio fitted for 0.08 ≲ k ≲ 0.2 arcmin−1, which corresponds to 0.03 ≲ l ≲ 0.1 pc.

Open with DEXTER
thumbnail Fig. 15

Gaussian and coherent reconstructed maps for 0.025 ≲ k ≲ 0.077 arcmin−1, which corresponds to 0.11 ≲ l ≲ 0.33 pc.

Open with DEXTER
thumbnail Fig. 16

Anisotropy measure as a function of the azimuthal angle θ and the wavenumber k as defined in Eq. (28).

Open with DEXTER

4.7 Anisotropies

Since the non-Gaussian segmentation is applied as a function of azimuthal angles, the segmentation should also reveal structure anisotropies as a function of spatial scales. Figure 16 shows the anisotropy measure as suggested by Farge (1992) for the totalset of wavelet coefficients, the Gaussian part and the coherent part. The anisotropy measure is defined as (28)

According to Fig. 16, non-Gaussianities start to appear in most directions at k ≈ 0.018 arcmin−1. Larger anisotropies are measured onlarge scales around 45° in the total set of coefficients and in the Gaussian part, whereas they are not detected as non-Gaussianities. It is important to remember that these large scales, the number of independent pixels is n ≈3, thus the sample for detecting outliers is very small. The next anisotropies are measured at k ≈ 0.033 arcmin−1, where n ≳20 in the Gaussian and in the coherent part of the field with a value of A(l, θ) ≈ 2.5, which is much lower than the intermittency measures shown in Fig. 5. Anisotropies for k ≳ 0.06 arcmin−1 fluctuate slightly around 1.0, which means that in the inertial range structures have no particular direction. This result is also in agreement with the Polaris filamentary structure analysis of Ossenkopf-Okada & Stepanov (2019). On the smallest scale, anisotropies are due to the pixelation effect.

5 Discussion

The MnGSeg method, via simple assumptions such as the self-similar nature of incompressible turbulence and the ergodicity inferred by the Fourier power spectrum analysis, is able to extract two fundamentally different statistical behaviours in molecular cloud gas density distributions. This dual nature of molecular clouds corresponds well to the previous description given by Falgarone et al. (2004). One fractal component, characterised by its self-similarity over the spatial scales and its unique power law (see Figs. 8 and 10), corresponds to the diffuse component of the Polaris flare. The other component, called the coherentcomponent, in contrast to the fractal component, displays a network of filaments of different sizes. These coherent filamentary structures are physically important because they are where dense cores are embedded and thus crucial for star formation. Ward-Thompson et al. (2010) identified the five densest cores within the MCLD 123.5+24.9 highest density structure (also called the saxophone). These starless cores were subsequently studied for their chemical composition and their mass by Shimoikura et al. (2012) and Wagle et al. (2015). Future analyses of the comparison between the physical properties and the spatial distribution of dense cores or young stellar objects and the coherent ISM structures will be important in order to shed the light of the impact of the ISM environment on the star formation processes.

The fact that filamentary structures affect the power spectrum analysis of a region was already found by Schneider et al. (2011) and Elia et al. (2014). However, in the present analysis, we demonstrate that the fluctuations associated with the coherent filamentary component of the cloud dominate over most of what is considered the inertial range of the power spectrum. This is demonstrated in the present analysis on the Polaris flare and was demonstrated previously on a larger area of the Galactic plane (Robitaille et al. 2014). This result challenges the recent study of Roy et al. (2019) suggesting that only filamentary structures with a sufficiently high area filling factor and/or a high column density contrasts have an impact on the scale-free power spectrum of dust continuum images. In the light of the present results, the important question to ask is not why filamentary structures do not produce any breaks in a power spectrum analysis, but how the coherent filamentary structures we extracted with the MnGSeg technique are related to the true scale-free component of molecular clouds.

The segmentation procedure presented in this paper corresponds well to the bifractal intermittency model described in the review on turbulence by Frisch (1995). From a statistical and geometrical point of view, the Gaussian “incompressible” component extracted using the MnGSeg technique is monofractal by construction, which means that a single power law for the inertial range is sufficient to describe the dynamics of this component. In the bifractal model described by Frisch (1995), the inertial range of an intermittent turbulent flow is in fact the result of a competition between two power laws, where the one with the smallest exponent dominates on small scales. In the case where the steeper power corresponds to the Kolmogorov incompressible turbulence, only structure functions of order higher than three would be affected by the second power law. This means that the Fourier power spectrum is not sensitive to the bifractal model of intermittency. In our case we show that the Gaussian part of the density map has a power law closer to the Kolmogorov incompressible turbulence; however, the coherent part of the density map still has an influence on the total Fourier power spectrum.

Frisch (1995) also proposed a more complete multifractal model of the intermittence. There is no reason to believe that the coherent component of the density map is monofractal. Theories on fully developed turbulence predict that local variations of the turbulent energy cascade can be defined as a multifractal system (Frisch 1995; Arnéodo et al. 1995). A full multifractal analysis of the Polaris density and velocity field is beyond the scope of this paper. However, the extraction of a monofractal component covering all spatial scales and filling the whole field is a first step towards a better comprehension of how turbulence is related to the ISM gas structures. An important comparaison of multifractal analysis based on a box-counting approach has been recently done by Elia et al. (2018) using Hi-GAL observations (a key programme of Herschel), fBm images, and numerical simulations. They found that all the investigated fields, which are located in the Galactic plane, exhibit a multifractal structure rather than a simple monofractal structure. According to this result, numerical simulations appear, depending on the specific model, more similar to observations than fBms. Multifractal analysis methods of images or signals based on wavelet transforms are already developed and are actively used in multiple domains from turbulence analysis to medical research (Arnéodo et al. 2000), and the MnGSeg analysis could be potentially adapted to these approaches.

It is important to remember that the Kolmogorov energy cascade was derived from the velocity structure function and that our analysis is based on the density power spectrum. The link between the density and velocity power spectrum has been investigated mainly in the case of the warm ionised medium (Armstrong et al. 1995; Terry & Smith 2007), where a power law of γ = 11 ∕ 3 for the electron density power spectrum has been measured over a broad range of spatial scales. In hydrodynamics and magnetohydrodynamics, equations show that fluctuations in turbulence velocity, magnetic field, and density are coupled. For subsonic to transonic turbulence, the velocity power spectrum depends slightly on the Mach number (Kritsuk et al. 2007). At higher Mach number, density fluctuations are more sensitive mainly due to the apparition of the compressive component of the turbulence, which is measured by the shallower energy power spectrum. The correlation between the density and velocitypower spectra has also been investigated through the velocity-channel analysis (VCA; Lazarian & Pogosyan 2000), which is based on the variations of the power spectra in velocity channels at changing velocity resolution.

As mentioned in Sect. 2, the scale-free nature of ISM power spectra is commonly associated with the energy cascade of the Kolmogorov theory of turbulence, which predicts a power law of γ = 11∕3 for a subsonic and incompressible flow. This value is rarely found in the ISM and the typical value is closer to that measured for the total Fourier power spectrum of the Polaris 250 μm map: γ ≈ 2.6 (Padoan et al. 1997a). The shallower slope is usually attributed to the presence of small-scale structures present in the compressible ISM. Our analysis has demonstrated that the non-Gaussianities associated with the coherent structures in the Polaris molecular cloud are present on a wide range of spatial scales and not only at small scales. These non-Gaussianities appear progressively towards the small scales along the inertial range of the fractal diffuse component of the cloud. The Fourier power spectrum being equivalent only to the second-order structure function, it is essentially not sensitive to these excesses of power as a function of spatial scales, referring here to the non-Gaussianities shown in the intermittence measure of Figs. 5 and 6. The local reconstruction of the fluctuations responsible for these excesses of power as a function of spatial scales shows that the non-Gaussianities are associated with the filamentary network of molecularclouds (see right panel of Fig. 9). Following our algorithm, the non-Gaussianities of the density distributions have been identified by calculating the third moment of the distributions (i.e. the skewness).

In the more recent literature, an unclear partition is made between the PDF analysis of density and velocity fields. While the global PDFs of star-forming regions are generally applied to column density maps (Federrath et al. 2010; Schneider et al. 2013), multiscale PDF analysis in order to detect intermittency is generally reserved to velocity field or centroid velocity increments (Hily-Blant et al. 2008; Bertram et al. 2015). For the turbulent ISM, it is generally assumed that the density PDF follows a lognormal distribution, where the standard deviation σ is related to the Mach number through the relation (29)

with b ≈ 0.5 determined by Padoan et al. (1997a) with supersonic magnetohydrodynamic (MHD) simulations. On the other hand, as shown by Lis et al. (1996), some normalised PDFs of centroid velocity increments have broad to near-exponential wings in PDFs, while global PDFs for the entire velocity field are approximately Gaussian. According to these observations, should the global PDF of turbulent density or velocity fluctuations still be considered a good measure of the turbulence regime? Orkisz et al. (2017) showed in observations, by applying a reconstruction method of the statistical properties of a vector velocity field (Brunt & Federrath 2014), that there can be a strong intra-cloud variability of the compressive and solenoidal fractions. In this case, a more local or scale-dependent approach for PDF analyses would more suitable.

In this paper we propose an alternative approach to identify the intermittent behaviour of the ISM, where the non-Gaussianities are identified by calculating the third moment of fluctuation distributions as a function of scales. One powerful aspect of our procedure is the use of complex-wavelet transforms, which, contrary to the structure function or Fourier power spectrum analysis, allows us to truly filter the spatial scales and extract the non-Gaussian component contributing normally to the shallower slope than ζp = p∕ 3 for structure functions with p ≥ 3. Another advantage of the MnGSeg analysis is that it allows us to calculate a non-biased power spectrum of a map, i.e. a power spectrum analysis not affected by the non-Gaussianities. For the column density map and for the intensity map at 250 μm, the Gaussian power spectrum has a power law closer to the Kolmogorov incompressible turbulence of 11 ∕ 3. This result could be interpreted, with caution, as a demonstration of the correlation between fluctuations in velocity and density and it shows that the segmentation procedure could extract the incompressible turbulent component. Further analyses applying the MnGSeg technique on centroid velocity maps and numerical simulations, where the intermittent behaviour of both the density and velocity field could be compared, are needed to confirm this interpretation.

We can also interpret from these results that the ~0.15 pc characteristic scale identified in the ratio of the coherent to the Gaussian power spectra in Fig. 14 is associated with the transition of regimes dominated by incompressible turbulence versus compressible modes and by other physical processes such as gravity. However, it has been seen in laboratory experiments that incompressible turbulence alone can also produce intermittence. For this reason, we cannot directly attribute the extracted non-Gaussianities to compressive modes only. Schneider et al. (2011) observed many characteristic scales by applying the Δ-variance method over several Galactic clouds. For a 13CO line-integrated map of Cygnus X they revealed two peaks, one at 4 pc and another one at 80 pc. For all the low-mass star formation regions they found a double-peak structure in their Δ-variance spectrum of their extinction map. These regions are localised around 0.4–1.5 pc and 2.9–4.6 pc. In contrast, they found no characteristic scale for Polaris. They stated three main possible causes for these Δ-variance peaks: preferred geometric scales, such as length and width of filaments or clumps; the decaying turbulence scales (Mac Low & Ossenkopf 2000); or scale energy injection from external or internal sources.

The origin of the characteristic scale identified with the power spectra ratio is thus uncertain. Many mechanisms and other forces in addition to the turbulence play a role in the structure formation in the ISM. Among the forces, gravity itself undoubtedly plays an important role. Coherent structures are also seen in velocity fields. Observations (Galván-Madrid et al. 2010; Hennemann et al. 2012; Hacar et al. 2017, 2018) and simulations (Smith et al. 2016; Zamora-Avilés et al. 2017) show that large filamentary structures represent a collection of smaller and coherent subfilaments, sometimes called fibres. According to the simulations, the turbulence being responsible initially for the large-scale density fluctuations, filaments and their substructures formation would then be gravity-driven by accretion. This scenario where dense substructures are no longer linked to the turbulence could be seen to contradict the results of this paper, where the hierarchical nature of the coherent structures seems intimately linked with the energy cascade of the turbulence. This could suggest, however, that coherent density-driven filaments originate from shocks directly associated with compressive turbulence. Our results are in agreement with the gravoturbulent models of Larson (1981), Hennebelle & Chabrier (2008), and more recently Lee et al. (2017) which take into account the scale dependance of the supporting thermal, turbulent, and magnetic energies. The magnetic field certainly also plays a role in the formation of large filaments, where it has been found oriented perpendicularly to nearby filaments, such as Musca and Taurus (Palmeirim et al. 2013; Planck Collaboration XXXIII 2016). In the case of the Polaris flare region, the projected magnetic field was found to be preferentially aligned with the dust filamentary structures (Panopoulou et al. 2015, 2016). However, even the denser filamentary structure of the Polaris cloud has a lower column density than Musca and Taurus, and in contrast with these two filaments, the Polaris flare does not show any sign of star-formation activity. These peculiar properties of Polaris suggest that the cloud is engaged in its initial phases of molecular cloud formation and it could also have an impact on the magnetic field configuration.

6 Conclusion

The MnGSeg technique is a new powerful analysis method which constitutes a robust alternative to the classical Fourier power spectrum. By coupling the multiscaled PDFs with the power spectrum analysis, this novel technique is sensitive to the progressive contribution of non-Gaussianities towards the small spatial scales. These contributions commonly attributed to turbulence intermittency are usually measured only with higher-order structure function (p ≥ 4). The great advantage of the MnGSeg technique over the Fourier power spectrum or the structure function is that it can easily expose the relationship between the self-similar Gaussian structures and the progressive contribution of non-Gaussianities, and that it allows the spatial reconstruction of both components.

Using the MnGSeg technique, we prove that the Fourier power spectrum of the Polaris flare is dominated by the power of its denser filamentary structures. The spatial combination of these non-Gaussianities with the fractal scale-free component produces no characteristic scale visible in the Fourier power spectrum. The origin of these non-Gaussianities appearing as a dense filamentary structure is likely diverse. Our results are in agreement with a global emerging scenario, also seen through numerical simulations, where turbulence plays a dominant role in the early stages of molecular cloud formation. Then the interplay of turbulence intermittency, gravity, and other mechanisms such as thermal instability, shocks, and the influence of magnetic fields are quickly responsible for the density enhancement of the medium and ultimately the formation of stellar objects. In this context the in-depth analysis of this transition between the two regimes, the incompressible random turbulent field (demonstrating a monofractal nature), and the coherent structures (most likely multifractal and intimately linked to the formation of stars) is fundamental for our global understanding of molecular cloud formation and star formation in the ISM. We tentatively measured this transition on a ~ 0.15 pc scale, but further tests need to be performed to confirm the association between this characteristic scale and the turbulent regime transition.

The MnGSeg technique deserves to be applied also to velocity or centroid velocity maps. However, the direct application to velocity quantities faces three main difficulties: the lower spatial resolution of radio line maps; the multiple velocity components along the line of sight, especially for dense filamentary structures; and the multiple tracers needed to observe the different gas phases and densities.

The future comparison of MnGSeg analyses of multiple column density maps of nearby molecular clouds will allow us to understand the link between the formation of coherent structures in the ISM and the emergence of star formation activity. The application of a fully multifractal analysis of the coherent structures will also help to expose more accurately the plural nature of molecular clouds and to validate whether the spatial variations of the turbulence dissipation is ultimately linked to the core mass function of star-forming regions.

Acknowledgements

We would like would like to thank Isabelle Joncour for the precious comments on this paper. This projecthas received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie Grant Agreement No. 750920. N.S. and S.B. acknowledge support by the french ANR and the german DFG through the project GENESIS (ANR-16-CE92-0035-01/DFG1591/2-1). This research has made use of data from the Herschel Gould Belt survey (HGBS) project ([Article]). The HGBS is a Herschel Key Programme jointly carried out by SPIRE Specialist Astronomy Group 3 (SAG 3), scientists of several institutes in the PACS Consortium (CEA Saclay, INAF-IFSI Rome and INAF-Arcetri, KU Leuven, MPIA Heidelberg), and scientists of the Herschel Science Center (HSC). This research has also made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration 2013) and its affiliated package APLpy (Robitaille & Bressert 2012).

References


1

The codes and tutorials applied on mock simulations are available at [Article]

2

Throughout the paper, bold variables denote vectorial quantities. For simplicity, when quantities are averaged over azimuthal angles, as in Eq. (2), the variable is then considered scalar.

3

The admissibility condition requires the zero mean value of the wavelet function, . Since without any correction for the Morlet wavelet (see Eq. (9)), this wavelet is only marginally admissible.

All Tables

Table 1

Fit values for the column density map Fourier power spectrum.

Table 2

Fit values for the column density map Gaussian and the non-Gaussian power spectra.

Table 3

Fitvalues for the total 250 μm map, the Gaussian, and the non-Gaussian power spectra.

All Figures

thumbnail Fig. 1

Schematic representation of 2D Fourier space, where u and v are the two dimensions, k is the wavenumber, and θ is the azimuthal angle.

Open with DEXTER
In the text
thumbnail Fig. 2

Top left panel: real part of the Morlet wavelet. Top right panel: Morlet wavelet rotated over an azimuthal angle of π, also called the Fan wavelet (Kirby 2005). Bottom left panel: Morlet wavelet in Fourier space. Bottom right panel: Morlet wavelet rotated over an azimuthal angle of π in Fourier space.

Open with DEXTER
In the text
thumbnail Fig. 3

Polaris column density map derived from Herschel imaging data taken as part of the HGBS and “Evolution of interstellar dust” key programs. The brightest structure with the highest column density value located on the south-western part of the field has been labelled as MCLD 123.5+24.9 and as the “saxophone” by Schneider et al. (2013).

Open with DEXTER
In the text
thumbnail Fig. 4

Fourier (solid lines) and wavelet (diamonds) power spectra of the Polaris flare region presented in Fig. 3. The black line and black symbols represent the original power spectrum, P(k) in Eq. (21), while the blue line and blue symbols represent Psky (k) in Eq. (22).

Open with DEXTER
In the text
thumbnail Fig. 5

Intermittency measure I(l, x) as a functionof the spatial frequency k. Frequencies from ~0.01 to ~ 0.75 arcmin−1 are plotted, and approximately correspond to scales where enough statistically independent points are available on large scales and where the signal is not dominated by noise on small scales.

Open with DEXTER
In the text
thumbnail Fig. 6

Intermittency measure in linear scale for three spatial frequencies. Compared to Fig. 5, the number of bins for each spatial scale respects the relation n = (Nx × Ny) ∕ l2, where k = |k0| ∕ l.

Open with DEXTER
In the text
thumbnail Fig. 7

Skewness, ς(l), for the intermittency measures defined in Eq. (23) and shown Figs. 5 and 6.

Open with DEXTER
In the text
thumbnail Fig. 8

Segmented power spectra for the Polaris flare column density map. The total Fourier power spectrum shown in Fig. 4 is represented by the solid black line. The red diamonds show the power spectrum for the Gaussian self-similar part of the map. The blue triangles show the power spectrum for the non-Gaussian coherent part of the map. The black dashed line represents the fits of the curves and the green dashed line represents the CIB and contribution level of the sources. The 3D Kolmogorov power law of − 11∕3 is plotted with the red dash-dotted line.

Open with DEXTER
In the text
thumbnail Fig. 9

Gaussian (left) and coherent (right) reconstructed column density maps following Eqs. (19) and (20); the correction factor Cδ has been multiplied to both maps, and the mean values of the initial map μ0 has been distributed to both maps.

Open with DEXTER
In the text
thumbnail Fig. 10

Top left panel: Gaussian reconstructed map filtered for scales k ≳ 0.3 arcmin−1. Top right panel: fractional Brownian motion simulation with the same power law as the Gaussian part of Polaris. The fBm is also filtered for k ≳ 0.3 arcmin−1. Bottom left panel: intermittency measure as defined in Eq. (23) for the Gaussian segmentation of wavelet coefficients. Bottom right panel: intermittency measure of the fBm simulation on the same spatial scales as Polaris.

Open with DEXTER
In the text
thumbnail Fig. 11

Intermittency measure as defined in Eq. (23) for the non-Gaussian wavelet coefficients.

Open with DEXTER
In the text
thumbnail Fig. 12

Segmented power spectra for the Polaris flare map at 250 μm. The total Fourier power spectrum is represented by the solid black line. The red diamonds show the power spectrum for the Gaussian self-similar part of the map. The blue triangles show the power spectrum for the non-Gaussian coherent part of the map. The black dashed line represents the fits of the curves and the green dashed line represents the CIB and sources contribution level. The fit values are listed in Table 3.

Open with DEXTER
In the text
thumbnail Fig. 13

Different segmentations applied on a subregion of the Polaris map at 250 μm. (a) Original map. (b) Gaussian map reconstructed from all spatial scales. (c) Coherent map reconstructed from all spatial frequencies. (d) Original map without the frequencies associated with the flattened and decreasing part of the Gaussian power spectrum related to the CIB (i.e. 0.9 ≲ k≲ 2.2 arcmin−1). (e) Gaussianmap without frequencies 0.9 ≲ k ≲ 2.2 arcmin−1. (f) Sum of the Gaussian frequencies 0.9 ≲ k ≲ 2.2 arcmin−1 dominated by the CIB signal.

Open with DEXTER
In the text
thumbnail Fig. 14

Ratio of the coherent to the Gaussian power spectra for the column density and the 250 μm maps. The dashed curve shows the exponential curve defined in Eq. 27 using the fitted power-law values for the column density map in Table 2. The dot-dashed curve shows the exponential for the power spectra ratio fitted for 0.08 ≲ k ≲ 0.2 arcmin−1, which corresponds to 0.03 ≲ l ≲ 0.1 pc.

Open with DEXTER
In the text
thumbnail Fig. 15

Gaussian and coherent reconstructed maps for 0.025 ≲ k ≲ 0.077 arcmin−1, which corresponds to 0.11 ≲ l ≲ 0.33 pc.

Open with DEXTER
In the text
thumbnail Fig. 16

Anisotropy measure as a function of the azimuthal angle θ and the wavenumber k as defined in Eq. (28).

Open with DEXTER
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.