Effects of dust absorption on spectroscopic studies of turbulence

TL;DR

This paper investigates how dust absorption influences the analysis of magnetohydrodynamic turbulence spectra and anisotropies, revealing new regimes and resolving previous degeneracies, with implications for interpreting astrophysical observations.

Contribution

It introduces a new understanding of dust absorption effects on turbulence spectra, including a novel asymptotic regime and the removal of spectral degeneracy, enhancing analysis methods like VCA and VCS.

Findings

Dust absorption limits the volume to optical depth unity.

A new asymptotic regime causes the spectrum to become shallower by one.

Dust absorption removes the universal $K^{-3}$ spectrum degeneracy.

Abstract

We study the effect of dust absorption on the recovery velocity and density spectra as well as on the anisotropies of magnetohydrodynamic turbulence using the Velocity Channel Analysis (VCA), Velocity Coordinate Spectrum (VCS) and Velocity Centroids. The dust limits volume up to an optical depth of unity. We show that in the case of the emissivity proportional to the density of emitters, the effects of random density get suppressed for strong dust absorption intensity variations arise from the velocity fluctuations only. However, for the emissivity proportional to squared density, both density and velocity fluctuations affect the observed intensities. We predict a new asymptotic regime for the spectrum of fluctuations for large scales exceeding the physical depths to unit optical depth. The spectrum gets shallower by unity in this regime. In addition, the dust absorption removes the degeneracy resulted in the universal $K^{-3}$ spectrum of intensity fluctuations of self-absorbing medium reported by Lazarian \& Pogosyan. We show that the predicted result is consistent with the available HII region emission data. We find that for sub-Alfv\'enic and trans-Alfv\'enic turbulence one can get the information about both the magnetic field direction and the fundamental Alfv\'en, fast and slow modes that constitute MHD turbulence.

Tables

Table 1. Table 1: Summary of scalings of various statistical techniques in the presence of dust absorption.

Technique	Scale	Tool	Scaling	Equation
thin slice VCA		structure function	$R^{1-\nu /2}$	Eq. (47)
	$R>{\Delta }_0$	correlation function	$R^{-\nu /2}$	Eq. (54)
thick slice (steep density)		structure function	$R^{1+{\nu }_{\rho }}$	Eq. (41)
	$R>R_c$	structure function	saturated	Eq. (41)
thick slice (shallow density)		correlation function	saturated	Eq. (37)
	$R>R_c$	correlation function	$R^{1-{\nu }_{\rho }}$	Eq. (37)
centroids		structure function	$R^{1+\nu }$	Eq. (82)
	$R>{\Delta }_0$	structure function	$R^{\nu }$	Eq. (81)
VCS	$k_v^{-1}>V_{\Delta B}$	power spectrum	$k_v^{-2\left(1-{\nu }_{\rho }\right)/\nu }$	Eq. (66)
		power spectrum	$k_v^{-2\left(3-{\nu }_{\rho }\right)/\nu }$	Eq. (67)

Equations

$$\mathrm{d}{I}_{\nu }=-g_{\nu }{I}_{\nu }\mathrm{d}z+j_{\nu }\mathrm{d}z-{\sigma }_{\nu }{I}_{\nu }\mathrm{d}z,$$

$$g_{\nu }$$

$$j_{\nu }$$

$${\sigma }_{\nu }$$

$$I_{\nu }=$$

$$\times \exp \left[-{\int }_0^z\mathrm{d}{z}^{\prime }\left(\alpha {\rho }_g(\mathbf{X},z^{\prime }){\Phi }_{\nu }(\mathbf{X},z^{\prime })+{\kappa }^{\prime }{\rho }_d(\mathbf{X},z^{\prime })\right)\right].$$

$$I_{v_1}(\mathbf{X})=\frac{ϵ}{\alpha }\left[1-{\mathrm{e}}^{-\alpha {\rho }_s(\mathbf{X},v_1)}\right],$$

$${\rho }_s(\mathbf{X},v_1)={\int }_0^S\mathrm{dz}\rho (\mathbf{x}){\Phi }_v(\mathbf{x}).$$

$$I_{\nu }(\mathbf{X})=ϵ{\int }_0^S\mathrm{dz}{\rho }_g(\mathbf{x}){\Phi }_v(\mathbf{x})\exp \left[-{\int }_0^z\mathrm{d}z{\kappa }^{\prime }{\rho }_d(\mathbf{X},z^{\prime })\right],$$

$${\Phi }_v(\mathbf{x})=\frac{1}{\sqrt{2\pi \beta }}\exp \left[-\frac{(v-u(\mathbf{x})-v_c(\mathbf{x}){)}^2}{2\beta }\right],$$

$${\kappa }^{\prime }{\rho }_d(\mathbf{x})=\kappa {\rho }_g(\mathbf{x})\equiv \kappa \rho (\mathbf{x}),$$

$${\xi }_{\mathcal{R}}({\mathbf{x}}_1,{\mathbf{x}}_2)=⟨\mathcal{R}({\mathbf{x}}_1)\mathcal{R}({\mathbf{x}}_2)⟩,$$

$$D_{\mathcal{R}}({\mathbf{x}}_1,{\mathbf{x}}_2)=⟨{\left(\mathcal{R}({\mathbf{x}}_1)-\mathcal{R}({\mathbf{x}}_2)\right)}^2⟩.$$

$$D_z(\mathbf{r})=D_z(S){\left(\frac{r}{S}\right)}^{\nu },r<S$$

$$D_{\rho }(\mathbf{r})\approx 2{\sigma }_{\rho }^2{\left(\frac{r_c}{r}\right)}^{-{\nu }_{\rho }},{\nu }_{\rho }<0,r<r_c,$$

$${\xi }_{\rho }(\mathbf{r})\approx {\rho }_0^2+{\sigma }_{\rho }^2{\left(\frac{r_c}{r}\right)}^{{\nu }_{\rho }},{\nu }_{\rho }>0,r>r_c,$$

$$I_{v_i}({\mathbf{X}}_1;\Delta V_i)={\int }_{-\infty }^{\infty }d{v}_1{I}_{v_1}(\mathbf{X})W(v_1-v_i;\Delta V_i)$$

$$\xi (\mathbf{R},v_i,\Delta V)=⟨I_{v_i}({\mathbf{X}}_1;\Delta V)I_{v_i}({\mathbf{X}}_2,\Delta V)⟩$$

$$\mathcal{D}(\mathbf{R},v_i,\Delta V)=⟨{\left(I_{v_i}({\mathbf{X}}_1,\Delta V)-I_{v_i}({\mathbf{X}}_2,\Delta V)\right)}^2⟩.$$

$$\xi (\mathbf{R})\sim {\int }_{-\infty }^{\infty }d{v}_1{\int }_{-\infty }^{\infty }d{v}_2⟨I_{v_1}({\mathbf{X}}_1)I_{v_2}({\mathbf{X}}_2)⟩$$

$$\xi (\mathbf{R})\sim \Delta V^2⟨I({\mathbf{X}}_1)I({\mathbf{X}}_2)⟩$$

$$\tau (\mathbf{X},z)\equiv \kappa {\int }_0^z\rho (\mathbf{X},z^{\prime })d{z}^{\prime },$$

$$\kappa {\int }_0^{\Delta (\mathbf{X})}\mathrm{d}z\rho (\mathbf{X},z)=1,\tau (\mathbf{X},\Delta (\mathbf{X}))\equiv 1,$$

$$I(\mathbf{X})\propto 1-{\mathrm{e}}^{-\tau (\mathbf{X},S)}.$$

$$e^{-\tau (\mathbf{X},z)}\approx \Theta \left(\Delta (\mathbf{X})-z\right),$$

$$I\approx ϵ{\int }_0^{\Delta }\mathrm{d}z\rho (\mathbf{x})=\frac{ϵ\tau (\mathbf{X},\Delta (\mathbf{X}))}{\kappa }\equiv \frac{ϵ}{\kappa }$$

$$I_v(\mathbf{X})\approx {\int }_0^{\Delta (\mathbf{X})}dz\rho (\mathbf{X},z){\Phi }_v(\mathbf{X},z)$$

$$I_v(\mathbf{X})\sim \Delta (\mathbf{X})\rho (\mathbf{X},0){\Phi }_v(\mathbf{X},0)\approx \frac{1}{\kappa }{\Phi }_v(\mathbf{X},0).$$

$$\displaystyle I_{v}(\bm{X})=\frac{1}{\kappa}\bigg{(}\left[\Phi_{v}(\bm{X},0)-\Phi_{v}(\bm{X},S)\mathrm{e}^{-\tau(\bm{X},S)}\right]$$

$$\displaystyle+\int_{0}^{S}\mathrm{d}z\,\frac{\mathrm{d}\Phi_{v}(\bm{x})}{\mathrm{d}z}\mathrm{e}^{-\kappa\int_{0}^{z}\mathrm{d}z^{\prime}\,\rho(z^{\prime})}\bigg{)}~{}.$$

$$I_v(\mathbf{X})\approx \frac{1}{\kappa }\left[{\Phi }_v(\mathbf{X},\Delta (\mathbf{X}))-{\Phi }_v(\mathbf{X},S){\mathrm{e}}^{-\tau (\mathbf{X},\mathbf{S})}\right].$$

$$I_v(\mathbf{X})\approx \frac{1}{\kappa }{\Phi }_v(\mathbf{X},\Delta (\mathbf{X})),\tau (\mathbf{X},S)>1$$

Full text

††pubyear: 2017††pagerange: Effects of dust absorption on spectroscopic studies of turbulence–C

1 Introduction

It is now an accepted fact that the interstellar medium (ISM) is magnetized and turbulent (see Elmegreen & Scalo 2004; McKee & Ostriker 2007; Mac Low & Klessen 2004). The Reynold’s number in the ISM may exceed $10^{8}$; in such case flows are expected to be turbulent. Moreover, analysis of scattering and scintillations observations shows turbulent spectra of electron density ranges many orders of magnitude, the so-called Big Power Law in the sky (see Armstrong et al. 1995; Chepurnov & Lazarian 2010). The fluctuations of synchrotron emissivity and polarization also reflect the underlying magnetised turbulence in the ISM (see Lazarian & Pogosyan 2012, 2016 and references therein).

In this paper, we are interested in using non-thermally broadened spectral lines in order to study the underlying turbulence statistics. The corresponding analytical description of statistics of spectral line intensity is provided in a series of our papers (Lazarian & Pogosyan 2000, 2004, 2006, 2008; Kandel et al. 2016b, 2017). These papers describe various aspects of the statistics of the position-position-velocity (PPV) cubes that can be available from observations111We mention parenthetically that the PPV cubes do not need to be fully sampled. For instance, the studies of turbulence with emission and absorption lines are possible for just a few lines passing through the volume (see Chepurnov & Lazarian 2010)..

The VCA was introduced in Lazarian & Pogosyan (2000, hereafter LP00) and was extended for self-absorbing media in Lazarian & Pogosyan (2004, hereafter LP04). The analytical expressions obtained predict how the intensity fluctuations in the velocity slice of PPV depend on the underlying turbulent velocity and density spectra with the change of the velocity slice thickness. This allows one to obtain the velocity and density spectra separately.

While the VCA technique has been successfully applied to study turbulence spectrum in our galaxy and the Small Magellanic Cloud (SMC) (see LP00; Stanimirović & Lazarian 2001; Lazarian 2009 and references therein), the studies of anisotropies using VCA have so far only been performed with the synthetic data. When combined with the study of synchrotron fluctuations as suggested in Lazarian & Pogosyan (2012, 2016), VCA allows studies of the variation of turbulence properties from the cold to hot media.

Another nobel way of obtaining information on turbulence was identified in Lazarian & Pogosyan (2006, hereafter LP06). There it was suggested to study fluctuation of intensity along the velocity axis, and analytical relations between these fluctuations and the underlying turbulence velocity and density spectra were established. The technique was termed Velocity Coordinate Spectrum (VCS), and it has been successfully applied for studies of turbulence in molecular clouds (Padoan et al., 2009), for the HI turbulence at high galactic latitudes (Chepurnov et al., 2010), and for the HI in the SMC.

The Velocity Centroid technique is perhaps the oldest technique for studying turbulence (see Münch 1958). Using the synthetic data obtained with simulation of MHD turbulence, Esquivel & Lazarian (2005) showed that centroids are able recover spectra of subsonic turbulence as well as to trace anisotropies arising from the anisotropic nature of MHD turbulence. The direction of the anisotropy was shown to trace the direction of magnetic field. Esquivel & Lazarian (2011) have shown that the degree of anisotropy is related to the media magnetization that can be characterized by the Alfvén Mach number $M_{A}=V_{L}/V_{A}$, where $V_{L}$ is the turbulence velocity at the injection scale $L$, while $V_{A}$ is the Alfvén velocity. The analytical treatment of velocity centroids that capitalized on the aforementioned PPV studies was presented in Kandel et al. (2017, hereafter KLP17), where it was shown how to use velocity centroids to study the Alfvén, slow and fast MHD modes observationally.

The aforementioned techniques are complementary. The VCA was shown to be efficient for studies of supersonic turbulence, or subsonic turbulence using the heavier species of a turbulent flow. On the contrary, velocity centroids are shown to be much less affected by the thermal Doppler broadening and therefore applicable to the studies of subsonic turbulence using the main ingredients of the flow, e.g. atomic hydrogen. However, studies (Esquivel & Lazarian, 2005) suggest that centroids are not reliable for the studies of supersonic turbulence.

All the studies that were mentioned above do not account for the absorption by dust, however. This is acceptable for studies of turbulence using radio lines, e.g. 21 cm HI line, CO lines. However, dust absorption affects optical and UV lines. In fact, this was a problem that the earlier researchers were aware of (see Münch 1958). This paper extends the description of the PPV statistics by accounting for the effect of the dust absorption for the PPV studies and shows how the interpretation of results using VCA, the VCS and velocity centroid technique change in the presence of dust absorption.

In what follows, we discuss the radiative transport in the presence of dust in Sec. 2, and introduce the power law statistics of density and velocity in Sec 3. The effects of dust on the recovery of the velocity and density spectra using the VCA, VCS and velocity centroids are studied in Sec. 4, 5 and 6, respectively. The studies of anisotropies using both the VCA and the centroids are described in Sec. 7. In Sec. 8, applicability of VCA, VCS and centroids for collisionally excited lines is discussed. Comparison of theory with observational data is carried out in Sec. 9, while the discussion of the results and the summary are presented in Sec. 10 and Sec. 11, respectively.

2 Radiative transport in the presence of dust

In the presence of dust, the standard equation of radiative transfer is (Spitzer Jr, 2008)

[TABLE]

where the term $-\sigma_{\nu}I_{\nu}\mathrm{d}z$ denotes extinction due to dust. In the case of self-absorbing emission in spectral lines that is proportional to first power of density

[TABLE]

where $\rho_{g}$ is the density of emitting gas, $\rho_{d}$ is the density of dust, $\Phi_{\nu}(z)$ is the frequency distribution of emitters along the LOS, and $\kappa^{\prime}$ is the dust extinction coefficient, which is in general frequency dependent. For the purpose of our analysis, we assume it to be constant within the line width of a Doppler broadened emission line. The solution to the Eq. (1) takes the following form

[TABLE]

In the case when dust absorption is negligible, i.e. when $\kappa^{\prime}\rightarrow 0$, the spectral intensity reduces to

[TABLE]

where $\rho_{s}$ is the density of emitters in PPV space, and is given by

[TABLE]

In writing the above equations, we have used that fact that the Doppler-shifted frequency information can be used to deduct LOS velocity of the emitters, and hence all frequency dependence in radiative transfer equation can be equivalently designated as velocity dependence $\nu\rightarrow v$.

The study of turbulence using emission from self-absorbing medium was studied in detail in Lazarian & Pogosyan (2004) and KLP17. In this paper, we focus on the studies of turbulence using emission from a dusty medium.

In the case when self-absorption is negligible, i.e. $\alpha\rightarrow 0$, Eq. (2) reduces to

[TABLE]

where $\Phi_{v}(\bm{x})$ is given

[TABLE]

where $u(\bm{x})$ is the LOS component of turbulent velocity and $v_{c}(\bm{x})$ is the LOS component of coherent velocity that is not part of turbulent cascade. In this paper, we assume that $v_{c}(\bm{x})$ is small, else it has to be modelled independently. Clearly, the main effect of dust is it cuts out signals from optical depths larger than 1, thus attenuating intensity signals. In the subsequent analysis, unless explicitly stated, we will consider the case when self-absorption is negligible. Before starting further derivations, we assume that the dust density is directly proportional to the gas density222 This assumption has been shown to hold good in various papers (Lada et al., 1994)., so that

[TABLE]

where $\kappa$ is the redefined dust absorption coefficient, which takes into account the proportionality factor between density of gas and density of dust.

3 Velocity and density fields

In the first part of this paper, we are interested in the spectra of density and velocity fields, which are fluctuating in a turbulent media. Two important statistical measure of fluctuations $\delta\mathcal{R}(\bm{x})$ of a random quantity $\mathcal{R}(\bm{x})$ are correlation function

[TABLE]

and structure function

[TABLE]

We assume turbulence to be homogeneous in which case, the correlation function and structure function can be represented purely as a function of spatial separation, so that $\xi_{\mathcal{R}}(\bm{x}_{1},\bm{x}_{2})=\xi_{\mathcal{R}}(\bm{r})$ and $D_{\mathcal{R}}(\bm{x}_{1},\bm{x}_{2})=D_{\mathcal{R}}(\bm{r})$, where $\bm{r}=\bm{x}_{2}-\bm{x}_{1}$.

Scaling of power distribution of random fields are of two types: steep and shallow. For a steep field, power of fluctuations is dominated by moderate to large scale fluctuations, while for shallow field, power is dominated by small scale fluctuations. Depending upon whether fluctuations are steep or shallow, it is appropriate to either use structure function or use correlation function. A major difference between correlation function and structure function is that, while the structure function at scale $r$ is determined by the integrated power of fluctuations over scales smaller than $r$, the correlation function is determined by the integral of the power over scales larger than $r$. Therefore, the correlation function is more appropriate for shallow spectra, while structure function is more appropriate for steep spectra.

Observations show that velocity spectra in the ISM is steep, while density can be either steep or shallow. Leaving discussion of angular dependence of velocity and density statistics to Sec. 7, we take LOS velocity structure function as

[TABLE]

where $\nu>0$. In the above equation, $S$ is the injection scale of turbulence, and for isolated clouds, we assume that $S$ is comparable to the size of the cloud. Modelling of the structure function at $r>S$ requires modelling of the injection process. The correlation length of velocity field is comparable to the energy injection scale.

Similarly, the structure function of steep density field is taken to be

[TABLE]

while for shallow density spectrum, the density correlation is modelled as

[TABLE]

where $\rho_{0}$ is the mean density, $\sigma_{\rho}$ the density dispersion, and $r_{c}$ the correlation length of the density field. In the case steep density field, $r_{c}$ is comparable to the size of injection scale, i.e. $r_{c}\sim S$, and the power law scaling regime lies at $r<r_{c}$, while for shallow density, $r_{c}$ is small, and the power law scaling is exhibited at $r>r_{c}$. In the subsequent parts of this paper, we use Eq. (14) when $\nu_{\rho}<0$ and Eq. (15) when $\nu_{\rho}>0$.

Here, as in our previous papers (see LP00), velocity fluctuations are assumed to be Gaussian. In LP00, no assumptions on density field were made in the formulation of VCA for optically thin emission lines. However, as we shall see in the next section, one needs an explicit model of density field in order to carry out further analysis in the presence of dust absorption. In the case when density dispersion is smaller than the mean density (which is so for the steep density), it is safe to assume density fluctuations as Gaussian random field. However, if the density dispersion is larger than the mean density (which can be so for the shallow spectra), Gaussian distribution is not a good model for density field, as it violates an important constraint that $\rho>0$ (see Fry 1986). Thus, for a shallow spectra, we take density fluctuations to obey log-normal distribution. This distribution is physically reasonable as it always satisfies the constraint $\rho>0$, and is computationally convenient (Barrow & Morgan, 1983). This assumption also fits with recent observations, which suggest that density probability distribution function of the diffuse gas in the ISM is lognormal (see Berkhuijsen & Fletcher 2008).

4 Effects of dust on VCA

Velocity channel analysis (see LP00) concerns with studies of the correlations of line intensities measured in velocity channels

[TABLE]

Here $v_{i}$ is the central velocity (frequency) of the channel “i”, and $W(v_{1}-v_{i};\Delta V_{i})$ is the channel profile of the width $\Delta V_{i}$ (such as $W(v_{1}-v_{i};\Delta V_{i})\sim 0$ when $|v_{1}-v_{i}|>\Delta V_{i}/2$). VCA’s focus is on intensity correlations between pairs of LOS as function of LOS separation and, importantly, as a function of velocity channel width. Relevant statistical measures are the correlation function

[TABLE]

and the structure function

[TABLE]

There are two important regimes: thin velocity slice regime, where information about velocity spectra can be obtained, and thick velocity slice regime, where density spectra can be obtained. Which regime is realized depends on the comparison between the characteristic difference of turbulence velocities between two LOS and the velocity channel width. The slice is thin if $\Delta V<D_{z}^{1/2}(R)$ (LP00). Thus, to measure velocity effects in VCA, one either needs a narrow velocity channel or the measurements of correlations at sufficiently separated LOS. In practice, the minimal width of channels is determined by the spectral resolution of the instrument, while wider channels are obtained synthetically by coadding signal of individual high resolution channels (assuming for simplicity that channels are co-adjacent). This gives VCA the ability to separate velocity and density effects, if thin channel regime is available, by varying the slice thickness in the analysis. Here we note that thermal broadening interferes with the ability to carry out thin slice measurements by effectively increasing the channel thickness as $\Delta V^{2}\to\Delta V^{2}+\beta_{\text{T}}$. Thus, VCA is most effective for supersonic turbulence or when used with lines from massive species with slow thermal velocities.

For asymptotic analysis, thick slice regime corresponds formally to the channel profile $W\sim 1$, i.e. integration over the whole line. The explicit expression for correlation function is then

[TABLE]

Thin slice asymptotics is reproduced for $W(v_{1}-v_{i})\approx\Delta V_{i}\delta_{D}(v_{1}-v_{i})$ and

[TABLE]

Previous papers on VCA did not consider interstellar dust extinction. In this section, we consider the effect of dust on VCA statistics, and try to answer the following question: can velocity and density spectra be recovered using VCA even in the presence of significant dust absorption? In Sec. 4.1 and 4.2, we will consider optically thin emission lines, and in Sec. 4.3 we will discuss the effect of self-absorption. The full expression for intensity correlation in the presence of dust absorption is presented in Appendix A.

We start with discussing main qualitative effects. The optical depth due to dust absorption

[TABLE]

and the LOS distance $\Delta(\bm{X})$ at which optical depth becomes unity,

[TABLE]

play the central role in determining the relevant statistical properties of intensity correlation. In the absence of the turbulent velocity effects, or more specifically, when spectral line is integrated over $\Phi_{v}$ distribution, Eq. (8) becomes

[TABLE]

This is the case applicable for thick velocity slices, further detailed in Sec. 4.1. Main observation here is that when optical depth through the volume exceeds unity along all LOS, $\tau(\bm{X},S)>1$, the intensity fluctuations due to density fluctuations are exponentially suppressed. This comes from compensation between density variations along the line of sight and corresponding variations of the physical LOS distance $\Delta(\bm{X})$ at which $\tau$ reaches unity. This can be seen by approximating exponential cutoff in LOS integral by a step $\Theta$ function

[TABLE]

after which Eq. (8) gives constant intensity

[TABLE]

This suppression of density fluctuation is not present when densities of emitters and dust do not follow each other (in particular, when distribution of dust is taken to be uniform), or when intensity is proportional to the square of the emitters density. It is also absent when absorption is weak, namely optical depth through the cloud is less than unity, $\tau(\bm{X},S)<1$.

Let us see what happens when velocity effects are present, as, for instance, when intensity is measured in thin velocity slices. Applying the step approximation to the full Eq. (8) we have

[TABLE]

which highlights two main effects of the strong absorption. Strong absorption means small $\Delta(\bm{X})$ which suggests an estimate

[TABLE]

Thus, firstly, the signal is reduced to the nearest 2D surface of the cloud rather than being integrated along LOS through 3D volume, and, secondly, the density fluctuations are compensated away, leaving only the velocity induced intensity fluctuations. However, correlation studies introduce additional scales besides $S$, that can be much shorter, namely the density correlation length $r_{c}$ and the separation between pairs of LOS, $R$ which can be taken arbitrarily small (limited by the instrument resultion). Thus we must consider whether $\Delta/R$ and $\Delta/r_{c}$ are small or large before concluding whether $\Delta(\bm{X})$ is sufficiently small or not.

The estimate Eq. (27) seemingly reguires density to be near constant through $[0,\Delta(\bm{X})]$ LOS interval, but to what extend it is essential to the described picture? We can develop a more precise estimate rewriting Eq. (8) by integrating by parts

[TABLE]

to eliminate density pre-factor before applying $\Theta$-function approximation for the exponential term. This gives

[TABLE]

This expression is valid only for strong absorption, $\tau(\bm{X},S)>1$, so the exponential term can be neglected. We leave this exponentially suppressed term here just to show that Eq. (29) is consistent with Eq. (23) when the velocities in the spectral line are integrated over. Otherwise, we may just write

[TABLE]

This expression explicitly neglects the width of the region to $\Delta(\bm{X})$, but does not, however, assume anything about the density uniformity. Thus, the conclusion about 2D surface nature of the signal at high absorption holds independently on density properties. Density fluctuations are introduced only through fluctuations of the LOS positions of this surface, and one can expect that if $\Delta(\bm{X})$ does not fluctuate much, contribution of density fluctuations is small. This is the manifestation of density fluctuation suppression discussed above.

In the following sub-sections, we discuss with more rigour two important VCA regimes: thin velocity slice and thick velocity slice regimes.

4.1 Thick velocity slice

We are first interested in a simpler case of thick velocity channel. In this case, all emission along the LOS is collected independent of the velocity of emitters, which can be seen using Eqs. (8) and (16) with $W(v)=1$, which gives

[TABLE]

where $\mathcal{N}(\bm{X})$ is the column-density, and is given by

[TABLE]

Eq. (31) shows that dust absorption leads to exponential suppression of density contribution. Thus, making measurements in thick slices in the presence of dust is informative only when absorption is not too large (more accurate criteria for extent absorption is discussed further on). Due to the appearance of density dependence in exponential fashion, in the presence of dust one has to carry out further analysis with an explicit model of density field, unlike in LP00, where no assumption about the model of density field had to be made.

Statistics of the column density $\mathcal{N}(\mathbf{X})$ is discussed in the Appendix B. For density field with steep spectrum of fluctuations, the most relevant case is when density dispersion is less than the mean density333For steep density, perturbations of largest amplitude have wavelength comparable to the size of the cloud. If we determine mean density over the size of the cloud, then small scale perturbations naturally have amplitude less than the mean.. Therefore, we can model density fluctuations $\delta\rho$ as a Gaussian random field, with $\langle\delta\rho\rangle=0$. This translates to the Gaussian behaviour of the column density as well.

On the other hand, for shallow density field, a more appropriate and general model for the density field is log-normal field. Numerical calculations present in Appendix B show that column density $\mathcal{N}$ for a shallow density still tends to be Gaussian if the integration length $S$ along the LOS is much larger than the correlation length $r_{c}$ of density field (which can be easily satisfied for shallow spectra, as $r_{c}$ in this case is very small). This is because, outside of the correlation length, the fluctuations of density are essentially independent (as correlation decays outside of correlation length), and the sum of uncorrelated random fluctuations should tend to Gaussian quantity according to the Central limit theorem. The near Gaussian nature of the column-density is clearly illustrated in Fig. 14. Appendix B also defines and discusses the column density variance $\sigma_{\mathcal{N}}^{2}$ and correlation $\xi_{\mathcal{N}}(\mathbf{R})$ and structure $d_{\mathcal{N}}(\mathbf{R})$ functions and their ranges of utility.

With this, we are now ready to proceed further on. In configuration space, to obtain correlation function of fluctuations of intensity, one needs to have knowledge about the mean intensity. From Eq. (31), we have the mean intensity as

[TABLE]

and correlation of fluctuations of intensity $I(\mathbf{X})-\bar{I}$ as

[TABLE]

With the assumption of homogeneity, we have $\left\langle\mathrm{e}^{-\kappa\mathcal{N}(\bm{X}_{1})}\right\rangle=\left\langle\mathrm{e}^{-\kappa\mathcal{N}(\bm{X}_{2})}\right\rangle$, which allows us to simplify Eq. (4.1) to

[TABLE]

For a Gaussian $\mathcal{N}$, the first term in the Eq. (35) can be written as

[TABLE]

where $\xi_{\mathcal{N}}(\bm{R})=\left\langle\left(\mathcal{N}(\bm{X}_{1})-\bar{\mathcal{N}}\right)\left(\mathcal{N}(\bm{X}_{2})-\bar{\mathcal{N}}\right)\right\rangle$. Thus, using Eqs. (35) and (36) we finally have

[TABLE]

For $\kappa^{2}\xi_{\mathcal{N}}(R)<1$, the above Eq. (37) reduces to

[TABLE]

Note that $\kappa^{2}\xi_{\mathcal{N}}(0)=\langle\delta\tau^{2}\rangle$, the variance of the optical depth through the volume. Thus, the expression Eq. (38) holds everywhere, and one can recover density statistics, if absorption due to dust is sufficiently uniform, $\langle\delta\tau^{2}\rangle<1$. However, it also holds for sufficiently large $R$ even if dust absorption is less uniform.

As discussed in Appendix B, correlation function approach is useful if density spectrum is shallow with $\nu_{\rho}>1$. In this case Eq. (38) leads to a very interesting conclusion that one can recover density scaling as in usual thick slice VCA even in the presence of dust if one studies correlation of intensity at large lags $R$. The useful scale range is given by the condition $\kappa^{2}\xi_{\mathcal{N}}(R)<1$, which for $\nu_{\rho}>1$ sets a correlation length for the observed intensities at

[TABLE]

Thus, for $\nu_{\rho}>1$ one can obtain density scaling measuring the correlation function at lags $R>\mathrm{max}(r_{c}^{*},r_{c})$. This condition is more and more limiting as fluctuations of the optical depth grow beyond $\langle\delta\tau^{2}\rangle=1$, but may still be satisfied for shallow density which has small correlation length $r_{c}$.

For $\nu_{\rho}<1$ one should focus on measuring the structure function of intensity. Using Eq. (37), it can be written as

[TABLE]

which can be further simplified to

[TABLE]

At sufficiently small lags $R$, the usual result of the structure function of optically thin emission lines is recovered, i.e.

[TABLE]

while at larger lags $R$, the structure function is saturated. The saturation scale, which is, in effect, the correlation length of intensity, is estimated by setting $\frac{\kappa^{2}}{2}d_{\mathcal{N}}(r_{c}^{*})\sim 1$, which now gives

[TABLE]

where we note that Eq. (39) and Eq. (43) are equivalent. Thus, by measuring structure function of density in thick slices, one can recover density spectra at scales $R<\mathrm{min}(r_{c}^{*},r_{c})$. For steep spectra with $\nu_{\rho}<0$ and $r_{c}\sim S$ the critical scale is less than the size of the cloud, i.e. $r_{c}^{*}<S$ only when dust optical fluctuations are exceeding unity, while for weak absorption fluctuations it is comparable to the size of the cloud, i.e. $r_{c}^{*}\sim S$. Therefore, in the case when dust inhomogeneity is strong, the recovery of steep density statistics in thick velocity slices may be limited by the spatial resolution of telescope. The situation is the least advantageous for $0<\nu_{\rho}<1$ range of shallow spectra where we should use the structure function measurements but need to resolve $r_{c}$ scale which now can be small.

We stress that the amplitude of correlations in both Eq. (38) and Eq. (41) will be exponentially suppressed if optical depth of the cloud is much larger than unity, and is challenging to measure if dust absorption is strong and instrumental noise is not properly factored out.

4.2 Thin velocity slice

To facilitate comparison between the case of presence and absence of dust absorption, we first review the main results of VCA in thin slices in the absence of dust.

In the absence of dust, Eq. (8) is reduced to

[TABLE]

i.e. the spectral intensity of optically thin emission lines is directly given by PPV density of emitters. In this case, the correlation of spectral intensity is given by

[TABLE]

where we have assumed uncorrelated density and velocity field. Assuming velocity fluctuations to be Gaussian random quantities, Eq. (4.2) can be used to finally obtain correlation of thin slice (see LP04)

[TABLE]

thus the structure function is

[TABLE]

For steep density spectra, the small $R$ asymptote of $\mathcal{D}(R)$ is $R^{1-\nu/2}$, while for shallow spectra, it is $R^{1-\nu_{\rho}-\nu/2}$. The latter asymptotics actually describes the scaling of $\mathcal{D}(R)$ when density index $\nu_{\rho}<1-\nu/2$. This presents a unique way to obtain velocity spectra of the turbulence.

We now study thin velocity slice limit in the presence of dust absorption. For that, we start with Eq. (30). We remind the reader that although the representation in Eq. (30) is able to describe regimes where dust absorption is strong, it cannot describe the regimes where dust absorption is weak. The weak absorption regime will be studied numerically later in this section. Strong absorption requirement for validity of Eq. (30) means short LOS distance $\Delta(\mathbf{X})$ to $\tau=1$ relative to all scales of the problem, namely $\Delta(\mathbf{X})<R\leq S$.

We now carry out explicit analytical calculations for the correlation of intensity in the presence of strong dust absorption. Using Eq. (30), the correlation of intensity in the thin velocity slice and in the strong dust absorption regime can be written as

[TABLE]

where $\langle..\rangle$ denotes statistical averaging in both velocity field and density field. Under the assumption of uncorrelated density and velocity field, this averaging can be performed independently, and thus we first carry out averaging with respect to velocity field by assuming Gaussian velocity field. This allows us to write Eq. (48) as

[TABLE]

where the remaining averaging is with respect to $\Delta(\bm{X})$ as a random fluctuating quantity. Using Eq. (13), correlation given by Eq. (49) can be explicitly written as

[TABLE]

where we have restored all the relevant dimensionalities. A crucial point that needs to be explained at this point is the regime of validity of Eq. (50). As a part of explanation, we take the case of constant density. In this case, on one hand, our approximate expression Eq. (50) gives a scaling of intensity correlation as $\sim R^{-\nu/2}$. On the other hand, an exact expression of intensity correlation can be obtained by using Eq. (46), with $\xi_{\rho}$ set to zero, and $S$ set to $\Delta_{0}$. This expression gives a correlation scaling of $R^{-\nu/2}$ only at $R>\Delta_{0}$. We stress that density effects are suppressed by dust absorption, and therefore it is reasonable to deduce the regime of validity of Eq. (50) by considering the constant density case. With all these discussion, we suggest the statistically averaged criterion for the validity of Eq. (50)

[TABLE]

Appendix C discusses statistics of $\Delta$ in various regimes of density perturbations, in particular gives expressions for $\langle\Delta^{2}\rangle$ in Eqs. (143, 159). For small perturbations, $\Delta_{0}\approx\langle\Delta\rangle\approx 1/(\kappa\bar{\rho})$. At lags $R$ satisfying Eq. (51), we obtain from Eq. (50) the density insensitive scaling of the intensity correlation function

[TABLE]

To evaluate averaging in Eq. (50) accurately, one needs to properly model the statistics of $\Delta$, in particular the probability distribution function of $\Delta$. To give specific example, let us take dust density to be Gaussian. This is reasonable when the density dispersion is less than the mean density. In the case when density fluctuations are small in comparison to the mean density (which is good for steep density spectra, but more restrictive for shallow density spectra), the above equation reduces to

[TABLE]

Performing averaging with Gaussian $\delta\rho$, one obtains

[TABLE]

where $U$ is the confluent Hypergeometric function, and $\bar{d}_{\rho}(R)$ is the normalized structure function of density fluctuations

[TABLE]

The asymptote of Eq. (54) is $R^{-\nu/2}$ at $R^{2}>\bar{d}(R)$. Since $\bar{d}(R)\leq 2\sigma_{\rho}^{2}/(\kappa^{2}\bar{\rho}^{4})\approx\langle\Delta^{2}\rangle\sigma_{\rho}^{2}/\bar{\rho}^{2}$, this asymptote matches the general result of Eq. (52) over all its range of validity.

At scales $R<\Delta_{0}$ the thickness of the transparent layer at the front of the cloud cannot be considered negligible, but more distant regions do not contribute. In other words, in place of Eq. (46) we can write

[TABLE]

Thus, dust absorption does not modify VCA predictions in this regime. This has been represented in Fig 2, where for $R\gtrsim\Delta_{0}$, the asymptote of correlation function is clearly $R^{-\nu/2}$, whereas for $R<\Delta_{0}$, the asymptote of structure function follows the usual VCA scaling. Eq. (56) also obviously covers the regime of weak absorption when $\Delta_{0}\gtrsim S$, in which case $\Delta_{0}$ should be replaced by $S$.

A crucial point that needs to stressed at this point is that Eq. (56) is applicable only when dust absorption is weak or moderate, i.e. $\Delta_{0}\sim S$. If dust absorption is strong, i.e. $\Delta_{0}\ll S$, density effects are erased, and the intensity correlation is given by

[TABLE]

To test this, we carry out numerical integration with exact analytical expression presented in Appendix A. Expression presented in the Appendix A has its own limitation that it is not applicable if the density dispersion is larger than the mean density, as this becomes inconsistent with the assumption of Gaussianity. This poses challenge in deducing whether or not the density effects are erased due to strong dust absorption, as for small $\sigma$, the pure velocity part dominates while the density effects are sub-dominant regardless of the strength of dust absorption. To overcome this difficulty, we carry out numerical integration of just the density contribution part in the correlation given by Eq. (A). As one might already see in Eq. (A), the $\xi_{\rho}$ part is partly compensated by the differential part that comes with the negative sign. As shown in Fig. 1, the contribution of density part to intensity structure function flattens as dust absorption becomes strong, suggesting that density effects are unimportant in this regime.

The validity of Eqs. (56) and (57) is further tested by carrying out numerical integration of the full expression of intensity correlation presented in Appendix A. The results are presented in Fig. 2 for steep density spectra. Left panel in Fig. 2 shows the correlation function of the intensity, confirming that large lag $R$ asymptotics follows Eq. (52) while right panel shows the structure function, emphasizing small $R$ scaling. Large $R$ correlation scaling is confirmed to be insensitive to the density fluctuations.

The dust cut-off depth $\Delta_{0}$ in the Eq. (51) is the scale $R$ where the structure function starts to saturate, so it plays the role of the correlation length of the intensity fluctuations. It is close to $1/(\kappa\bar{\rho})$, though an accurate analytical determination of this scale is not possible within our asymptotic analysis. However, our study shows that measuring the intensity correlation length in observations of dusty regions provides direct information about the level of dust absorption.

At small $R$ such that $D_{z}(R)<\Delta V^{2}$ the velocity slice thickness $\Delta V$ can not be neglected and thin slice approximation in Eq. (56) is not applicable anymore. To cover this regime, we should restore the velocity integral in the correlation function, which for strong dust absorption leads to

[TABLE]

An apparent complication that might arise with our analysis is that while at small $R$ structure function is applicable, correlation function is used at larger $R$. An unified picture of the regimes with strong dust absorption and that with weak dust absorption can be obtained with the language of power spectrum. The power spectrum $P(K)$ of structure function of intensity is defined as

[TABLE]

which signifies that the power spectrum of structure function and correlation function are the same upto some constant factor.

In the language of power spectrum, one can see the break in the slope of power spectrum as one goes from high $K$ to low $K$ in the thin slice limit. This has been clearly shown for both steep and shallow density in Fig. 3. Notice in Fig. 3 that the slope changes from $-2+\nu/2$ to $-3+\nu/2$ for steep density and to $-3+\nu_{\rho}+\nu/2$ for shallow density with weak dust absorption. Thus, it is useful to use correlation function at low $K$ and structure function at high $K$.We stress that the main effect is switching from integrated signal at small $R$ to just two dimensional slice in configuration space at large $R$. However, the non-trivial effect is at large $R$, or for strong dust absorption density effects are erased.

4.3 Dust effect on self absorbing medium

The study in the previous sections were carried out for optically thin emission lines, i.e. for the case when self-absorption is negligible. However, effect of self-absorption can be important in various interstellar environments, for instance in molecular clouds. An important study on the effects of self-absorption in the intensity statistics was carried out in LP04. Their study suggests that self-absorption introduces non linear effects on the statistics of intensity, and that the recovery of velocity and density spectrum might not be always possible.

With the result of previous section, one can extend the results of previous section to account for the effect of both self-absorption and dust-absorption. The main result of LP04 for structure function of intensity fluctuations can be expressed using $\xi_{s}(\bm{R},v)$ as:

[TABLE]

The only difference between this result and the result of LP04 is the replacement of $S$ by $\Delta_{0}$ in all the integrals that are carried out over the LOS. To quantify the effect of dust on self-absorption, we first start with the expression of $\tilde{d}_{s}(0,v)$ in Eq. (60) (see LP04)

[TABLE]

The difference between this expression and the one in LP04 is the replacement of $S$ by $\Delta_{0}$. Eqs. (60) and (61) physically imply that dust absorption decorrelates signals in physical space, while self-absorption decorrelates signals in velocity space. From Eqs. (60) and (61), one can immediately see that the effective window introduced by self-absorption decreases for any $\Delta_{0}<S$, thus diluting the effect of self-absorption. This is an interesting result: stronger dust absorption results in weaker self-absorption.

The results of LP04 suggest that in the case of negligible dust absorption, self-absorption leads to a universal regime where structure function or power spectrum of intensity fluctuations are independent of underlying velocity spectra, i.e. $\mathcal{D}(R)\sim R$. This extent of this regime is set by the strength of self-absorption, and is seen at scales

[TABLE]

where $v_{ab}$ is the velocity scale at which $\alpha^{2}\tilde{d}_{s}(0,v_{\text{ab}})=1$ (c.f Eq. (61)). As a simple approximation, in the presence of dust absorption, the scale at which this universal regime is seen is given by

[TABLE]

Clearly, the scales at which this regime is seen in the presence of dust absorption is reduced by a factor $\Delta_{0}/S<1$. Physically, this can be understood in the following sense: strong dust absorption decorrelates signals from separations larger than the dust cut-off, and therefore the extent of self absorption a signal experiences is decreased.

5 Effect of dust on VCS

Velocity coordinate spectrum (VCS, LP06) is another powerful technique to study supersonic turbulence. This technique measures power spectrum of intensity fluctuations along the velocity axis in PPV space. With VCS, one can study turbulence even in a spatially unresolved cloud. In this section, we discuss how dust absorption might affect studies of turbulence using VCS.

As it was shown in LP00, the power spectrum of intensity fluctuations along the velocity axis is given by

[TABLE]

which is

[TABLE]

where $B(\bm{R})$ is the window of the telescope beam, which describes the resolution of the instrument. The two regime of interest are high spatial resolution (narrow beam) and poor spatial resolution. For a narrow beam, $B(\bm{R})=\delta(R)$, which gives asymptotic form of $P(k_{v})$ at high $k_{v}$ for optically thin emission lines with negligible dust emission (see LP06)

[TABLE]

while for a wide beam the window is $B(\bm{R})=1$ (which upon integration over the entire image is equivalent to setting $\bm{K}=0$), and the asymptote is

[TABLE]

As was shown in previous sections, the main effect of dust absorption is the introduction of cut-off scales beyond which intensity signals are decorrelated. Thus, the dust cut-off scale $\Delta$ effectively acts as LOS extent $S$ of the turbulent cloud, and Eq. (5) can be written as

[TABLE]

where $\Delta_{\text{av}}$ is some average of dust cut-off, and comes from the fact that the dust absorption cut-off along different LOS are different. In general, $\Delta_{\text{av}}\sim\Delta_{0}$. In the case when $k_{v}^{2}D_{z}(\Delta_{0})\gg 1$, the LOS integration limit in Eq. (5) can be extended to infinity, fluctuations of $\Delta$ along different LOS does not affect the averaging, allowing one to recover the usual results of VCS. Formally, this condition corresponds to

[TABLE]

If Eq. (69) is satisfied, the asymptotic forms represented by Eqs. (66) and (67) for narrow and wide beam respectively will be valid even with the presence of dust absorption, as shown in Fig 4, where the asymptote is established at some $k_{v}$ regardless of what the dust cut-off is. However, it is important to keep in mind that these asymptotic forms are applicable only at large wavenumbers, which depends on the value of $\Delta_{\text{av}}$ is. This is what Fig. 4 demonstrates, where it is shown how different dust cut-off $\Delta_{\text{av}}$ lead to asymptote at different $k_{v}$.

An important scale is the velocity scale $V_{\Delta B}$ where the transition from poor to high resolution occurs

[TABLE]

It is clear from Eq. (70) that due to dust absorption the width of the telescope beam effectively increases from $\Delta B$ to $\Delta B(S/\Delta_{\text{av}})$, which implies that dust absorption makes spatial resolution poorer. At scales

[TABLE]

the beam is narrow, while on shorter scales its width is important. In the case of negligible thermal broadening, the maximum value of $V_{\Delta_{B}}$ is $\sqrt{D_{z}(S)}$. Thus, a narrow beam requires $\Delta B<\Delta_{0}$, implying that it might not always be possible to achieve narrow beam especially if dust absorption is strong (i.e. $\Delta_{0}$ is small).

The condition required for the asymptotes represented by narrow beam is more stringent than what Eq. (71) gives. Combining Eqs. (69) and (71), one can see that there is a range of wavenumbers $k_{v}$ for which one can achieve narrow beam:

[TABLE]

An important result that can be inferred from the above equation is that while the above equation can be easily satisfied by large $k_{v}$ for large $\Delta_{\text{av}}$, the range of $k_{v}$ becomes narrower as $\Delta_{\text{av}}$ decreases. In particular, if $\Delta_{\text{av}}\lesssim\sqrt{S\Delta B}$, one cannot achieve the narrow beam.

It is also important to note that if dust absorption is strong, the power spectrum of fluctuations of intensity in velocity space is purely due to velocity effects, and the contribution due to density fluctuations are suppressed. In the weak dust absorption regime, the usual VCS statistics is restored.

6 Centroids

Centroid, which is velocity moment of intensity, is another powerful technique to study velocity field in a turbulent medium using Doppler broadened emission lines. Traditionally, centroid is defined as

[TABLE]

The normalisation in the definition of centroid makes analytical construction difficult. To overcome this, LE03 introduced a notion of unnormalized velocity centroids (UVC), which is defined as

[TABLE]

For optically thin emission lines with negligible dust-absorption, the above definition gives

[TABLE]

and the structure function in the case of constant density is given by

[TABLE]

which at $R\ll S$ has an asymptotic form of

[TABLE]

Thus, by measuring the structure function of centroids, one can obtain the spectral slope of velocity field.

We now focus on effects of dust absorption for the studies of turbulence using velocity centroids. The formalism developed in Sec. 4.2 can be easily extended to centroids. We start with the definition of centroids

[TABLE]

In the case of weak dust absorption, we do not expect much changes in the centroid statistics. However, for strong dust absorption, using similar idea as in Sec. 4.2, one can easily obtain

[TABLE]

Therefore, the structure function of centroids is given by

[TABLE]

where the averaging is over $\Delta$ as a random fluctuation. The above approximation is valid for $R>\Delta_{0}$. For $R<\Delta_{0}$, one should obtain the usual centroid structure function. With this, one can write

[TABLE]

The above scalings are also represented in Fig. 5.

Notice that unlike in VCA, dust does not change the statistics from steep at short scales to shallow at large scales.

7 Anisotropic MHD Turbulence

Our analysis in the previous sections focused on the spectra of density and velocity field, and thus we worked with the assumptions that the fields are isotropic and homogeneous. The interstellar medium is magnetised and thus there exists a preferred direction (along the direction of magnetic field), and therefore isotropy is not a good assumption any more. The properties of MHD turbulence depend on the degree of magnetization. An useful measure of the degree of magnetization is the Alfvén Mach number $M_{\text{A}}=V_{L}/V_{A}$, where $V_{L}$ is the injection velocity at the scale $L$ and $V_{A}$ is the Alfvén velocity. For $M_{\text{A}}\gg 1$ magnetic forces are not important at large scales and the cascade should be similar to ordinary hydrodynamic cascade in the vicinity of the injection scale. The modern theory of MHD turbulence was proposed in Goldreich & Sridhar (1995, hereafter GS95), where it was argued that turbulence eddies are strongly anisotropic, and are elongated towards the direction of local magnetic field. This theory was elaborated in further studies, in particular the concept of local system of reference (Lazarian & Vishniac 1999; Cho & Vishniac 2000; Maron & Goldreich 2001; Cho et al. 2002) was introduced. According to this concept, turbulent motions should not be viewed in the system of reference of the mean magnetic field as all earlier theories of MHD turbulence attempted to do, but in the system of reference of magnetic field comparable with the size of the eddies. However, from the point of view of the observational studies of the turbulence in a volume when the only available statistics are those averaged along the line of sight, the measurements should be carried out in the system of mean magnetic field, rather than the local system of reference. Therefore, one has to describe Alfvénic turbulence in the global system of reference (see the discussions in Cho & Lazarian 2002; Esquivel & Lazarian 2005; LP12). This modifies the available statistics. For instance, in the local system of reference GS95 predicts the existence of two different energy spectra, namely, the parallel and perpendicular, in the global system of reference only the spectrum of dominant perpendicular fluctuations is available. Similarly, while in the local system of reference the anisotropy increases with the decrease of size of the eddies, the anisotropy stays constant in the global system of reference.

MHD turbulence can be presented as a superposition of interacting fundamental modes, i.e. Alfvén, slow and fast. The first theoretical considerations in favour of this were given in GS95 (see also Lithwick & Goldreich 2001), which were extended and numerically tested in Cho & Lazarian (2002, 2003) and in Kowal & Lazarian (2010). Because the compressible and incompressible modes weakly exchange their energy (Cho & Lazarian 2002), it is possible to consider the modes separately. The details of velocity structure function of these fundamental modes are presented in Kandel et al. (2016b, hereafter KLP16).

In this section, we aim to discuss effects of dust absorption on the study of anisotropy. For that we will explicitly use some of the important results from KLP16. We first briefly review different fundamental MHD modes.

7.1 Different MHD modes

In a magnetised medium, turbulence can be described by decomposition of basic MHD modes : Alfvén, slow and fast. This decomposition is carried out with respect to the mean magnetic field, and is considered to be reasonable as long as perturbations of magnetic field are less than the mean magnetic field. In general, the Fourier component of velocity in a mode is given by $v(\bm{k})=a_{k}\hat{\xi}_{k}$, where $\bm{k}$ is the wavevector, $a_{k}$ is the random complex amplitude of a mode and $\hat{\xi}$ is the direction of allowed displacement. With this, the power spectrum can be written as

[TABLE]

where $\mathcal{A}(k,\hat{k}\cdot\hat{\lambda})\equiv\langle\hat{a}_{k}\hat{a}_{k}^{*}\rangle$ is the power spectrum of amplitude fluctuations, and is in general anisotropic (as its angular dependence reflects), and $(\hat{\xi}_{k}\otimes\hat{\xi}_{k}^{*})_{ij}$ is the tensor structure of the mode, which is also in general anisotropic.

In the subsequent subsections, we briefly review general properties and power spectrum of each of the fundamental MHD modes.

7.1.1 Alfvén mode

Alfvén modes are incompressible modes with most of the energy residing in them at around the wavenumbers where eddy turnover time is equal to the period of the Alfvén wave. The displacement of this mode in a plasma is orthogonal to the wave-vector and the direction of magnetic field, i.e. $\hat{\xi}_{k}\propto\hat{k}\times\hat{\lambda}$, which leads to anisotropic tensor structure of Alfvén mode as (see KLP16)

[TABLE]

In the global system of reference (which is the only reference system available for observations), the anisotropy of power spectrum is determined by the anisotropy at outer scale, and the power spectrum of amplitude fluctuations is given by

[TABLE]

where $\hat{\lambda}$ is the unit vector pointing towards the direction of mean magnetic field, and $M_{A}$ is the Alfvén Mach number. Clearly, if the LOS direction is aligned with the direction of magnetic field, then one can not see any anisotropic effect.

7.1.2 Fast mode

Fast modes are compressible modes and in a plasma with $\beta\equiv P_{\text{gas}}/P_{\text{mag}}\gg 1$, they behave like sonic waves, while in a low-$\beta$ plasma, they propagate with Alfvén speed due to the compressions of magnetic field. In the absence of collisionless damping, turbulent cascade of fast mode are expected to persist over time spans longer than that of Alfvén or slow modes, and they are expected to marginally interact with Alfvén modes (Cho & Lazarian, 2002). The power spectrum of amplitude fluctuations of this mode is isotropic, i.e.

[TABLE]

The tensor structure of this mode can be derived by considering the allowed displacement of this mode in a plasma. In the high-$\beta$ regime, the displacement of fast mode is purely radial, i.e. parallel to wavevector $\hat{k}$, which leads to isotropic tensor structure

[TABLE]

while in the low-$\beta$ regime, the displacement of this mode is orthogonal to the direction of magnetic field, i.e. $\hat{\xi}_{k}\propto\hat{\lambda}\times(\hat{\lambda}\times\hat{\bm{k}})$, and leads to tensor structure

[TABLE]

which is clearly anisotropic.

7.1.3 Slow mode

Slow modes in high-$\beta$ plasma are similar to pseudo-Alfvén modes, while in low-$\beta$ plasma they are density perturbations propagating with sonic speed parallel to magnetic field (see Cho & Lazarian 2003). The power spectrum of this mode is anisotropic, and is in fact same as that of Alfvén mode. The tensor structure of this mode is derived by considering the allowed displacement of this mode in a plasma. In the high-$\beta$ regime, the displacement of this mode is perpendicular to the wavevector $\hat{k}$, i.e. $\hat{\xi}_{k}\propto\hat{\bm{k}}\times(\hat{\bm{k}}\times\hat{\lambda})$, leading to anisotropic tensor structure (see KLP16)

[TABLE]

On the other hand, in the low-$\beta$ regime, the displacement of this mode is parallel to the direction of magnetic field, i.e. $\hat{\xi}_{k}\propto\hat{\lambda}$, and leads to tensor structure

[TABLE]

which is also clearly anisotropic.

7.2 Effects of dust absorption on VCA and centroid anisotropy

The intensity statistics of emission lines is dependent on the LOS velocity. In a magnetised medium, motions, and hence intensity statistics, are affected by the magnetic field. These effects are manifested through the anisotropic nature of correlation and structure function of intensity and centroids. The study of intensity anisotropy was carried out in KLP16 and the study of centroids anisotropy was carried out in KLP17. In this section, we will discuss how the results in these paper would be affected by dust-absorption.

Due to the cut-off introduced by dust absorption, it was shown in Section 4 that the statistics of correlation function and structure function is modified at different lags $R$. In particular, for thin slice case, it was shown that intensity statistics changes from steep at small lags $R$ to shallow at large lags $R$. Similarly, it was shown that in the thick slice regime, intensity structure function saturates for steep density at $R$ larger some critical value, which depends on dust absorption coefficient. Thus, it is natural to expect different anisotropic behaviour of intensity correlation and structure function depending upon the lags $R$ under consideration.

Let us first start with thick velocity slice regime of optically thin emission lines. From Eq. (41), it is clear that for steep spectra at large scales $R$, the intensity structure function saturates, while at small scales the usual VCA structure function is recovered. Thus, one can conclude that at small scales, the anisotropy will be the same as that predicted by VCA (see KLP16), while at large scales the statistics is isotropic. Thus, without a good spatial resolution, one cannot study anisotropies in the presence of dust absorption. On the other hand, shallow density is usually isotropic, and thus we do not expect anisotropies in correlation function of intensity in thick slice regime if density is shallow.

We now consider thin slice case. In this case, at short lags $R$, dust absorption does not affect power law scaling, and level of anisotropy is similar to that predicted in KLP16. On the other hand, results from Sec. 4.2 show that the structure function is saturated at scales $R>\Delta_{0}$, and thus correlation function is better measure at these scales. Therefore, one can expect isotropization of intensity statistics at these scales if one uses structure function. On the other hand, anisotropies are manifested in the intensity correlation function at these scales, and the anisotropy of correlation function of intensity is directly related to the anisotropy of velocity structure function. It is clear from this discussion that the study of anisotropy in the thin slice regime is complicated by the fact that one has to use two different statistics at two different scale ranges. A way to unify the study of anisotropy is to study anisotropy of power spectrum of intensity fluctuations. An anisotropic power spectrum can be decomposed into multipoles

[TABLE]

where $P_{n}$ is the multipole moments of power spectrum. An useful measure of level of anisotropy in the language of power spectrum is the quadrupole to monopole ratio $\mathcal{R}$ defined as

[TABLE]

which depends on the scale $K$ and magnetisation of media.

The anisotropy level of different MHD modes at small lags $R$ has already been presented in KLP16. Here, we proceed to compute this anisotropy level for lags $R>\Delta_{0}$. From Eq. (50), and for sufficiently large lags, it is enough for our discussion to write the correlation function as

[TABLE]

First important conclusion is that at large lags $R$, the level of anisotropy is independent of $R$. However, level of anisotropy is highly dependent on Alfvén Mach number $M_{A}$ and $\gamma$. Here, we use their result to predict the level of anisotropy of intensity fluctuations due to various MHD modes. The quadrupole to monopole ratio of power spectrum of intensity fluctuation in the thin slice regime of VCA has been illustrated for various MHD modes in figures 7, 8, 9 and 10. These figures clearly show that the quadrupole to monopole ratio of each of these MHD modes change their sign from small wavenumber $K$ to large $K$. As Alfvén modes and low-$\beta$ fast modes do not contribute at $\gamma\sim 0$, while slow modes do not contribute at $\gamma\sim\pi/2$, in a realistic situation involving mixture of different modes, the anisotropy level decreases at large $K$ in comparison to that at small $K$. It is important to stress that the change in anisotropy level from small $K$ to large $K$ is simultaneously accompanied by the change is the spectral slope of the power spectrum. In the presence of self-absorption, it was shown in KLP16 that the anisotropy level in the universal regime is similar to that of the thin slice regime. As our finding in Section 4.3 suggests, dust absorption restricts the range of universal regime, but the anisotropy level in this universal regime is still similar to that of the thin slice regime of negligible dust absorption.

The discussion carried above applies to anisotropy of centroids as well. At scales $R<\Delta_{0}$, the anisotropy level of centroids in the presence of dust-absorption is the same as that with the presence of dust-absorption. At $R>\Delta_{0}$, the centroid structure function is to a good approximation

[TABLE]

and combining this expression with the results of KLP16, one can easily obtain the anisotropy level of different MHD modes. The anisotropy level of centroid structure function of different MHD modes at $R\gtrsim\Delta_{0}$, where dust absorption is significant, is shown in figures 11 and 12. If one compares these anisotropy levels with that in the absence of dust absorption (see KLP17), it is easy to see that the anisotropy level decreases at these scales in comparison to the one in the absence of dust or at lags $R\lesssim\Delta_{0}$. It is again important to note that the change of anisotropy level in going from $R\lesssim\Delta_{0}$ to $R\gtrsim\Delta_{0}$ is accompanied by the change in the spectral index of centroid structure function.

8 Effect of dust on collisionally excited emission lines

In some interstellar environments, like H II region, the ultra-violet, visible and infra-red spectra of emission lines are primarily collisionally excited lines of metal ions and recombination lines of Hydrogen and Helium. In this section, we study how collisionally excited lines can be used to obtain spectra of density and velocity field in a turbulent medium.

8.1 VCA for collisionally excited emission lines

The spectral intensity of a collisionally excited emission line in the presence of negligible dust absorption is given by

[TABLE]

where $\rho_{g}$ is the density of emitters and $\rho_{g}^{2}$ characterises the collision rate. Under the assumption of zero correlation between density and velocity, all the major equations given by VCA is still applicable for collisionally excited emission lines, but the density correlation has to be properly modelled. The main equation for the correlation of intensity is given by

[TABLE]

Eq. (8.1) differs from the usual expression only through the structure of density correlation $\langle\rho_{g}^{2}(\bm{x}_{1})\rho_{g}^{2}(\bm{x}_{2})\rangle$. For sub-sonic turbulence, which has steep density spectra, density fluctuations are smaller in comparison to the mean density and therefore breaking density into mean and fluctuating part $\rho_{g}=\rho_{0}+\delta\rho_{g}$ gives

[TABLE]

Thus, at this level of approximation, one returns to the usual predictions of VCA (with emissivity proportional to density of emitters rather than the square of their density).

The subsequent part of section is focused specifically for supersonic turbulence, which has shallow density spectra, and thus has density fluctuations which are comparable to the mean density. For such density spectra, linear assumption is not reasonable. One way to proceed is assume a particular form of density PDF. Various studies have shown that the density PDF of molecular clouds (where turbulence is supersonic) is lognormal. Therefore, we model dust density $\rho$ as

[TABLE]

where $\mu$ is a random Gaussian field with zero mean and variance $\sigma_{\mu}^{2}$ and is described by the two-point probability distribution function given by Eq. (133). Using Eqs. (139) and (139), it can be shown that

[TABLE]

We illustrate the above correlation with a specific model. For a shallow density, the correlation

[TABLE]

for $r>r_{c}(\sigma/\bar{\rho})^{-2/\nu_{\rho}}$. Eq. (100) represents a form of correlation that is similar to that of density correlation represented in Eq. (15), with an overall factor and different weight on amplitudes of correlation. Thus, for scales $R>r_{c}(\sigma/\bar{\rho})^{-2/\nu_{\rho}}$, all analytical expression derived in VCA are applicable to collisionally excited lines.

In the presence of dust absorption, the intensity of collisionally excited line is given by

[TABLE]

If density fluctuations are smaller than the mean, i.e. $\delta\rho<\bar{\rho}$, the factor $\rho^{2}$ can be approximated as $\rho^{2}\approx\bar{\rho}^{2}+2\bar{\rho}\delta\rho$. At this level of approximation, all the results that were derived in Sec. 4 are applicable to collisionally excited emission lines as well. Therefore, for sub-sonic turbulence (for which density fluctuations are smaller than the mean), the results of Sec. 4 are applicable.

For supersonic turbulence, density fluctuations are usually comparable to the mean density. To obtain an appropriate model of intensity correlation can be worked out for strong dust absorption regime, we carry out integration by parts of the spectral intensity. We do this for thick and thin velocity slice regimes separately.

In the thick velocity slice regime, by integrating over entire range of velocity, one obtains

[TABLE]

At this point, it is instructive to rewrite Eq. (102) as

[TABLE]

which upon integration by parts gives

[TABLE]

where $\Delta$ is the physical depth where optical depth reaches unity. The statistics of $\Delta$, including its mean $\bar{\Delta}$ and variance $\sigma_{\Delta}$, is presented in Appendix C. Comparing Eqs. (8.1) and (25), one can see that even in the case of strong dust absorption, density effects do not get washed away for collisionally excited lines.

Using Eq. (8.1), the correlation function of intensity correlation can be written as

[TABLE]

The evaluation of the above correlation is tricky due to the fact that $\Delta$ is itself a function of density. For fluctuating density, $\Delta$ across different lines of sight are different. The overall effect of fluctuating density is therefore smearing of the correlation. This effect can be captured by evaluating the average of density correlation function over a range of $\Delta$ given by $[\bar{\Delta}-\sigma_{\Delta}/2,\bar{\Delta}-\sigma_{\Delta}/2]$ such that

[TABLE]

Clearly, if $\sigma_{\Delta}$ is large, the correlation gets smeared out, while for small $\sigma_{\Delta}$, the intensity correlation function directly gives the density correlation function of approximately a two-dimensional thin slab. In particular, one can obtain useful information about the correlation function of density by studying correlation function of intensity at lags $R$ given by

[TABLE]

while at shorter lags, the correlation is saturated.

Using similar idea as presented above, one can readily show that under the assumption of uncorrelated density and velocity field, the intensity correlation function in the thin-slice regime in the presence of strong dust absorption is

[TABLE]

where $D_{z}(R,z)$ is the velocity structure function. Unlike our previous discussion of the case which excluded collisionally excited lines, we see that in this case the density effects are not erased even at large lags. Again, for small $\sigma_{\Delta}$, one obtains correlation of intensity of approximately two-dimensional thin physical slab.

8.2 VCS for collisionally excited emission lines

For collisionally excited emission lines, the power spectrum of intensity fluctuations along velocity coordinate in the presence of dust-absorption is given by

[TABLE]

where $\xi_{\rho^{2}}(\bm{r})\equiv\langle\rho^{2}(\bm{x}_{1})\rho^{2}(\bm{x}_{2})\rangle$. The difference between this case and the usual VCS comes through the replacement of $\xi_{\rho}(\bm{r})$ by $\xi_{\rho^{2}}(\bm{r})$. However, as it was shown in the previous section, the form of $\xi_{\rho^{2}}(\bm{r})$ looks structurally similar to that of $\xi_{\rho}(\bm{r})$ for both steep and shallow spectra. Therefore, all the results of VCS from Sec. 5 apply to collisionally excited emission lines. In particular, the extent of dust cut-off sets the range of wavenumbers where one can achieve narrow beam, and this range is given by Eq. (72). For small dust cut-off scales, achieving narrow beam is difficult. All the relevant asymptotics are the same as in the usual VCS, and are given by Eqs. 66 and 67.

8.3 Centroids for collisionally excited emission lines

The centroids were shown to be useful for studies of subsonic turbulence. As it was mentioned earlier, density fluctuations are smaller than mean density for subsonic turbulence. Therefore, all the results derived in Sec. 6 is applicable even for collisionally excited emission lines.

9 Comparison with observations

Several studies (see Münch & Persson (1971)) suggest that the Orion Nebula is dusty. The study of the reddening of the brightest region of the Orion Nebula as well as the behaviour of the He I lines $2^{3}P\rightarrow 2^{3}S$ and $3^{3}S\rightarrow 2^{3}P$ in Münch & Persson (1971) suggests that the dust and gas in the Orion are well mixed. Thus, the results of this paper are applicable to this nebula. More importantly, ionised gas usually leads to collisionally excited emission lines. In this paper, we have discussed how VCA works with collisionally excited emission lines and shown that the previous analytical results are applicable even for the collisionally excited emission lines. The study of turbulence in ionised Orion Nebula has been carried out by Arthur et al. (2016, hereafter AMH16) by analysing observational data through VCA and centroids. Thus, our results validate the calculations performed by AMH16 in order to find spectral index of velocity and density fields. However, one important point to note is that for strong dust absorption, density effects are suppressed. Therefore, this changes the result for shallow density spectra.

Of the results presented in AMH16, the most relevant ones for comparison with our prediction in the power spectrum of intensity fluctuations of various spectral lines in the thin velocity slice regime. In this regime, our prediction was that the power spectrum of intensity fluctuations exhibits a break in the slope, and the breakpoint is characteristic to the mean strength of dust absorption. Moreover, the change in slope between the two regimes should be of unity. We used data points obtained from AMH16 to perform fit of the largest scales (smallest $K$), which was missing in AMH16. The so obtained plots for various spectral lines are shown in Fig. 13. Note that we do not include any dispersion in the slopes of power spectrum, unlike in AMH16, where dispersion in slopes were from the different samplings of the PPV slits. Also notice the difference in Fig. 13, where $K^{3}P(K)$ is plotted against $K$ and our schematic representation in Fig. 3, where $P(K)$ is plotted against $K$.

As shown in Fig. 13, clearly, there is a break in the spectral slope, and the difference in the slope is close to unity for all the spectral lines. Moreover, the characteristic break point for all the spectral lines is of order $\sim 20$ pc*-1*, which corresponds to the dust-absorption cut-off $\Delta_{0}\sim 0.05$ pc. For a medium ionization zone of Orion nebula, whose markers are [O III], H II, He I and [Cl III], O’dell 2001 quote that the optical depths (which should be of similar to dust-absorption cut-off) is $\sim 0.06$ pc, totally consistent to our $\Delta_{0}\sim 0.05$ pc. Thus, one possible reason for this break is the dust extinction. Notice that, in Fig. 13, different spectral lines have different break points, which is naturally expected as different spectral lines have different $\kappa$ values. Also notice that some of the markers in Fig. 13, [N II] marks low ionization zone and [S II] marks ionization front. O’dell (1994, 2001); O’dell & Wen (1994) find that the ionization front has an optical depth of $10^{-4}$ pc and low ionization zone $2\times 10^{-3}$ pc. Our $\Delta$ for these spectral lines seems to be different from the results of O’dell (1994, 2001); and O’dell & Wen (1994), thus it is difficult to assert if the break in the slopes for these spectral lines is due to dust-absorption.

Attribution of the break to dust absoprtion can also be affected by the presence of foreground extinction, and the finite thickness of the ionization layer (see AMH16 for more details on this). In particular, if the thickness of the ionization layer is less than the physical depth at which optical depth reaches unity, then it would not be correct to conclude that the break is due to dust-absorption.

10 Discussion

10.1 Important results and their interpretation

In this paper, we showed how dust absorption affects different statistical methods to study turbulence. Our results suggest that VCA, VCS and velocity centroids are able to trace velocity spectra even in the presence of dust absorption. However, it might not be possible to obtain density spectra with these techniques, especially if dust absorption is strong and the resolution of telescope is poor.

In the presence of dust absorption, power spectrum of intensity fluctuations in the thin-slice regime has spectral index greater than -2 at small $K$, while at large $K$, it has spectral index less than -2. Thus it is appropriate to use correlation function at large $R$, while structure function is appropriate at small $R$. Our results suggest that density fluctuation is unimportant if dust absorption is strong as one effectively samples signals from only optical optical depth of unity, or equivalently from a fixed column-density. This allows one to obtain cleaner velocity statistics at larger scales. To study velocity and density spectra, one can also use power spectrum of intensity fluctuations. In the language of power spectrum of the intensity fluctuations, the spectral slope changes by unity444This is true for only if density spectra is steep or if dust absorption is strong. For shallow density in the weak absorption limit, spectral slope changes by $1-\nu_{\rho}$ (see Sec. 4.2)., at the critical break point $K_{c}=\Delta_{0}^{-1}$. This is a manifestation of the fact that at the lags $R$ associated with these wavenumbers, the depth along LOS is much smaller than the separation between the lines of sight, thus statistics are effectively two dimensional. The break in the slope near dust cut-off scale is also exhibited by centroids. Thus, study of velocity spectra is possible using both centroids and thin velocity slice regime of VCA even in the presence of strong dust absorption, albeit intensity signal will be diminished if the dust absorption is strong.

On the other hand, one can only study density spectra if absorption is not strong by either using structure function at short scales $R$ or by using correlation function at large scales $R$. Our calculations show that the correlation length set by dust absorption is short for shallow density and long for steep density. Thus, it is appropriate to study steep density using structure function, while correlation function is appropriate for shallow density. An important limitation is that the amplitude of correlation function is exponentially suppressed, and thus it is challenging to measure this correlation especially if instrumental noise is not properly factored out. We stress that in the thin slice regime of VCA, density effects are not important if dust absorption is strong, and only velocity statistics is manifested in the intensity statistics. This allows one to cleanly obtain velocity spectra even if the information of density spectra is not available.

Our results in this paper show that dust absorption effectively decreases the extent of self-absorption. This means that the universal regime $K^{-3}$ of intensity statistics, originally obtained in LP04, is suppressed by strong dust absorption; in particular, the range of scale where one can see universal regime is decreased by a factor of $\Delta_{0}/S$. This is due to the fact that the cut-off in physical space introduced by dust absorption decreases the velocity window introduced by self-absorption, thus suppressing the universal regime.

Our results also show that one can study intensity anisotropies even in the presence of dust absorption. However, in the context of VCA, dust absorption might lead to isotropization of statistics if proper statistical measure is not used. In particular, dust absorption leads to isotropic structure function at lags $R>\Delta_{0}$, where one needs to use correlation function in order to study both spectrum and anisotropy of intensity fluctuations. Similarly, one needs to study both spectrum and anisotropy at small $R$ for thick velocity slice limit of steep density. This ambiguity of choice of statistical measure is resolved if one studies anisotropy using power spectrum of intensity fluctuations. We stress that our description of anisotropy of interstellar turbulence uses the theory of strong non-linear MHD turbulence. The theory, including the issues of anisotropies, has been studied extensively both theoretically and numerically, in particular in the papers by Cho & Lazarian (2002, 2003) and Kowal & Lazarian (2010). In these papers the anisotropies of MHD turbulence were studied for a wide range of sonic and Alfvén Mach numbers. These anisotropies are not expected to depend on the ionization of the media, provided that the scales under study are larger than the scale of the decoupling of ions and neutrals (see Xu et al. 2016). Therefore, our description of anisotropies applies to both turbulence within ionized gas in massive star forming regions and diffuse multiphase interstellar medium. Moreover, we study the anisotropies arising from the sub-Alfvénic and trans-Alfvénic turbulence.

In the context of VCS, we showed that both narrow beam and wide beam asymptote remain the same even in the presence of dust absorption. However, an important consequence of dust absorption is that the effective width of telescope beam $\Delta B$ increases, and one will be able to obtain narrow beam only if this width is less than the dust cut-off, i.e. $\Delta B<\Delta_{0}$. This might be difficult to achieve if dust absorption is strong.

Finally, an important piece of work that was carried out was understanding how different analytical techniques work with collisionally excited emission lines. Our results show that overall, all the analytical expressions and asymptotic predictions obtained in the previous studies are applicable to collisionally excited emission lines. However, an interesting regime is when dust absorption is present. In this case, density effects are not suppressed for collisionally excited emission lines even when dust absorption is strong.

In conclusion, the previous techniques of studying turbulence applied to mostly radio lines. Our study suggests that these techniques can be applied even to the optical and ultra-violet (UV) lines for which dust absorption can not be neglected.

10.2 Model Assumptions

In carrying out the analysis in this paper, we used several assumptions. Firstly, we considered the case when emission is proportional to the density of emitters, which excludes the case of collisionally excited emission, where emission proportional to the square of density of emitters. Secondly, we assumed dust density to be proportional to the gas density. Correlation between HI and 100 micron (see Boulanger et al. 1996) justifies this assumption for the gas clouds in the ISM. Thirdly, we assumed that the dust extinction coefficient $\kappa$ is independent of frequency within the frequency range of Doppler broadened spectral line. While this is not strictly true, it does not affect our analysis in any significant manner, especially in the case when thermal broadening of intensity profile is not significant. It is important to note that the dust extinction coefficient is different for different emission lines.

Although density is better modelled as log-normal field, we argued that the column density can be approximated with reasonable accuracy to be Gaussian. For steep density, the density dispersion is small, therefore density and thus column-density can be modelled as Gaussian. For shallow density, we showed in Appendix B that if the LOS depth is larger than the correlation length of the field, column-density tends to be Gaussian.

To study anisotropy of intensity correlations, we explicitly used models of different MHD modes, and obtained anisotropy for each of these modes. Description of turbulence in terms of MHD cascades is only applicable if these MHD modes approximately form independent cascade. Numerical study carried out by Cho & Lazarian (2002) shows that the coupling between various MHD modes is marginal as long as sonic Mach number is not too high, thus validating the applicability of this decomposition. Some possibilities of the validity the MHD modes in describing turbulence in the diffuse ISM has been questioned by Caldwell et al. (2017) by considering the recent dust polarization result obtained by Planck observations. By assuming MHD model, Caldwell et al. (2017) have shown that there is a very narrow range of parameters that could explain the Planck’s result. Thus, they raise the possibility that at large scales ($\sim$ pc), turbulence could either be unimportant, or the ISM is not described by MHD model that is used. However, our ongoing work shows that the description of dust polarization using anisotropic MHD model is not inconsistent with the Planck result.

10.3 Other techniques to study turbulence

In this paper, we mainly focused on how dust absorption affects VCA, VCS and centroids. However, there are other techniques to study turbulence, and these are also likely to be affected by dust absorption. An important technique to study turbulence using velocity slice of PPV space is the spectral correlation function (SCF, see Rosolowsky et al. 1999). SCF is very similar to VCA if one removes the adjustable parameters from it, as both of these technique measure correlations of intensity in velocity slices of PPV. However, SCF treats outcomes empirically. The description of how intensity is affected by dust absorption that was presented in this paper is likely to help us understand the effect of dust absorption on SCF.

There also exist numerous techniques to identify and analyze clumps and shells in PPV (see Houlahan & Scalo 1992; Williams et al. 1994; Stutzki & Guesten 1990; Pineda et al. 2006; Ikeda et al. 2007). Moreover, a more advanced technique to study a hierarchical structure of the PPV, namely dendrogram technique (Goodman et al., 2009), can provide a complementary insight to the values of sonic and Alfvén Mach numbers (Burkhart et al., 2013). Synergy of these different available techniques is very important to obtain a comprehensive picture of MHD turbulence.

10.4 Spectroscopic and synchrotron studies of MHD turbulence

In LP16, a number of techniques were suggested aimed at obtaining the information about magnetic field and the density of cosmic electrons using synchrotron intensity and polarization measurements. These techniques use position-position-frequency (PPF) cube (an analogue of PPV cube used in the VCA). Such studies provide a complimentary method to study the mean magnetic field direction, the degree of magnetization of the media and the contribution of Alfvén, slow and fast modes in a turbulent cloud.

Understanding the properties of turbulent cascade in different interstellar phases is essential to understand the dynamics of the ISM as well as the transport process of heat and cosmic rays. While the turbulence in HI and molecular gas is well studied via spectroscopic measurements, synchrotron emission is useful in studying diffuse hot and warm turbulent media. The correspondence of the different important parameters characterizing MHD turbulence in different interstellar phases would testify a single turbulent cascade on the galactic scale, which can be a finding with important consequences for different branches of astrophysical research, e.g. for cosmic ray physics (see Schlickeiser 2002). Since dust absorption is important at shorter wavelengths, it is less likely to affect synchrotron studies. Thus, synergy of synchrotron studies with spectroscopic studies can be beneficial to understand magnetic properties of a turbulent medium, even when it is significantly dusty.

11 Summary

On the basis of the analytical theory of the intensity fluctuations in the PPV space developed in the earlier works (LP00; LP04; LP06; KLP16) we provided the analytical description of how the results of those works are affected by dust absorption. The results of our study within the assumptions that were used can be summarized as follows:

• •

Simple model of dust absorption was used to show how dust absorption affects study of VCA. In particular, it was shown that it is possible to obtain velocity spectra in the thin slice regime of VCA even if dust absorption is present. However, the spectral index of power spectrum of intensity changes by one at low wavenumbers to high wavenumbers. This change corresponds to studying the statistics at the scales comparable to the thickness of the dust layer with the optical depth of unity.

• •

It was shown that the velocity spectra can be obtained using centroids in the presence of dust absorption. Similar to the VCA, the centroids, due to dust absorption, exhibit the break in the spectral slope of centroid structure function. The reason for that is the same as in the case of the VCA, namely the transfer to studying 2D surface layer of turbulence.

• •

The study of turbulence using VCS in the presence of dust absorption shows that the asymptotes at different wavenumbers $k_{v}$ for both poor and high resolution beam remains the same as in the case of negligible dust absorption. One can state therefore, that the technique is robust in relation to the dust absorption.

• •

Anisotropies measured by both VCA and centroids at small 2D lags $R$ are unaffected by dust absorption, whereas at larger lags one needs to use proper statistical measure (either correlation function or power spectrum) in order to study anisotropies.

• •

In the presence of both dust absorption and self-absorption of the gas, the statistics of the observed total intensities and velocity centroids change. In particular, the universal regime of fluctuations $\sim K^{-3}$ that is independent from both underlying velocity and density fluctuations disappears when dust cut-off is much smaller than the LOS extent of the cloud. This allows study of velocity spectra even when self-absorption is present.

• •

All the analytical formulations of VCA, VCS and centroids are valid even with collisionally excited emission lines. However, in the regime of strong dust absorption, density effects are not suppressed for collisionally excited emission lines.

• •

Our theoretical expectations correspond to the observations of turbulence in HII regions and explain the change of the spectral slope that is reported for studies of fluctuations at large scales.

12 Acknowledgements

D. K. and D. P. thank the Institut Lagrange de Paris, a LABEX funded by the ANR (under reference ANR-10-LABX-63) within the Investissements d’Avenir programme under reference ANR-11-IDEX-0004-02, where this work was started. A. L. acknowledges the NSF grant AST 1212096 and Center for Magnetic Self Organization (CMSO). He also acknowledges Institut D’Astrophysique de Paris and Institut Lagrange de Paris for hospitality during his visit.

Appendix A Full expression of intensity correlation

We now focus on a simple case for which self-absorption is negligible. In this situation, Eq. (2) gives

[TABLE]

which can be equivalently written as

[TABLE]

Assuming uncorrelated density and velocity field, the mean spectral intensity is given by

[TABLE]

Assuming $\int_{0}^{z}\mathop{\mathrm{d}z^{\prime}}\delta\rho(\bm{X},z^{\prime})$ to be Gaussian random quantity (which holds good if the density fluctuations are small), Eq. (112) can be written as

[TABLE]

The correlation of spectral intensity is given by

[TABLE]

We assume fluctuations of column density to be Gaussian, so that

[TABLE]

where

[TABLE]

With this, one can immediately write

[TABLE]

and

[TABLE]

Using Eqs. (A), (115) and (A), one can immediately write

[TABLE]

where $\langle\rho(\bm{x}_{1})\rho(\bm{x_{2}})\rangle\equiv\rho_{0}^{2}+\xi_{\rho}(\bm{x}_{1},\bm{x}_{2})$.

Appendix B Statistics of column-density

The column density

[TABLE]

is an important parameter, especially at thick slice limit of the VCA, where it is proportional to the intensity of emission. We first write several important statistical descriptors of $\mathcal{N}$. The mean column-density is given by

[TABLE]

while the variance of the column-density fluctuations $\mathcal{N-\bar{N}}$, $\sigma_{\mathcal{N}}^{2}=\int_{0}^{S}\mathrm{d}z_{1}\int_{0}^{S}\mathrm{d}z_{2}\,\xi_{\rho}(0,z_{1}-z_{2})$ is evaluated to

[TABLE]

The correlation function of column-density fluctuations $\mathcal{N-\bar{N}}$ is given by

[TABLE]

which for shallow density at $R>r_{c}$ evaluates to

[TABLE]

The structure function of the column density, which can be computed as

[TABLE]

for density with $\nu_{\rho}<1$ behaves at $R<r_{c}$ as

[TABLE]

Note that since the column density is a projected 2D field, its behaviour reflects the shifted value $\nu_{\rho}-1$ of the 3D density index. Correspondingly, the correlation function $\xi_{\mathcal{N}}(\mathbf{R})$ is applicable for large separation analysis only for $\nu_{\rho}>1$ rather than $\nu_{\rho}>0$ as for 3D density field. Similarly, one should use the structure function $d_{\mathcal{N}}(\mathbf{R})$ for $\nu_{\rho}<1$ to study scaling regime through small $R$, while too steep underlying 3D density with $\nu_{\rho}<-1$ will lead to saturated structure function statistics and needs to be studied separately. For instance, Kolmogorov index for 3D density, $\nu_{\rho}=-2/3$ will result in steep $5/3$ scaling of the column density structure function.

We can estimate the correlaton length as the scale where $d_{\mathcal{N}}(R)$ or $\xi_{\mathcal{N}}(R)$ reaches the value of the variance, extrapolating Eqs. (125) and (127). The second form of Eqs. (125) and (127) shows that the correletion length of the column density is given by the underlying density scale $r_{c}$ for $\nu_{\rho}>1$ but is increased to $S$ due to accumulation of slow falling long wave correlations for $\nu_{\rho}<1$. Seemingly abrupt jump in column density correlation scale at $\nu_{\rho}=1$ is in reality a smooth transition which is not recovered by our asymptotic analysis. While correlation lengths similar to the size of the cloud $S$ are expected for $\nu_{\rho}<0$, note that in the $0<\nu_{\rho}<1$ range, the column density may have considerably larger correlation scale than underlying 3D density which is still shallow and is expected to have $r_{c}<S$.

Let us turn to the kind of statistical model appropriate for column density. For a steep density spectra, fluctuations are smaller than the mean density, and Gaussian approximation for density fluctuations works well which translates to Gaussian behaviour of the column densty. However, for shallow density, fluctuations can be larger than the mean density, and thus density fluctuations follow lognormal distribution better than the Gaussian distribution. If so, what kind of distribution does column density follow? We answer this question numerically. For that we generate random correlated densities in one dimension. The most convenient way to do this is to start in spectral domain, and generate Gaussian random amplitudes $a(k)$ in spectral domain with dispersion given by their power spectra, which is related to the correlation function $\xi(\bm{r})$ as

[TABLE]

where $n$ is the dimensionality of space, and in our test case $n=1$. We define $P(k)$ to satisfy power-law

[TABLE]

where $r_{c}$ is the correlation length of the density field, and is in general very small for a shallow density. The advantage of working in spectral domain is that these random amplitudes are uncorrelated. The Fourier transform of these random amplitudes corresponds to $\log\rho$ in real space:

[TABLE]

With the obtained $\rho(r)$, we find out $\mathcal{N}$, and this is done for several realisation, so as to obtain distribution of $\mathcal{N}$. As shown in Fig. 14, our results show that $\mathcal{N}$ is to a very good extent Gaussian. One can understand this Gaussian behaviour by considering the fact that if the integration length $S$ is much larger than the correlation length $r_{c}$, then at LOS lengths larger than $r_{c}$, the fluctuations are essentially independent. Thus, from Central limit theorem, the sum of fluctuations should tend to be Gaussian.

Appendix C Statistics of dust cutoff Delta

The main effect of dust absorption is the introduction of cut-off scale $\Delta(\bm{X})$ beyond which signals are effectively cut out. This cut-off scale is given by

[TABLE]

$\Delta(\bm{X})$ is fluctuating random quantity, since $\rho(\bm{r})$ is fluctuating. In this appendix we derive several important statistical properties for $\Delta$ assuming lognormal distribution for $\rho$. We model dust density $\rho$ as

[TABLE]

where $\mu(r)$ is a random Gaussian field with zero mean and variance $\sigma_{\mu}^{2}$, and is described by the two-point probability function

[TABLE]

where $\mu_{1}=\mu(\bm{x}_{1})$, $\mu_{2}=\mu(\bm{x}_{2})$ and $\xi_{\mu}(|\bm{x}_{1}-x_{2}|)$ is the correlation function of $\mu$ field.

Using the following relations

[TABLE]

as well as

[TABLE]

we obtain for the moments of the dust density distribution

[TABLE]

and for the density correlation function

[TABLE]

Non-Gaussianity of the density field is reflected in (higher order moments).

It can also be shown that the correlation of square of density is

[TABLE]

For dust absorption we shall consider two regimes, depending on comparison of the absorption length $\Delta({\bm{X}})$ and the density perturbations correlation length $r_{c}$.

First regime is $\Delta({\bm{X}})<r_{c}$. A weak, but explicit necessary condition for this to hold is is $\kappa\bar{\rho}r_{c}\gg 1$. This regime is primarily relevant for the case of steep spectrum of density, when $r_{c}$ is associated with the largest scale of the inertial range. Since density can be thought as almost constant within $r_{c}$ length scale, we have

[TABLE]

With the help of Eqs. (132), (133) and (141), statistical properties of the absorption length are immediately available. The mean and the second moment of $\Delta$ are

[TABLE]

therefore relative fluctuations $\delta\Delta/\langle\Delta\rangle$ have the variance equal to that of density fluctuations

[TABLE]

The correlation function is

[TABLE]

and

[TABLE]

An interesting observation that can be made from Eq. (142) is that the mean depth $\bar{\Delta}$ is longer than $1/(\kappa\bar{\rho})$.

Second regime is $\Delta(\bm{X})>r_{c}$, which is particularly relevant for shallow spectra, which has $r_{c}$ much smaller than the size of the cloud. For shallow spectra density is uncorrelated at scales larger than the correlation length. Formally, given a distribution of $\rho$, we are interested in the distribution of $\Delta$ such that the following holds

[TABLE]

In this situation, the accumulation of optical depth along the LOS can be described by Brownian random walks with a net drift.

A formal way to obtain the distribution of physical depths $\Delta$ is via diffusion approximation. Let $\Pi(\tau_{n},\Delta_{n})$ be the distribution of paths for which the threshold optical depth of 1 has not been met. In this picture $\delta\tau\equiv\tau_{n}-\tau_{n-1}=\kappa\rho r_{c}$ and $\delta\Delta=\Delta_{n}-\Delta_{n-1}=r_{c}$, where $r_{c}$ is the minimum length at which correlation of random jumps at each step can be neglected. The case of constant density implies $r_{c}\rightarrow 0$. The distribution function $\Pi$ at any optical depth depends upon the distribution at the immediately prior optical depth, and hence satisfies Chapman-Kolmogorov equation

[TABLE]

where $P(\tau_{n},\Delta_{n}|\tau_{n-1},\Delta_{n-1})$ is the transition probability to have a field value of $\tau_{n}$ at $\Delta_{n}$, given the field value of $\tau_{n-1}$ at $\Delta_{n-1}$. If higher moments of $\Pi$ are ignored, Eq. (148) leads to Fokker-Plank equation

[TABLE]

The initial condition that $\Pi$ satisfies is $\Pi(\tau,0)=\delta_{D}(\tau)$, where $\delta_{D}(\tau)$ is the Dirac delta function, and the boundary condition satisfied by $\Pi$ is $\Pi(1,\Delta)=0$.

Instead of solving for complete form of probability distribution $\Pi$, it is convenient to introduce first passage length (an analogue of first passage time in random walks), which is the physical length required to reach optical depth of unity for the first time. We denote this length to be $\Delta(\tau)$, and the probability density associated with it to be $f(\tau,\Delta)$, which gives the probability that optical depth $\tau$ is reached at physical depth $\Delta$. In order to calculate $f(\tau,\Delta)$, we need to solve Eq. (149) with proper boundary conditions.

We denote the first passage distance from $\tau_{0}$ to $\tau=1$ by $\Delta(\tau_{0})$. In this scenario, it is useful to introduce the notion of survival probability $S(\tau_{0},\Delta)$, which is the probability that at physical depth $\Delta$, optical depth of unity has not been reached yet, assuming that the motion started at optical depth of $\tau_{0}$. The survival probability satisfies the following relations

[TABLE]

We are now interested in moments of physical depths, i.e. $\langle\Delta^{n}\rangle$, which is given by

[TABLE]

which was obtained by integration by parts and assuming $S(\infty)-S(0)=0$ for a well-defined probability distribution function $S(\Delta)$. Eq. (C) can also be written as

[TABLE]

where

[TABLE]

Integrating by parts Eq. (153) $n$ times and using Fokker Plank equation (Eq. (149)) finally gives differential equation governing moments of $\Delta$

[TABLE]

satisfying the initial and boundary conditions The moments satisfies the following boundary condition

[TABLE]

[TABLE]

Eq. (155) can be solved for $\Delta_{1}(0)$ and $\Delta_{2}(0)$, to obtain $\langle\Delta\rangle$ and $\langle\Delta^{2}\rangle$. Here, we quote the final result

[TABLE]

and

[TABLE]

In the above expression, $r_{c}$ is the minimal step size at which correlation between subsequent steps can be neglected. Clearly in the case when $r_{c}\rightarrow 0$, the variance of $\Delta$ is square of its mean.