A&A 404, 557-567 (2003)
DOI: 10.1051/0004-6361:20030480

Measurements of the interstellar turbulent plasma spectrum of PSR B0329+54 using multi-frequency observations of interstellar scintillation

V. I. Shishov1 - T. V. Smirnova1 - W. Sieber2 - V. M. Malofeev1,5 - V. A. Potapov1 - D. Stinebring3 - M. Kramer4 - A. Jessner5 - R. Wielebinski5


1 - Pushchino Radioastronomy Observatory of Lebedev Physical Institute, 142290 Pushchino, Russia, and Isaac Newton Institute of Chile, Pushchino Branch
2 - Hochschule Niederrhein, Reinarzstr. 49, 47805 Krefeld, Germany
3 - Oberlin College, OH 44074, Oberlin, USA
4 - University of Manchester, Jodrell Bank Observatory, Macclesfield SK11 9DL, UK
5 - Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany

Received 2 July 2002 / Accepted 21 March 2003

Abstract
Interstellar scintillation multi-frequency observations of PSR 0329+54 in the frequency range from 102 MHz to 5 GHz were analysed to estimate the spectrum of interstellar plasma inhomogeneities in the direction of this pulsar. Based on the theory of diffractive scintillation, the composite structure function of phase fluctuations covering a large range of turbulence scales was constructed. We found that the spectrum is well described by a power law with n = 3.5 for scales from 106 to 109 m, which differs from the value known for a Kolmogorov spectrum. We can, however, within the accuracy of our data not exclude a Kolmogorov spectrum. It became also clear that angular refraction of emission must be taken into account to fit the data points at all observing frequencies. The size of the irregularities responsible for the angular refraction is estimated to be about  $3\times 10^{13}$ m. They can be identified with clouds of neutral hydrogen that can be considered as holes of electron density.

Key words: stars: pulsars: general - turbulence - ISM: structure - stars: pulsars: individual: PSR B0329+54

1 Introduction

A generally accepted point of view at the present time is that scattering as well as diffractive and refractive interstellar scintillation effects are caused by electron density fluctuations in the interstellar medium (ISM) (Scheuer 1968; Rickett 1969; Sieber 1982; Rickett et al. 1984; Armstrong et al. 1995; Stinebring et al. 2000). It could be shown that for the different kinds of pulsar observations connected with the propagation through the interstellar plasma, a composite structure function of phase fluctuations can be constructed and that this structure function follows a power law over a very wide range of scales (106 to 1013 m). This function fits the experimental data quite well for the near ISM ($\le$1 kpc) (Armstrong et al. 1995) as well as for the distant ISM (>1 kpc) (Shishov & Smirnova 2002). The exponent of the structure function is about 1.7. This means that the 3-dimensional spatial spectrum of electron density can be described by a Kolmogorov spectrum

 
$\displaystyle % \Phi_{N_{\rm e}} (q) = C_{N_{\rm e}}^2 \vert\vec{q}\vert^{-n} , \quad n = 11/3$     (1)

with the power law index n = 11/3, where $\Phi_{N_{\rm e}}$ denotes the squared Fourier spectrum of the electron density fluctuations. $C_{N_{\rm e}}^2$ characterizes the level of plasma turbulence and $\vec{q}$ is the spatial frequency.

Although the Kolmogorov spectrum describes the data sufficiently well in a statistical sense, the dispersion of points is large and in particular directions the spectrum can differ from a Kolmogorov one. Smirnova et al. (1998) showed that two types of turbulent spectra can exist in the ISM as measurements of the flux density variations of 21 pulsars at 610 MHz prove, which show two different types of behaviour for the dependence of the modulation index on the dispersion measure DM. The data can be explained if one group follows a pure power law spectrum (Eq. (1)) with $n\simeq 11/3$, whereas the other group corresponds to a medium with a power law spectrum pieced together by an exponent of $n\simeq 11/3$ and an inner scale of length $l_i \simeq 3\times 10^{10}$ cm and an exponent n1 > 4 at higher space frequencies:

 
$\displaystyle % \begin{array}{lcl} \Phi_{N_{\rm e}} (q) \propto \vert\vec {q}\v... ... e}} (q) \propto \vert\vec {q}\vert^{-n_1}, & n_1 > 4, & q > 1/l_i. \end{array}$     (2)

It should be noted that the accuracy of the determination of the value of li is low and that the shape of the turbulence spectrum in the region around li is unknown. This means that the shape of the spectrum can be different for different regions in the sky and that the variations of the turbulence spectrum must be investigated in dependence on these directions.

In general, there are two types of power law spectra known for turbulent plasmas: in addition to the aforementioned Kolmogorov spectrum with $n_{\rm K} = 11/3$ (Tu et al. 1984) the spectrum of weak plasma turbulence with $n_{\rm wpt} = 3.5$(Iroshnikov 1963; Kraichnan 1965). The difference between $n_{\rm K}$ and $n_{\rm wpt}$ is unfortunately small, about 5%, whereas the accuracy of the measurements of the value of n is not better than 10%. New measurements with much higher accuracy are therefore urgently needed.

In this paper we propose a new method for the determination of the turbulent spectrum in the direction of a given pulsar using multi-frequency scintillation observations covering a wide frequency range from centimeter to meter wavelengths. We will then address the question whether the exponent of the power law spectrum equals a Kolmogorov spectrum or not.

2 Basic equations

In the theory of wave propagation through random media the fundamental function is the gradient of the phase structure function computed in the geometrical optics approximation (Prokhorov et al. 1975). This function is given for a plasma by the equation (Shishov & Tokumaru 1996):

 
$\displaystyle % D (\vec{\rho}) = 4 \pi (\lambda r_{\rm e})^2 \int {\rm d}^2 \ve... ...ho}\right)\right] \Phi_{N_{\rm e}} \left(\vec{q_{\perp}}, q_{\parallel}\right),$     (3)

where $\vec{\rho}$ denotes a two-dimensional vector between two points in the plane normal to the line of sight, $\lambda$ is the wavelength, $r_{\rm e}$ the classical electron radius, $\vec{q}$the spatial frequency, $q_{\parallel}$ its component along the line of sight, and $\vec{q_{\perp}}$ its normal component. For a power law spectrum (1) $D(\vec{\rho})$ is described by the equation (Smirnova et al. 1998):
 
$\displaystyle % \begin{array}{lcl} D(\vec{\rho}) & = & A (n) (\lambda r_{\rm e}... ...splaystyle\frac{2^{4-n} \pi^3}{\Gamma^2(n/2) \sin (\pi n / 2)}\cdot \end{array}$     (4)

$D(\vec{\rho})$ describes the turbulence spectrum $\Phi_{N_{\rm e}} (\vec{q})$ for the spatial frequency region near $q = 1/\rho$. The equation is correct when $C_{N_{\rm e}}^2$ and n varies smoothly with $\vec{q}$; it describes also spectra determined by Eq. (2) with $\rho > l_i$. If n1 is however greater than 4, $D(\vec{\rho})$ becomes proportional to $\rho^2$for $\rho > l_i$ so that it is impossible to obtain information on the turbulence spectrum in the region q > 1/ li.

When the displacement of the line of sight is determined by a movement of the source with the velocity $V_{\perp}$, which is assumed to be greater than internal motions within the medium or the observer velocity, the structure phase function is given by

 
$\displaystyle % D_{\rm S}(t) = \int\limits_o^R {\rm d}r D[V_{\perp} t(1-r/R)],$     (5)

where R is the distance from the observer to the source and r is the variable distance from the source. $D_{\rm S}(t)$ can be measured directly using pulsar timing data (Cognard & Lestrade 1997) for large time intervals t. The structure function of the variations of the residuals $\delta \tau$ due to propagation in the turbulent medium, $D_{\rm\delta\tau, turb} (t)$, is reduced to $D_{\rm S}(t)$ by the relation
 
$\displaystyle % D_{\rm S}(t) = (2\pi f)^2 D_{\rm\delta \tau, turb} (t),$     (6)

where f denotes the radio wave frequency. $D_{\rm S}(t)$ can be measured for small values of t using data on intensity variations in the saturated scintillation regime. For this regime the correlation function of flux variations I(t) is described by the equation
 
$\displaystyle % B_{\rm I}(t) = <I>^2 \exp[- D_{\rm S} (t)].$     (7)

If t0 is the characteristic scale of intensity variations we have
 
$\displaystyle % D_{\rm S} (t_0) = 1$     (8)

and Eq. (7) can be reduced for $t \ll t_0$ to the equation
 
$\displaystyle % B_{\rm I}(t) \simeq B_{\rm I}(t=0) - <I>^2 D_{\rm S} (t),$     (9)

and one obtains
 
$\displaystyle % D_{\rm S}(t) \simeq \frac{B_{\rm I}(0) - B_{\rm I}(t)}{<I>^2} \simeq \frac{1}{2} \frac{D_{\rm I}(t)}{<I>^2}, \quad t \ll t_0,$     (10)

where $D_{\rm I}(t)$ is the structure function of the intensity fluctuations.

The phase structure function can be determined from observations of intensity scintillation in the weak scintillation regime. In this case the intensity correlation function of intensity fluctuations is determined by the equation (Prokhorov et al. 1975; Malofeev et al. 1996)

 
$\displaystyle % \begin{array}{lcl} B_{\rm I}(t)\; & \simeq & \;8 \pi (\lambda r... ...vec{q_{\perp}} {\displaystyle{\frac{R}{r}}},q_{\parallel}=0\right), \end{array}$     (11)

where $k = 2\pi/\lambda$ is the wave number. For small values of t we may write
 
$\displaystyle % t \ll t_{\rm Fr} = \frac{(R/k)^{1/2}}{V_{\perp}}\cdot$     (12)

A small value of t corresponds to a large value of $q_{\perp} \gg (t_{{\rm Fr}} V_{\perp})^{-1}$ so that we may replace $\sin^2 [q_{\perp}^{2} R (R-r)/2kr]$ by its average value 1/2; it follows that 0.5 times the structure function of relative intensity fluctuations (1/2) DI (t) / <I>2 equals the phase structure function DS (t) (see, for example, Tatarskii 1971). Equation (10) may therefore be used for the determination of the phase structure function in the weak scintillation regime.

$D_{\rm S}(t)$ may be measured at any frequency f; if a fixed frequency f0 is given, the result may be converted by use of the universal factor $\left( \frac{\lambda_0}{\lambda} \right)^2 = \left( \frac{f}{f_0} \right)^2$

 
$\displaystyle % D_{\rm S}(t,f_0) = \left( \frac{f}{f_0} \right)^2 D_{\rm S}(t,f).$     (13)

It should be mentioned that $D_{\rm S}(t,f)$ is valid for a given time t, which must be kept constant when $D_{\rm S}(t,f)$ is converted to  $D_{\rm S}(t, f_0)$.

An estimation of the value of the phase structure function can be obtained from the scintillation index of refractive scintillation $m_{\rm ref}$ (Smirnova et al. 1998). The refractive scintillation index of strong scintillation is described in case of isotropic turbulence by

 
$\displaystyle % m_{\rm ref}^2 = \lambda^4 r_{\rm e}^2 R^2 \int \limits_0^{\inft... ... \Phi_{N_{\rm e}}\left[q_{\perp} \frac{R}{r}, q_{\parallel} = 0\right] \exp(-F)$     (14)


 
$\displaystyle % F = \int\limits_0^r {\rm d}r' D\left[\frac{q_{\perp} r'}{k} \le... ...right] + \int\limits_r^R {\rm d}r' D\left[\frac{q_{\perp}}{k} (R - r')\right] .$     (15)

The refractive scintillation index is determined by the value of the turbulence spectrum $\Phi_{N_{\rm e}}$ in the region of spatial frequencies near to $1/b_{\rm ref}$, where $b_{\rm ref}$ is the characteristic spatial scale of refractive scintillation. $m_{\rm ref}^2$ can be expressed as (Shishov & Smirnova 2002)
 
$\displaystyle % m_{\rm ref}^2 = A(n) D_{\rm S} (V_{\perp} T_{\rm ref}) \left( \frac{t_{\rm dif}}{T_{\rm ref}} \right)^2.$     (16)

Here $T_{\rm ref}$ and $t_{\rm dif}$ are the characteristic temporal scales of refractive and diffractive scintillation respectively, A(n) is a numerical coefficient, which depends on the power law index n and on the distribution of the medium along line of sight. In the case of a statistically homogeneous medium with $4 - n \ll 1$, A (n) is given by
 
$\displaystyle % A(n) \simeq 6\times (4-n), \quad 4-n \ll 1.$     (17)

The phase structure function can also be measured using observations of the frequency structure of intensity variations.

Ostashov & Shishov (1977) showed that the two-frequency field coherence function of a spherical wave $B_{\rm E}(\Delta f)$ in the case of a power law turbulence spectrum and for small values of the frequency difference

 
$\displaystyle % \Delta f \ll \Delta f_{\rm dif}$     (18)

can be described by the equation

\begin{eqnarray*}B_{\rm E} (\Delta f, \rho=0) \simeq <I> \biggl[1 + \int\limits_... ...r, \vec {\rho'}) D\left(\frac{r}{R} \vec {\rho'}\right) \biggl] \end{eqnarray*}


with
 
$\displaystyle % G (R,r,\vec {\rho'}) = \displaystyle\frac{k (r/R)}{2 {\rm i} \p... ...e\frac{k (r/R) (\vec {\rho'})^2}{2{\rm i} (\Delta f / f) (R - r)} \right\}\cdot$     (19)

Equation (19) is the first order approximation of an iteration series resulting from an expansion of the successive scattering. The first term is much smaller than the zero term, which equals to 1. In Eq. (19) i is the imaginary unit and $\Delta f_{\rm dif}$ the characteristic scale of the frequency structure of diffractive scintillation (for more details see Smirnova et al. 1998). Using Eq. (19) one may derive the relation which determines the frequency correlation function of intensity fluctuations
 
$\displaystyle % \begin{array}{l} B_{\rm I,dif} (\Delta f, \rho=0) = B_{\rm E} (... ... G (R,r, \vec {\rho'}) D\left(\frac{r}{R} \vec {\rho'}\right)\big]. \end{array}$     (20)

Equation (20) can be reduced to the following equation in case of a power law turbulence spectrum
 
                                       $\displaystyle B_{\rm I,dif} (\Delta f, \rho=0) \simeq~ <I>^2 \left[1 - D_{\rm S,1} \left(\rho_{{\rm dif}} (\Delta f)\right)\right]$  
    $\displaystyle D_{\rm S,1} \left(\rho_{{\rm dif}}(\Delta f)\right) = \int\limits... ...)^{1/2} \left( 1 - \frac{r}{R} \right)^{1/2} \rho_{{\rm dif}} (\Delta f)\right]$  
    $\displaystyle \rho_{{\rm dif}}(\Delta f) = \left[\Gamma \left( \frac{n}{2} \rig... ...{4 \pi}\right]^{1/(n-2)} \left( \frac{\Delta f}{f} \right)^{1/2} \rho_{\rm Fr},$  
    $\displaystyle \rho_{\rm Fr} = \left(\frac{R}{k} \right)^{1/2}\cdot$ (21)

Using these relations one can estimate the phase structure function
 
$\displaystyle % D_{\rm S,1} (\Delta f) = D_{\rm S,1} \left(\rho_{{\rm dif}} (\D... ...f)\right) \simeq \frac{B_{\rm I,dif}(0) - B_{\rm I,dif} (\Delta f)}{<I>^2}\cdot$     (22)

Equations (21) and  (22) are correct for $D_{\rm S,1} (\Delta f) \ll 1$. The difference between $D_{\rm S}$and $D_{\rm S,1}$ depends on the distribution of the turbulent medium along the line of sight. In case of a statistically uniform distribution one obtains
 
$\displaystyle % \begin{array}{rcl} D_{\rm S} \left(\rho_{{\rm dif}}(\Delta f)\r... ...&\simeq &(n-2) D_{\rm S,1} \left(\rho_{{\rm dif}}(\Delta f)\right). \end{array}$     (23)

One can reduce the measured value $D_{\rm S,1} (\Delta f)$ to a given frequency f0 by use of the relations

 \begin{displaymath}% D_{\rm S,1} \left(\rho_{{\rm dif}}(\Delta f), f_0\right) = ... ...ight)^2 D_{\rm S,1} \left(\rho_{{\rm dif}}(\Delta f), f\right) \end{displaymath} (24)


 \begin{displaymath}% \Delta f (f_0) = \left( \frac{f_0}{f} \right)^2 \Delta f (f). \end{displaymath} (25)

It should be pointed out here that the value of $D_{\rm S,1}$ as well as the value of the frequency scale $\Delta f$ must be reduced to the same frequency f0.

Another type of relation between the frequency correlation function $B_{\rm I}(\Delta f)$ of the intensity fluctuations and the spatial phase structure function $D_{\rm S}(\rho)$ is realised in case of strong angular refraction. According to (Shishov 1973) we may introduce the accumulated angle of refraction $\Theta_{\rm ref}$ at distance r, which is a random function of r and depends only weakly on the coordinates in the plane perpendicular to the line of sight. If the characteristic value of the refraction angle $\Theta_{\rm ref}$ is much bigger than the characteristic value of the diffractive or scattering angle

 \begin{displaymath}% \Theta_{\rm ref} \gg \Theta_{\rm dif} \end{displaymath} (26)

then the frequency structure of the scintillation is formed by the frequency dependence of the displacement of the beam path (Shishov 1973). Using the equation for the two-frequency field coherence function derived by Shishov (1973) and taking into account inequality (26) we obtain

 \begin{displaymath}% \frac{\partial B_{\rm E}}{\partial r} + \frac{2}{r} \frac{\... ...la_{\eta} B_{\rm E} = -\frac{1}{2} D (r \vec{\eta}) B_{\rm E}, \end{displaymath} (27)

where we have used a spherical co-ordinate system, in which $\vec{\eta}$denotes a two-dimensional angle. The solution is

 \begin{displaymath}% \begin{array}{lcl} B_{\rm E} (\vec{\eta},\Delta f) & = & <I... ... D \left(r \vec{\eta} + \vec{\rho_{\rm ref}}\right) \end{array}\end{displaymath} (28)


 \begin{displaymath}% \vec {\rho_{\rm ref}} = 2 \frac{\Delta f}{f} \int\limits_0^R {\rm d}r'\frac{r}{r'} \vec{\Theta_{\rm ref}}. \end{displaymath} (29)

This solution describes the spatial pattern of diffractive scintillation for $\Delta f = 0$ and corresponds to Eq. (7). The fine frequency structure of scintillation is determined by the refractive relative displacement of beam paths at two frequencies that causes rearrangement of speckles of the diffractive scintillation pattern with a changing of frequency.

For the frequency correlation function of intensity fluctuations one can show that for small values of $\Delta f$ and in case of a power law turbulence spectrum the following relations hold

 \begin{displaymath}% B_{\rm I,ref} (\Delta f) = <I>^2 \left[1- D_{\rm S,1} \left(\rho_{\rm ref,eff} (\Delta f)\right)\right] \end{displaymath} (30)


 \begin{displaymath}% \rho_{\rm ref,eff} (\Delta f) \simeq \frac{\Delta f}{f} \frac{1}{3} R \Theta_{\rm ref,eff} \propto f^{-3}, \end{displaymath} (31)

where $\Theta_{\rm ref,eff}$ is the effective refraction angle

 \begin{displaymath}% \Theta_{\rm ref,eff} = \frac{\left[ \int\limits_0^R {\rm d}... ...r^R {\rm d}r'\frac{r}{r'}\right\}^{n-2}\right]^{1/(n-2)}}\cdot \end{displaymath} (32)

The numerical coefficient in Eq. (31) was calculated for the case of a statistically homogeneous medium. Using Eq. (31) one can estimate the phase structure function
 
$\displaystyle % D_{\rm S,1} (\rho_{\rm ref,eff} (\Delta f))$ $\textstyle \simeq$ $\displaystyle \displaystyle\frac{[B_{\rm I,ref} (0) - B_{\rm I,ref} (\Delta f)]}{<I>^2}$  
  $\textstyle \sim$ $\displaystyle \left( \displaystyle\frac{\Delta f}{\Delta f_{\rm ref}} \right)^{n-2}\cdot$ (33)

The spatial scale $\rho_{\rm ref,eff}$ is determined by Eq. (31). $\Delta f_{\rm ref}$ is the characteristic frequency scale of diffractive scintillation in the presence of strong angle refraction (Shishov 1973)

 \begin{displaymath}% \Delta f_{\rm ref} = \frac{B(n)c}{(R \Theta_{\rm ref,0} \Theta_{\rm dif,0})}, \end{displaymath} (34)

where c is the velocity of light, $\Theta_{\rm ref,0}$ the characteristic value of the refractive angle, and $\Theta_{\rm dif,0}$ the scattering angle. The numerical coefficient B(n)depends on the power law index n and on the distribution of the medium along the line of sight.

The weak dependence of $\rho_{\rm\Delta f,dif}$ on $\Delta f$ in comparison with that of $\rho_{\rm\Delta f,ref}$ leads to dominant diffraction effects for very small values of $\Delta f$

\begin{eqnarray*}\Delta f \ll \Delta f_0, \end{eqnarray*}


with

 \begin{displaymath}% \Delta f_0 \simeq f \left( \frac{3 \varphi_{\rm Fr}}{\Theta... ...quad \varphi_{\rm Fr} = \left( \frac{1}{kR} \right)^{1/2}\cdot \end{displaymath} (35)

One can use Eqs. (20) and (33) also if the scintillation is weak because these equations correspond to single or weak scattering on inhomogeneities with scales of the order of $\rho_{\rm dif} (\Delta f)$ or  $\rho_{\rm ref} (\Delta f)$.

It is important for practical observations that the diffractive and refractive model gives indeed for small values of $\Delta f$different functional dependences of the intensity correlation function on $\Delta f$. The power law indices of the temporal and frequency structure functions of intensity fluctuations are different: (n-2) and (n-2)/2 in case of the diffractive model. The indices are equal and have the value (n-2) for the refractive model.

The dynamic spectra of pulsars must show a strong frequency drift in case of a phase screen model, if refractive effects dominate over diffractive effects. The effect of the frequency drift is however weak as shown by Shishov (1973) in case of an extended random medium. It is possible that the frequency drift can in this case be the result of an interaction of refractive scintillation with strong angular refraction. This problem must be analyzed in more detail in future.

3 Observations

A comprehensive collection of observations of PSR B0329+54, covering the wide frequency range from 102 MHz to 5 GHz, was available for this paper to compare the theoretical results of Sect. 2 with measurements. The collection of data was used to construct the phase-structure function as outlined above.

3.1 Weak scintillation at 5 GHz

The emission of PSR B0329+54 is fortunately strong at high frequencies - our observations were taken at 4.7 and 10.6 GHz - showing clearly visible scintillation (Malofeev et al. 1996). We completed and extended the available observational material with specific measurements in July 1997 at the 100-m radio telescope of the Max-Planck-Institut für Radioastronomie, when we used a filterbank receiver with four channels of 60 MHz bandwidth each. The system temperature was in average about 60 K and the channels were centered at 5060, 4940, 4760 and 4640 MHz. Individual pulses were sampled at 1/1024 of the observing period and integrated subsequently in a data logger over ten periods to improve the signal-to-noise ratio. These blocks of data were used in the following analysis separately for each channel. More details on the observing method, the receiver and the pulsar backend are presented by Kijak et al. (1997).


  \begin{figure} \par\includegraphics[width=14cm,clip]{shishov_fig1.eps}\end{figure} Figure 1: Flux density variations over 400 min at four frequencies.
Open with DEXTER

Examples of intensity variations in the four channels are shown in Fig. 1. The measurements were smoothed over 20 data points by a running mean to improve the signal-to-noise ratio. The intensity fluctuations contain fast noisy oscillations and variations with a temporal scale order of the order of 10 min due to scintillation. The sharp peaks in the first 110 min are due to interference. They were excluded before we started the correlation analysis. It is interesting to note that the scintillation patterns at 5060 and 4640 MHz are partly decorrelated. To establish the degree of decorrelation, normalized cross-correlation functions were calculated between frequency channels. We introduced also an artificial time shift of 20 data points (142.6 s) between the channels to minimize the influence of noise. The dependence of the cross-correlation coefficient on frequency separation is shown in Fig. 2. The effect of decorrelation is obvious for frequency lag $\Delta f$ of about 8% of the observing frequency f.

At these high observing frequencies we are in the weak scintillation regime and the flux density variations should be well correlated in all four channels, which is obviously not the case. The decorrelation indicates that we see the influence of strong refraction, changing with frequency. We will discuss this later.


  \begin{figure} \par\includegraphics[width=6.5cm,clip]{shishov_fig2.eps}\end{figure} Figure 2: Cross-correlation coefficients of the flux variations in dependence on frequency separation corrected for noise.
Open with DEXTER

The temporal structure function of intensity fluctuations $D_{\rm I}(t)$ for time shifts t>200 s (calculated for f=4640 MHz) is shown in Fig. 3. It was corrected for noise by subtracting the mean value of the structure function calculated over the first 10 points.


  \begin{figure} \par\includegraphics[width=6.5cm,clip]{shishov_fig3.eps}\end{figure} Figure 3: Structure function of intensity variations normalized by the mean intensity in the time domain at a frequency of 4640 MHz corrected for noise.
Open with DEXTER

3.2 Diffractive scintillation at 610 MHz

At 610 MHz we used part of the data of a two year monitoring program of PSR B0329+54 at the NRAO 42-m telescope in 1994-1995. The observations were made using 1024 channels of the NRAO Spectral Processor covering a total bandwidth of 20 MHz in two orthogonal polarisations from which the total intensity was computed as the sum of the two signals. The spectra were accumulated during 59 s and written to magnetic tape for subsequent off-line analysis. The off-pulse spectra were subtracted from the on-pulse spectra. Each observation lasted for about 70 min.

We computed for all the 70-min dynamic spectra the normalized autocorrelation function (ACF) versus frequency and time to visualise the diffractive scintillation pattern, from which we calculated the characteristic frequency and time scales, $\Delta f_{\rm d}$ and  $\Delta t_{\rm d}$. $\Delta f_{\rm d}$ and $\Delta t_{\rm d}$ were defined as half of the ACF width at the level of 0.5 along the frequency axis and at the level of 1/e along the time axis after removing the spike at zero lag due to noise. Both, $\Delta f_{\rm d}$ and $\Delta t_{\rm d}$, are changing with time up to a factor of 3, an effect which is in accordance with the observations of Bhat Ramesh et al. (1999) at 327 MHz and which we will consider in a forthcoming paper.


  \begin{figure} \par\includegraphics[width=8.8cm,clip]{shishov_fig4.eps}\end{figure} Figure 4: Mean auto- (top) and cross-correlation (bottom) function between successive spectra in dependence on frequency lag. One frequency lag corresponds to 19.5 kHz. The data points are averaged over 70 min.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.8cm,clip]{shishov_fig5.eps}\end{figure} Figure 5: Mean auto- (top) and cross-correlation (bottom) function at f = 610 MHz in dependence on time lag. One time lag corresponds to 59 s. The data points are averaged over 1024 frequency channels.
Open with DEXTER

We computed in addition average autocorrelation functions for each day. The mean frequency ACF was calculated from the normalized ACFs averaged in time over all 59-s time bins whereas the mean time ACF was calculated by averaging over all normalised ACFs in each frequency channel. As an example the average time and frequency ACFs are shown for the observations of 12 March 1994 in the upper panels of Figs. 4 and 5. Figure 4 presents the frequency dependence, where one lag in frequency corresponds to 19.5 kHz, and Fig. 5 the time dependence with one lag representing 59 s in time (the spike at zero lag is due to noise). The corresponding crosscorrelation functions between intensity variations in neighbouring time bins ($\tau = 1$) and frequency channels ( $\Delta f = 1$) averaged over all pairs are shown at the bottom of these figures.

More informative for a comparison with theory are the structure functions (SF), since their slopes on a log-log scale contain information on the spectral distribution of the inhomogeneities in the ISM. We calculated therefore time and frequency structure functions separately for each channel and time bin and averaged them afterwards. The structure functions were normalized by the square of the mean intensity. The off-spectra, containing information on the noise, were treated in the same way. A special treatment is necessary for the first points of the SFs, $D_1 (D (\tau = 1),\ D (\Delta f = 1))$, since these points include noise and correlated signal. We subtracted from the SFs a value of $D_1 \times k_i$, where $k_1 = CCF~ (\tau = 0,\ \Delta f = 1) / ACF~ (\tau = 1)$ and $k_2 = CCF~ (\Delta f = 0,\ \tau = 1)/ACF ~ (\Delta f = 1)$ for a correction of the time and frequency structure functions. The ACFs and CCFs in the time and frequency domain were calculated as described earlier.


  \begin{figure} \par\includegraphics[width=8.8cm,clip]{shishov_fig7.eps}\end{figure} Figure 6: Mean time structure function of intensity variations at f = 610 MHz corrected for noise as described in the text and normalized by 2 m2, where m is the modulation index. The line is a fit to the first points.
Open with DEXTER


  \begin{figure} \par\includegraphics[width=8.8cm,clip]{shishov_fig6.eps}\end{figure} Figure 7: Mean frequency structure function of intensity variations at f = 610 MHz corrected for noise as described in the text and normalized by 2 m2, where m is the modulation index. The line is a fit to the first points.
Open with DEXTER

Figures 6 and 7 present the mean time and frequency structure functions at 610 MHz. The lines delineate fits to the first points of the SFs. To present the data of different frequencies in the same scale, the mean SFs were normalized by $2\times m^2$, where mis the modulation index computed on the basis of a two-dimensional frequency time array (70 min $\times$ 20 MHz).

3.3 Diffractive scintillation at 102 MHz

At 102 MHz we used observational material published by Popov & Soglasnov (1984); the observations were made with the Large Phased Array (BSA) at Pushchino in 1978. The analysis followed a different normalization procedure since the noise contribution in their correlation functions is large with the result that the correlation for small lags in time and frequency is significantly less then 1. We decided therefore to renormalize the time and frequency correlation functions $B_{\rm I}(t)$ and $B_{\rm I}(\Delta f)$ so that the correlation at zero lag becomes indeed 1: $B_{\rm I}(t=0) = B_{\rm I}(t=1) = 1$ and $B_{\rm I}(\Delta f=0) = B_{\rm I}(\Delta f=1) = 1$. These modified functions were used for a comparison with our data at cm and dm wavelengths as shown in Figs. 8 to 10.

Popov & Soglasnov (1984) mention two frequency scales for the intensity variations: 100 Hz and 750 Hz, the latter feature amounting to about 20% of the 100 Hz component in the correlation function. The 100 Hz structure corresponds quite well to observations at other frequencies and the estimated scintillation time scale, whereas the weaker 750 Hz feature might very well be due to statistical fluctuations. We consider therefore the 100 Hz fluctuations to be real and caused by ISS, whereas the 750 Hz feature is doubtful.

4 Compilation of the observational data at different frequencies


  \begin{figure} \par\includegraphics[width=6.5cm,clip]{shishov_fig8.eps}\end{figure} Figure 8: Time structure function of phase fluctuations at f0 = 1 GHz as compiled from the following data points 1.) 102 MHz observations: filled circles; 2.) 610 MHz: open circles; 3.) 5 GHz: stars. The solid line corresponds to the best fit to the first points of the structure functions. The slope is 1.50. The dashed line corresponds to n = 1.67 (Kolmogorov model).
Open with DEXTER

Using the theory developed in Sect. 2 one can reduce the measurements obtained at many different frequencies to one given frequency f0. Figure 8 shows the resulting time structure function $D_{\rm S}(t)$ at the frequency f0 = 1 GHz as compiled from Eqs. (10) and (13). Data points at different frequencies are marked by different symbols. The first points obtained at 5 GHz (stars) correspond very well to those at 610 MHz proving that we used the right correction procedure for measurements at different frequencies. The uncertainty increases, however, at large time lags (>2000 s) for the 5 GHz observations because the ratio between the time duration of the observations and the characteristic scintillation time scale becomes smaller and smaller. For a fit only the first three points were used. The measurements follow a power law:

 \begin{displaymath}% D_{\rm S} (t, f_0) = \left( \frac{t}{t_0} \right)^{\alpha}\cdot \end{displaymath} (36)

The best fit corresponds to the parameters: t0 = 1000 s and $\alpha = 1.50$ which is shown in Fig. 8 by the solid line. The accuracy of the values of $D_{\rm S}$ in Fig. 8 can be estimated to be about 15% of $D_{\rm S}$, which means that the error bars, are comparable in size to the dots in the figure. The most objective estimation of the total error may be obtained from the scatter of the $D_{\rm S}(t)$ values measured at different frequencies. The resulting accuracy in the determination of the exponent is $\alpha = 1.50 \pm 0.05$. The dashed line corresponds to: t0 = 800 s and the Kolmogorov value $\alpha = 1.67$, which describes the data points within $3 \sigma$.

For the computation of the frequency dependence, i.e. for the computation of $D_{\rm S}(\Delta f)$ at f0 = 1 GHz, one may use two different approaches. One approach is to base the analysis on the diffractive scintillation model, described by Eqs. (20), (22), and (23). The result of such a compilation is shown at the top of Fig. 9. One recognises immediately a strong discrepancy between the data sets at different frequencies.


  \begin{figure} \par\includegraphics[width=7.7cm,clip]{shishov_fig9.eps}\end{figure} Figure 9: Frequency structure function of phase fluctuations at f0 = 1 GHz based on data at 102 MHz, 610 MHz and 5 GHz (the different frequencies are characterized by different sizes of the circles) for two models: 1.) diffractive scintillation model; 2.) model with angular refraction. The solid line corresponds to the best fit to the first data points of the structure function. The slope is 1.47.
Open with DEXTER

Alternatively one may assume the strong refraction model (model 2) using Eqs. (23), (31), and (33). It should be remembered that $\Delta f$ is reduced to f0 in this case according to

 \begin{displaymath}% \Delta f (f_0) = \Delta f (f) \left( \frac{f_0}{f} \right)^3\cdot \end{displaymath} (37)

The results for the refractive model are shown at the bottom of Fig. 9 demonstrating now a good agreement between the data sets at different frequencies. The plotted power law corresponds to

 \begin{displaymath}% D_{\rm S,1} (\Delta f, f_0) = \left( \frac{\Delta f (f_0)}{\Delta f_1} \right)^{\alpha} \end{displaymath} (38)

with $\Delta f_1 = 3.6 \times 10^3$ kHz and $\alpha =1.47$. The accuracy of the values of $D_{\rm S}$ in Fig. 9 is about 15% so that the error bars a comparable in size to the shown points. No change in the local value of the slope can now be recognized down to the smallest value of $\Delta f$. We may conclude therefore that the refractive effect dominates over all ranges of the variation of $\Delta f$.

It is interesing to compare the angular refractive effect discussed above with the slanting features which are sometimes observed in the dynamical spectra of pulsars (see, for example, Rickett 1969; Gupta et al. 1994; Stinebring et al. 1996; Bhat et al. 1999). The angular refractive effect, which is evident in the frequency correlation function of intensity fluctuations, is stable over more than 15 years, because the time lag between the observations at 102 and 610 MHz as well as 5 GHz is about 15 years. The observations do, however, not show corresponding stable slanting features. This is an additional argument in favour of an extended model of the interstellar turbulent plasma near to the Sun because the frequency drift of diffraction patterns is strong for a phase screen and weak for an extended medium (Shishov 1973).

In case of a power law turbulence spectrum, the theory predicts small angle refraction of the order of more or less $\delta \Theta_{\rm ref} \approx m_{\rm ref} \Theta_{\rm dif,0}$, where $m_{\rm ref}$ is the scintillation index of refactive scintillation and $\Theta_{\rm dif,0}$ the scattering angle. The observations show however sometimes strong frequency drifts, that correspond to very large values of the refraction angle $\Theta_{\rm ref} \approx (5-10) \Theta_{\rm dif,0}$. It is difficult to explain such values of the refraction angle in terms of a power law turbulence spectrum only. It may be necessary to introduce into the model additional large scale inhomogeneities that give strong angular refraction, a modification which may open up new possibilities. The problem of refractive scintillation in presence of strong angle refraction is difficult and needs obviously special consideration.

One may reduce the temporal structure function to a spatial structure function $D_{\rm S}(\rho)$ using the pulsar velocity $V = 139~ {\rm km}~ {\rm s}^{-1}$ and the equation

 \begin{displaymath}% \rho = V t. \end{displaymath} (39)

An alternative would be to reduce the frequency structure function $D_{\rm S,1} (\Delta f)$ to the spatial structure function $D_{\rm S}(\rho)$ by use of Eqs. (23) and (31); one has then to specify the distance R and the refractive angle $\Theta_{\rm ref}$. Good agreement between $D_{\rm S}(t, f_0)$ and $D_{\rm S,1} (\Delta f (f_0),~ f_0)$ can be obtained if we suppose that R =1.4 kpc and the refractive angle is equal to

 \begin{displaymath}% \begin{array}{lcl} \Theta_{\rm ref} & \simeq & 3V_{\perp} \... ...\simeq 10 \varphi_{\rm Fr}, \quad f = 5~ {\rm GHz}. \end{array}\end{displaymath} (40)

Notice that the evaluation of $\theta_{{\rm ref}}$ in Eq. (40) is independent on the value of R because the ratio $V_{\perp}/R$ designates an angular velocity representing a measured quantity.

The whole compilation for $D_{\rm S}(\rho)$ is shown in Fig. 10, where the data points for $D_{\rm S}(t)$ are taken from Fig. 8 and converted to the spatial scale using Eq. (39).


  \begin{figure} \par\includegraphics[width=6.3cm,clip]{shishov_fig10.eps}\end{figure} Figure 10: Structure function of phase fluctuations at f0 = 1 GHzversus spatial scale of the inhomogeneities. The signs are the same as in Fig. 8; the cross designates the value computed from an analysis of refractive scintillation; the point at 1014 m is an upper limit from timing observations at 102 MHz. The points at the highest spatial scales are from DM variations of pulsar pairs in globular clusters. The line 1 corresponds to the best fit to the first points of the structure function, the line 2 to the Kolmogorov spectrum, and the line 3 is computed from Eqs. (41) and (44) as explained in the text.
Open with DEXTER

Further data points for the determination of $D_{\rm S}(\rho)$ can be obtained from measurements of refractive scintillation at 610 MHz as given by Smirnova et al. (1998) and from timing data at 102 MHz as published by Shabanova (1995). These points have been added to Fig. 10. The 610 MHz measurements allow one to calculate a value of $D_{\rm S}(\rho)$ at $\rho = T_{\rm ref} V = 2 \times 10^{11}$ m by use of Eqs. (16) and (17) with n=3.5 and assuming a scintillation index m = 0.37 and a characteristic time scale of $T_{\rm ref} = 17$ days. In addition, the variations of the pulse arrival times at 102 MHz with $\sigma_{\tau} = 0.7$ ms provide an upper limit for $D_{\rm S}$ at $\rho =10^{14}$ m if we use Eq. (6).

Data of dispersion measure variations for pulsars in globular clusters (Shishov & Smirnova 2002) were used to estimate $D_{\rm S}(\rho)$ for large values of $\rho$ in the range of one parsec (see Fig. 10). These data were reduced to our reference value of the dispersion measure, i.e. $DM_0 = 30~ {\rm pc~ cm}^{-3}$. The errors for these values are about 50% of $D_{\rm S}$. A detailed discussion of this data material can be found at Shishov & Smirnova (2002).

All the experimental data points were fitted by two lines where line 1 has a slope of $\alpha =1.5$ which corresponds to an exponent of n = 3.5 in the spectrum and where the line 2 has a slope of $\alpha = 1.67$ corresponding to a Kolmogorov spectrum. The point characterising the refractive scintillation falls above line 1 with $\alpha =1.5$ but we know that this point depends on the used theoretical model for the spectrum and the distribution of inhomogeneities along the line of sight to PSR 0329+54 and can be reduced.

We show in addition in Fig. 10 by line 3 the relation

 \begin{displaymath}% 2[\Delta S (\rho)]^2 = (k \Theta_{\rm ref} \rho)^2 \simeq 3.6 \times 10^{-15} (\rho [m])^2. \end{displaymath} (41)

This line corresponds to the case when the refraction dominates and $\Delta S (\rho)$ designates the phase difference measured at points separated by a distance $\rho$. $2[\Delta S (\rho)]^2$ can be considered then as a rough estimate of the phase structure function $D_{\rm S,ref}(\rho)$ influenced mainly by refractive inhomogeneities. The exact value of $D_{\rm S,ref}(\rho)$ can be described by Eq. (41), if $\Theta_{\rm ref}^2$ is substituted by  $<\Theta_{\rm ref}^2>$.

5 Discussion

The observations show that the interstellar plasma near to the Sun along the line of sight to PSR 0329+54 consists of two types of inhomogeneities: turbulent irregularities, which produce the scattering and scintillation effects, and much larger scale irregularities, which are responsible for angular refraction effects only.

The turbulence spectrum can be described by the Karman model (Tatarskii 1971)

 \begin{displaymath}% \Phi_{N_{\rm e}} (q) = K (n) <(\Delta_{N_{\rm e}})^2> \frac{L_0^3}{(1 + \vert{\vec q}\vert^2 L_0^2)^{n/2}} \end{displaymath} (42)

with n = 3.5 and

\begin{eqnarray*}K(n) = \frac{\Gamma (n/2)}{[\pi^{3/2} \Gamma ((n-3)/2)]} = 5 \times 10^{-2}. \end{eqnarray*}


L0 designates here the outer scale of the turbulence spectrum and $<(\Delta {N_{\rm e}})^2>$ the average of the squared electron density fluctuations. $\Gamma$ is the gamma function.

Using this model one obtains for small values $\rho \ll L_0$

\begin{eqnarray*}D_{\rm S,sph} (\rho) \simeq 2 M(n) <(\Delta S)^2> (\rho / L_0)^{(n-2)} \end{eqnarray*}


with

 \begin{displaymath}% M (n) = \displaystyle\frac{\Gamma ((n-2)/2)}{\left[2^{n-2} (n-1) \Gamma (n/2)\right]} \simeq 0.19, \quad n = 3.5. \end{displaymath} (43)

The observational data are described very well by Eq. (43), if n = 3.5. This value of the power law index indicates conclusively that the data correspond to the model of weak plasma turbulence (Iroshnikov 1963; Kraichnan 1965).

For large values of $\rho > L_0$ $D_{\rm S}(\rho)$ is given by

\begin{eqnarray*}D_{\rm S}(\rho) \simeq 2 <(\Delta S)^2> \end{eqnarray*}


with
 
                                         $\displaystyle <(\Delta S)^2> = G (n) (\lambda r_{\rm e})^2 <(\Delta N_{\rm e})^2> L_0 R$  
    $\displaystyle G(n) = 2 \pi^{3/2} \displaystyle\frac{\Gamma((n-2)/2)}{\Gamma((n-3)/2)} \simeq 3.7, \quad n = 3.5.$ (44)

Using these equations one can estimate the values of $<(\Delta {N_{\rm e}})^2>$ and L0. Unfortunately, there are no reliable data points for large values of $\rho$ available. One can therefore only estimate the limits of L0

 \begin{displaymath}% 2 \times 10^{11}~ {\rm m} < L_0 < 10^{17}~ {\rm m}. \end{displaymath} (45)

If we take the upper limit of L0 = 1017 m and $2<(\Delta S)^2>\; = 3.5 \times 10^{12}$ for f=1000 MHz, the value of the density fluctuations becomes a maximum and equal to

 \begin{displaymath}% \delta N_{\rm e} = [<(\Delta N_{\rm e})^2>]^{1/2} \simeq 5 \times 10^{-4} ~ {\rm cm}^{-3} \end{displaymath} (46)

and the level of turbulence will be

 \begin{displaymath}% \delta^2 = \frac{<(\Delta N_{\rm e})^2>}{N_{\rm e}^2} \simeq 3 \times 10^{-4}. \end{displaymath} (47)

The nonlinear decrement of absorption will be in this case

 \begin{displaymath}% \gamma \approx \delta^2 \omega_0 = \frac{\delta^2 \nu_{\rm s}}{L_0} \approx 3 \times 10^{-17}~ {\rm s}^{-1} \end{displaymath} (48)

and waves can propagate through a distance of about

 \begin{displaymath}% R_0 \approx \frac{\nu_{\rm s}}{\gamma} \approx 10~ {\rm kpc}, \end{displaymath} (49)

where $\nu_{\rm s}$ is the sound velocity. The value of R0 is important for the problem of a localization of the turbulence energy sources. If L0 is greater or equal to 0.3 pc then R0 must be greater or equal than 1 kpc and waves generated in the spiral arms of the Galaxy can be the source of the energy of the turbulence in the space between the spiral arms. For smaller values of L0, i.e. L0 < 0.3 pc, the energy sources of the turbulence must be distributed in the interarm space.

The type of irregularities which are responsible for the angular refraction are characterised by much smaller values of L0

 \begin{displaymath}% L_{\rm0,ref} \approx 3 \times 10^{13}~ {\rm m}. \end{displaymath} (50)

The value of $L_{\rm0,ref}$ cannot be less than $3\times 10^{13}$ m because the frequency correlation functions of scintillation obtained at different frequencies and at epochs stretching over about 10 years which we used still show good correspondence.

Using the value (50) in the Eq. (44) one obtains

 \begin{displaymath}% \delta N_{\rm e,ref} \approx 10^{-2}~ {\rm cm}^{-3}. \end{displaymath} (51)

Such irregularities can be from the radio physics viewpoint dense clouds or holes in the interstellar plasma. Dense plasma clouds are however - as plasma physics tells - unstable and will expand and generate MHD turbulence. We see, however, not yet such a turbulent high-frequency tail. The irregularities can be identified with clouds of neutral hydrogen which can be considered as holes of the electron density. Dense cold clouds can indeed be stable in a surrounding hot plasma. If we suppose that the local variation of the electron density in clouds is $N_{\rm e}$, then the filling factor of such clouds is of the order of

 \begin{displaymath}% F = \frac{(\delta N_{\rm e,ref})^2}{N_{\rm e}^2} \approx 0.1. \end{displaymath} (52)

The observational data may therefore be explained by the following scenario. The main part of the inter-arm space is filled by a hot plasma with weak plasma turbulence. 10% of the space shows, however, inhomogeneities in the sense that there exists an electron deficit on an characteristic spatial scale of the order of $3\times 10^{13}$ m. These volumes of low electron density can be identified with clouds of neutral hydrogen. It is possible that there exists a distribution of cloud sizes and that we detect only the small-size end of the distribution. We could show in addition that the three-dimensional spectrum of electron density fluctuations follows a power law with index $\alpha + 2 = 3.50 \pm 0.05$. Reliable data points are available for spatial scales between 106 and 108 m. It appears very important to investigate the density fluctuation spectrum also in the region of the outer scale of turbulence. To solve this problem it seems desirable to obtain multi-frequency scintillation observations of pulsars together with good timing measurements.

This scenario is not universal. There is evidence that different types of turbulent spectra exist in different regions of the Galaxy (see, for example, Smirnova et al. 1998; Lambert and Rickett 2000). Multifrequency scintillation observations of different pulsars should allow one to investigate the real change of interstellar plasma turbulence parameters in different regions of the Galaxy.

Acknowledgements
The authors thank the referee B. J. Rickett for useful comments, which improved the paper considerably. This work was supported by INTAS grant No. 00-00849, NSF grant No. AST 0098685, the Russian Foundation for Basic Research, project code 00-02-17850, and the Russian Federal Science and Technology Program in Astronomy. We thank the NRAO operated by Associated Universities for support with the 610 MHz observations and L. B. Potapova and G. Breuer for technical assistance.

References


Copyright ESO 2003