Free Access
Issue
A&A
Volume 531, July 2011
Article Number A126
Number of page(s) 8
Section Astronomical instrumentation
DOI https://doi.org/10.1051/0004-6361/201015890
Published online 28 June 2011

© ESO, 2011

1. Introduction

The intrinsic polarization of a synchrotron emitting source together with knowledge of propagation effects through intervening media provide critical diagnostics for magnetic field orientation and fluctuations in a wide range of astrophysical contexts. Faraday rotation is a physical phenomenon where the position angle of linearly polarized radiation propagating through a magneto-ionic medium is rotated as a function of frequency. As introduced in Brentjens & de Bruyn (2005) and Heald (2009), Faraday rotation measure synthesis is an important tool for analysing radio polarization data where multiple emitting regions are present along a single line of sight. Observations of extragalactic sources, which by necessity must be viewed through the Faraday rotating and emitting Galactic interstellar-medium (de Bruyn et al. 2006; Brown & Rudnick 2009; Schnitzeler et al. 2007, 2009), are an obvious example of this regime. Burn (1966) introduced the Faraday dispersion function F(φ), which describes the intrinsic polarized flux per unit Faraday depth φ (in rad m-2), and its relationship with the complex polarized emission P(λ2) as (1)where λ is the wavelength. Note that P can also be written as P = Q + iU, where Q and U represent the emission of Stokes Q and Stokes U, respectively.

To study multiple emitting and Faraday rotating regions along each line of sight, we need to reconstruct the Faraday dispersion function, which is, in general, a complex-valued function of the Faraday depth φ. From Eq. (1), we can invert the expression to yield: (2)However, the problem is that we can not observe the polarized emission at wavelengths where λ2 < 0. Even for the wavelength range λ2 > 0, it is impossible to observe all wavelengths or frequencies. Brentjens & de Bruyn (2005) propose a synthesis method by first introducing an observing window function M(λ2). The observed complex polarized emission can then be described as (3)In this paper, the tilde denotes the observed quantities.

If the observing window function is M(λ2) with m channels, the RM spread function (RMSF) is be defined by (4)where the parameter is the mean of the sampled values between  and  within the observation window M(λ2); i is the ith channel in the observation window, and K is a normalising constant of the window function M(λ2). In this paper, we assume as a simplification that all channels have uniform weights for the m channels in the observing window function.

In Brentjens & de Bruyn (2005), the reconstructed Faraday rotation measure synthesis can be written in discrete form as (5)where is the reconstructed Faraday dispersion function. From Eq. (5), we can see that the Faraday dispersion function can be reconstructed provided that the spectral coverage is sufficient.

However, the reconstructed results generally include some side lobes. Using the terminology of radio interferometry, the result of Brentjens & de Bruyn’ method is a dirty version of the Faraday dispersion function and is abbreviated as “the dirty curve”. It is the convolution of F(φ) and the RMSF, and a deconvolution step may be used to clean it up. By borrowing the cleaning procedure in the image deconvolution method of Högbom CLEAN (Högbom 1974). Heald (2009) proposes the RM-CLEAN method which deconvolves  with the RMSF to remove the sidelobe response.

Recently, Frick et al. (2010) proposed a wavelet-based Faraday RM synthesis method. In that approach, the authors assume specific magnetic field symmetries in order to project the observed polarization emissions onto λ2 < 0.

Compressive sensing/sampling (CS) (Candès & Wakin 2008; Candès 2006; Wakin 2008) has been one of the most active areas in signal and image processing over the last few years. Since CS was proposed, it has attracted very substantial interest, and has been applied in many research areas (Wakin et al. 2006; Lustig et al. 2007; Puy et al. 2010; Mishali et al. 2009; Bobin & Starck 2009). In radio astronomy, CS has attracted attention as a tool for image deconvolution. Wiaux et al. (2009) compare the CS-based deconvolution methods with the Högbom CLEAN method (Högbom 1974) on simulated uniform random sensing matrices with different coverage rates. They apply compressive sampling for deconvolution by assuming the target signal is sparse. Wiaux et al. (2009) proposed a new spread spectrum technique for radio interferometry by using the non-negligible and constant component of the antenna separation in the pointing direction. Recently, a new CS-based image deconvolution method was introduced in Li et al. (2011) in which an isotropic undecimated wavelet transform is adopted as a dictionary for sparse representation for sky images.

In this paper, we propose three new CS-based RM synthesis methods. In Sect. 2, the three CS-based RM synthesis methods are proposed. The implementation details of the general experiment layout is given in Sect. 3. Simulation results from the traditional methods are compared with those from CS-based methods in Sect. 4. The final conclusions are given in Sect. 5.

2. CS-based RM synthesis

CS is primarily a sampling theory for sparse signals. A sensing matrix (Candès et al. 2006b) is used to sample a signal with sparsity (few non-zero terms) or a sparse representation with respect to a given basis function dictionary. Given a limited number of measurements, generally less than the number of unknowns in the target signal, the target signal can be reconstructed by optimisation of an L1 norm. More information on the key concepts (such as sparsity, incoherence, the restricted isometry property, and the L1 norm reconstruction) and results can be found in Candès & Romberg (2007); Candès et al. (2006a); Candès & Wakin (2008); Candès (2006).

CS includes two steps: sensing/sampling and reconstruction. This is in contrast to Nyquist-Shannon theory which measures the target signal directly without the reconstruction step. In this paper, we will focus on the reconstruction step (calculating the Faraday dispersion function given an observing window) rather than the sensing step (the selection of the observing window), because the observing frequency range and the bandwidth for each channel are usually fixed for a given telescope array.

To proceed with the CS approach, we rewrite the Fourier relationship as a matrix equation. The projection of the Faraday dispersion function to the polarized emission can be described as a matrix Y of size m  ×  N(6)The inverse of the projection is the conjugate transpose of Y (7)where  ∗  denotes the conjugate transpose. Suppose f denotes the original Faraday dispersion function F(φ) in a vector format of length N, then the relationship between the Faraday dispersion function and the observed radio emission is: (8)where denotes the observed polarized emission in a vector format of length m.

Because we can only measure a limited number of observations with the limited number of channels, i.e. m ≪ N, there are many different potential Faraday dispersion functions consistent with the measurements. To resolve these ambiguities, the usual approach is to use some prior information to select a solution. The prior information can be: the Faraday dispersion function is real; the Faraday dispersion function has only point like signals which are sparse in the Faraday depth domain or the Faraday dispersion function has a sparse presentation with respect to a dictionary of basis functions, to name just a few. Our three synthesis methods are based upon the last two structural assumptions.

Before introducing our new RM synthesis methods, we need to review two technical terms: Faraday thin and Faraday thick. A source can be either Faraday thin if λ2 △ φ ≪ 1, or Faraday thick if λ2 △ φ ≫ 1, where  △ φ is the extent of the source along the axis of Faraday depth φ. Faraday thin sources can be well described by Dirac δ function of φ, while Faraday thick sources have extensive support on the Faraday depth axis (Brentjens & de Bruyn 2005). Note that the definition of Faraday thin or thick is wavelength dependent.

2.1. RM synthesis for Faraday thin sources

The relationship between the Faraday dispersion function and the observed polarized radio emission is a Fourier pair if λ2 = πu, where u is a wavelength related parameter. Since the space and Fourier domain are perfectly incoherent (Candès & Romberg 2007), we can apply CS for RM synthesis in a straightforward manner provided there are Faraday thin sources along the line of sight since the screen is necessarily sparse.

In this context, CS recommends solving for the Faraday dispersion function by minimising the L1 norm (summed absolute value) of the dispersion function as inimising the L1 norm optimises the sparsity of the reconstruction. There remains one further obstacle – the dispersion function is complex. We handle this by summing the L1 norm of the real and imaginary parts: (9)where Re(•) and Im(•) denote the real and the imaginary parts, respectively. By forming a real-valued vector of double length (comprising of the real part and the imaginary part) of the complex-valued vector, almost all L1 norm optimization solvers can be used for Eq. (9). This CS-based rotation measure synthesis for Faraday thin sources is abbreviated as CS-RM-Thin. This is similar in concept to RM-CLEAN because the assumption for RM-CLEAN is that the Faraday dispersion function comprising of spike like signals. However, results in Sect. 4 show that CS-RM-Thin can provides superior results to RM-CLEAN.

2.2. RM synthesis for Faraday thick sources

CS-RM-Thin can work effectively when the Faraday dispersion function includes Faraday thin sources along the line of sight. This limits its application for the case when there are some Faraday thick sources along the line of sight. However, CS can still reconstruct the Faraday dispersion function efficiently provided that we can find a suitable dictionary of basis functions that can decompose the extended sources into a sparse representation as described in Candès (2006). In this paper, we adopt the Daubechies D8 wavelet transforms (Daubechies 1992) as the dictionary. Other wavelet transforms can also be adopted, the selection depends on the property of the Faraday dispersion function. We choose the D8 wavelet transform, because we assume that the Faraday dispersion function with thick sources is a sinc-like signal.

We can rewrite Eq. (8) as (10)where W-1 is the inverse wavelet transform matrix of size N  ×  N; α is the wavelet coefficient of the Faraday dispersion function f. The wavelet transform matrix is denoted as W, therefore, α = Wf. Other symbols follow the definitions in Eq. (8). Under the condition that , we adopt the following assumption: both the real part and the imaginary part of the Faraday thick sources will have a sparse representation in the wavelet domain independently. Then the wavelet based CS RM synthesis method for Faraday thick sources can be written as (11)This CS-based rotation measure synthesis for Faraday thick sources is abbreviated as CS-RM-Thick.

2.3. RM synthesis for Faraday mixed sources

So far, we have proposed two RM synthesis methods: CS-RM-Thin and CS-RM-Thick for solving Faraday thin sources and thick sources, respectively. However, this begs the question: which method should be selected if there are both Faraday thin sources and thick sources along the line of sight? Moreover, how can we make a selection if we have no prior information about the Faraday dispersion function, i.e. we are not sure what it looks like? Clearly, neither of them is suitable, we therefore need another solution for solving the above problems. Let us assume that there are both Faraday thin sources and Faraday thick sources in F(φ). Suppose fthin denotes the Faraday thin sources in F(φ) in a vector format of length N; fthick denotes the Faraday thick sources in F(φ) in a vector format of length N, then fthin + fthick = f. Equation (8) can be rewritten as (12)Since we know that the L1 norm can preserve sparsity; Faraday thin sources show sparsity in the Faraday depth domain; Faraday thick sources show sparsity in the wavelet domain, we propose the following solution for the mixed circumstance (13)where the definition of W is the same as the above subsection. The above solution for Faraday mixed sources is still based on the spirit of CS by preserving the sparsity in the Faraday depth domain for Faraday thin sources and in the wavelet domain for Faraday thick sources, simultaneously. This CS-based rotation measure synthesis for Faraday mixed sources is abbreviated as CS-RM-Mix.

3. Implementation details

In this section, the implementation details of the above three proposed CS-based RM synthesis methods will be given.

3.1. Preparation

To create a general experiment layout, we borrow some definitions and conclusions from Brentjens & de Bruyn (2005). Both the diagram of the wavelength square λ2 and φ are displayed in Fig. 1. The maximum observable Faraday depth  ∥ φmax ∥  is given in Brentjens & de Bruyn (2005)(14)where δλ2 is the width of an observing channel. The full width at half maximum (FWHM) of the main peak of the RMSF can be estimated by (15)where  △ λ2 is the width of the total λ2 distribution.

Before using CS-based RM synthesis, the following two steps are needed:

  • 1.

    Select the resolution of the Faraday depth φR. Since we know the maximum observable Faraday depth φmax from Eq. (14) and the FWHM of the main peak of the rotation measure spread function from Eq. (15), we can select a grid resolution parameter φR in φ space, which should be four or five times less than δφ, to achieve Nyquist sampling. However, for some observational window functions M(λ2), the maximum scale (Faraday thickness) that one is sensitive to, estimated by , is actually smaller than δφ. In these cases, it might be practical to Nyquist sample this smaller scale in order to calculate φR. Based on φmax and  △ φ, we can calculate the number of grid points N as (16)

  • 2.

    Constructing the two matrices Y and Y.

thumbnail Fig. 1

This diagram shows the relationship between parameters in λ2 domain and φ domain, respectively.

The selection of these CS-based RM synthesis methods depends on the prior knowledge about the Faraday dispersion function. If we assume it includes Faraday thin sources only along the line of sight, we should select CS-RM-Thin. On the other hand, CS-RM-Thick should be used if we know that there are Faraday thick sources only. When we know that there are both Faraday thin sources and thick sources along the line of sight, CS-RM-Mix should be used. In most circumstances, we have no prior information about the Faraday dispersion function, CS-RM-Mix can always be used to reconstruct a reliable result as a compromise.

3.2. L1 norm solvers for CS-Based RM synthesis methods

For CS-RM-Thin and CS-RM-Thick, many optimization methods (Beck & Teboulle 2009; Becker et al. 2011; Boyd & Vandenberghe 2004) can be used to solve the L1 norm minimization problem in a straightforward manner. There are many solvers or toolboxes, for example, L1-Magic Matlab toolbox which can be download from http://www.acm.caltech.edu/l1magic/. In this paper, L1-Magic is adopted for solving Eqs. (9) and (11). Fast Iterative Shrinkage-Thresholding Algorithm (FISTA) (Beck & Teboulle 2009) can also be used for solving this problem if we rewrite Eqs. (9) or (11) in a Lagrangian form.

As far as CS-RM-Mix is concerned, the solvers or toolboxes introduced above can be used for solving Eq. (13) indirectly. Suppose αthick denotes the wavelet coefficients of the thick sources fthick in the Faraday dispersion function in a vector format, and W-1 is the inverse wavelet transform matrix, then we have (17)Bring the above equation into Eq. (12), we have (18)Furthermore, Eq. (18) can be rewritten as (19)where I denotes the identity matrix of size N × N, and O is the matrix of all zeros with the size of N × N. If we denotes Ymix =  [YYm × 2N, T =  and , almost all L1 norm minimization solvers can be used to solve Eq. (13) with: (20)To help readers who are unfamiliar with L1 norm minimization to use or implement our proposed CS-based RM synthesis methods, we have developed a simple algorithm based on the iterative soft-thresholding algorithm (ISTA) (Beck & Teboulle 2009) for CS-RM-Mix as an example. The algorithm is as follows:

  • 1.

    Initialization:

    • (a)

      choose parameters: the soft-threshold τ (this can be set by 1 for most circumstances), the stopping-threshold δ (this can be set by the noise level);

    • (b)

      set the number of iteration l = floor(τ/δ);

    • (c)

      ; fthick = 0.

  • 2.

    Within l iterations:

    • (a)

      Reconstructing the Faraday thin sources

      • i.

        calculate the residual;

      • ii.

        calculate the gradient d = Y ∗ r;

      • iii.

        update fthin = fthin + d;

      • iv.

        soft threshold Re(fthin) and Im(fthin), respectively. Set any values below τ to zero and update fthin.

    • (b)

      Reconstructing the Faraday thick sources

      • i.

        calculate the residual ;

      • ii.

        calculate the gradient d = Y ∗ r;

      • iii.

        update fthick = fthick + d;

      • iv.

        calculate the wavelet coefficients for both the real part and imagery part of fthick, i.e. W·Re(fthick) and W·Im(fthick);

      • v.

        soft threshold the wavelet coefficients of W·Re(fthick) and W·Im(fthick). Set any values below τ to zero and update W·Re(fthick) and W·Im(fthick);

      • vi.

        calculate the inverse wavelet transform for both the real part and imaginary part, respectively, then update fthick.

    • (c)

      τ = τ − δ.

  • 3.

    Reconstructed Faraday dispersion function .

Note the Faraday thin sources and thick sources are reconstructed separately. This can be helpful when astronomers focus on either Faraday thin sources or thick sources only.

If we ignore the step 2b and set fthick = 0, the above algorithm will degenerate to CS-RM-Thin. On the contrary, if we ignore the step 2a and set fthin = 0, the algorithm will degenerate to CS-RM-Thick. In this paper, these CS-based rotation measure synthesis methods are implemented in MATLAB. Our code may be found at http://code.google.com/p/csra/downloads1.

4. Experimental results

We have adopted the standard test platform in Brentjens & de Bruyn (2005). In this platform there 126 observing channels within Window 1 (0.036 to 0.5 m) evenly distributed in λ2. Three different Faraday dispersion functions are simulated to test these CS-Based RM synthesis methods for Faraday thin sources, Faraday thick sources and mixed sources cases.

thumbnail Fig. 2

We have tested our methods on a Faraday dispersion function with four Faraday thin sources. From left to right in the first row are: a) original F(φ); b) dirty curve; c) RM-CLEAN. From left to right in the second row are: d) CS-RM-Thin; e) CS-RM-Thick; f) CS-RM-Mix. The thin solid line shows the real value, the dashed line the imaginary part, and the thick solid line the amplitude. All horizontal axis units are rad m-2, i.e. φ, and all vertical axis units are Jy m2 rad-1.

4.1. Simulation results for Faraday thin sources

We simulate a Faraday dispersion function containing four Faraday thin sources. See Fig. 2 for the function. From left to right, these sources are: F(−10) = 10−4i Jy m2 rad-1, F(−17) =  −7 + 5i Jy m2 rad-1, F(40) = 9−7i Jy m2 rad-1 and F(88) =  −4 + 3i Jy m2 rad-1. The true Faraday dispersion function is shown in the top left corner of Fig. 2. The thin solid line shows the real value, the dashed line the imaginary part, and the thick solid line the amplitude. The Faraday dispersion function in this test is complex valued, i.e. the intrinsic polarization angles are non-zero. We have selected this simulated dispersion function rather than the standard test (real valued Faraday dispersion function) in Brentjens & de Bruyn (2005) to investigate the behaviour when the intrinsic polarization angles are non-zero. From Eq. (15), we can calculate that the FWHM of the RMSF is around 14 rad m2. The dirty curve which is calculated by assuming all unmeasured emission to be zeros, is shown in Fig. 2b. The result of RM-CLEAN is shown in Fig. 2c. Note that RM-CLEAN cannot correctly reconstruct the magnitude or phase of the Faraday dispersion function in this case. It has been demonstrated previously that RM-CLEAN has difficulty with the separation of sources below the FWHM (Farnsworth et al. 2011). In Fig. 2a, from left to right, the distance between the first and the second sources is smaller than the FWHM. In Fig. 2c, the first source and second source are merged to form an unrealistic source. Results of CS-RM-Thin, CS-RM-Thick and CS-RM-Mix are shown in the Figs. 2d − f, respectively. As one might expect, the result of CS-RM-Thick is poor, because the Faraday dispersion function has no Faraday thick sources. In Fig. 2, though CS-RM-Mix does a better job than CS-RM-Thick and RM-CLEAN in terms of magnitude and phase of the reconstructed disperse functions, the result is far from good enough. The result of CS-RM-Thin is shown in Fig. 2d. We can see that CS-RM-Thin gives the best result, reconstructing the Faraday dispersion function without any error. This is consistent with CS theory which says that we can reconstruct the sparse signal exactly with “overwhelming probability” (Candès & Wakin 2008; Candès 2006; Wakin 2008). In this test, we select φR = 3.6 rad m2 which is around one quarter of FWHM, so N = 480.

To carry out a numerical comparison, we use the root mean square (rms) error to characterise the difference between the reconstructed  and the Faraday dispersion function f: (21)

The rms error is calculated for all the candidate methods, and the results are listed in Table 1. Results for this test can be found in the first row of the table. CS-RM-Thin gives the best result.

thumbnail Fig. 3

Reconstructed results of a Faraday dispersion function with two Faraday thick sources. From left to right in the first row are: a) original F(φ); b) dirty curve; c) RM-CLEAN. From left to right in the second row are: d) CS-RM-Thin; e) CS-RM-Thick; f) CS-RM-Mix.

Table 1

Numerical comparison results (rms error).

4.2. Simulation results for Faraday thick sources

We now test our CS-based methods for a Faraday dispersion function with Faraday thick sources. Here we assume that the Faraday dispersion function includes two sources F(φ) =  2−2i Jy m2 rad-1 where  −120 ≤ φ ≤ 40 and F(φ) =   −6−3i Jy m2 rad-1 where 30 ≤ φ ≤ 70. The simulated Faraday dispersion function is shown in Fig. 3a. We adopt the previous observing window for this test and the dirty Faraday dispersion function is shown in Fig. 3b. Note that this only provides us with the approximate shape of the Faraday dispersion function. As mentioned in Frick et al. (2010), the magnitude of F(φ) indicates the polarized emission of the region with Faraday depth φ and its phase defines the intrinsic position angle. For the study of polarized emission of galaxies, the magnitude of F(φ) is very important, and for the study of orientation of the magnetic field perpendicular to the line of sight, the phase information of F(φ) is crucial. Unfortunately, Brentjens & de Bruyn’ method can not reconstruct reliable phase information for F(φ). The cleaned version is shown in Fig. 3c. RM-CLEAN also failed to reconstruct the phase information. CS-RM-Thin does not work well for this case, because F(φ) is not sparse. The result of CS-RM-Thin is shown in Fig. 3d. CS-RM-Mix is also used for this test by assuming that there are both Faraday thin sources and thick sources. The result of CS-RM-Mix is shown in Fig. 3f. The reconstructed result from CS-RM-Thick is shown in Fig. 3e. Even though CS-RM-Mix gives a much better result than CS-RM-Thin and RM-CLEAN, the result is not as good as that of CS-RM-Thick. CS-RM-Thick provides the best approximation to the original Faraday dispersion function in terms of both magnitude and phase. This is also supported by the numerical comparison in Table 1. From the second row of the table, we can see that CS-RM-Thick gives the smallest rms error 0.72 which is slightly better than that of CS-RM-Mix 0.77.

4.3. Simulation results for Faraday mixed sources

So far we have tested our methods for both Faraday thin sources and Faraday thick sources, then we will test our CS-based methods for the mixed circumstances – both Faraday thin sources and thick sources. The simulated Faraday dispersion function for this test is shown in the top left corner of Fig. 4. There are three sources along the line of sight. From left to right, these sources are: F(−58) =  −4 + 3i Jy m2 rad-1, F(−30) = 10−4i Jy m2 rad-1, F(φ) = 2−6i Jy m2 rad-1 where 41 ≤ φ ≤ 100. The original Faraday dispersion function is shown in the top left corner of Fig. 2. We use the previous observing window for this test. Based on the observing window (with 126 observing channels distributed between 0.036 m to 0.5 m), the RMSF is calculated and shown in Fig. 4b. the dirty curve is shown in Fig. 4c. The cleaned version by RM-CLEAN is shown in Fig. 4d. RM-CLEAN performs badly in this test, because there is a Faraday thick source. In general, RM-CLEAN can only work well when there are Faraday thin sources along the line of sight. The results of CS-RM-Thin and CS-RM-Thick are shown in Figs. 4e and f, respectively. We can see that CS-RM-Thin reconstructs the two Faraday thin sources nicely, but failed to reconstruct the Faraday thick source. On the contrary, CS-RM-Thick can properly reconstruct the Faraday thick source, but expends the two Faraday thin sources by mistake. As introduced above, CS-RM-Mix can separate the Faraday thin components fthin and the Faraday thick components fthick during the reconstruction. In this test, the soft-threshold τ = 1 and δ = 0.001 in the proposed algorithm for CS-RM-Mix. The results of CS-RM-Mix are shown in the last row of Fig. 4. The separated Faraday components fthin and fthick are shown in Figs. 4g and h, respectively. The sum of the separated components is the reconstructed Faraday dispersion function which is shown in Fig. 4i. We can see that the result of CS-RM-Mix takes advantage of the results of both CS-RM-Thin and CS-RM-Thick, and it gives the closest approximation to the original F(φ). From objective evaluation point of view, CS-RM-Mix gives the minimum rms error 0.80 which can be seen from the third row of Table 1.

thumbnail Fig. 4

Reconstructed results of a Faraday dispersion function with two Faraday thin sources and a thick source. From left to right in the first row are: a) original F(φ); b) RMSF of the observing window with 126 observing channels distributed between 0.036 m to 0.5 m; c) dirty curve. From left to right in the second row are: d) RM-CLEAN; e) CS-RM-Thin; f) CS-RM-Thick. From left to right in the third row are: g) thin components fthin by using CS-RM-Mix; h) thick components fthick by using CS-RM-Mix; i) CS-RM-Mix i.e. fthin + fthick. All horizontal axis units are rad m-2, and all vertical axis units are Jy m2 rad-1.

4.4. Discussion

In the figures above, we show the reconstructed dispersion function without smoothing and addition of the residuals – a step commonly know as restoring. For real applications, restoration is an option if the robustness is insufficient.

From the above three tests, we can see that there is no single CS-based RM synthesis method with an outstanding performance for any circumstances. The best reconstruction can only be achieved when we have some prior knowledge about the Faraday dispersion function and select the relevant CS-based RM synthesis method. If such information is not available, CS-RM-Mix can always be used as a compromise. Another option is that we can either select CS-RM-Thin with a large φR or CS-RM-Thick with a small φR. If a large φR is selected, f is likely to be a sparse vector, so CS-RM-Thin should be selected for the reconstruction. On the other hand, a small φR can expand compact sources into extended sources, so CS-RM-Thick will become suitable for the reconstruction. It does not mean that CS-RM-Thick can be used for any cases with a small φR, because a smaller φR (which means a larger N from Eq. (16)) brings more unknowns in f and more uncertainty. We have to balance these factors.

For RM synthesis, the observing window is quite similar to the frequency filters. For example, the previously introduced observing Window 1 is like a low pass filter in the wavelength squared domain. Therefore, we should bear in mind, to observe radio sources with Faraday thick sources under the same restriction of m and δλ2, the higher frequency observing band the better.

These CS-based RM synthesis methods are not limited to the optimisation methods L1-Magic solver, FISTA and ISTA. Other L1 norm optimization solvers can also be adopted for solving Eqs. (9), (11) and (13). Though the wavelet transform is used as the sparse representations dictionary in this paper, there are some other potential basis functions can be used to achieve sparsity for the Faraday thick sources.

In summary, the performance of CS-based RM synthesis methods depend on the observing window, the resolution of φ, the number of measurements, and the sparsity of the Faraday dispersion function. The reconstruction of these CS-based RM synthesis methods take less time than RM-CLEAN in general. For example, CS-RM-Thin takes 3 s for the above tests in contrast with 5 s of RM-CLEAN. The calculation time really depends on the construction of the matrix Y, the larger N and m (N = 480, m = 126 for the above tests), the more time it takes. The computer is a 2.53-GHz Core 2 Duo MacBook Pro with 4GB RAM.

5. Conclusions

Faraday rotation measure synthesis is a very useful tool to study astrophysical magnetic fields. The problem in RM synthesis is to reconstruct the Faraday dispersion function given incomplete observations. From CS, we know that a signal with sparsity can be well reconstructed based on few measurements. We propose three CS-based RM synthesis methods by finding sparse representations of the Faraday dispersion functions F(φ) for different circumstances. CS-RM-Thin, CS-RM-Thick and CS-RM-Mix can be used for Faraday thin sources only, thick sources only and mixed sources, respectively. In general, Faraday thin sources show sparsity in the Faraday depth domain φ, therefore, we apply the CS reconstruction methods (L1 norm optimization solvers) in a straightforward manner i.e. CS-RM-Thin. Although Faraday thick sources are not sparse in the Faraday depth domain, they are sparse in the wavelet domain for a suitably chosen basis wavelet. Therefore, we apply the L1 norm optimization solvers in the wavelet domain i.e. CS-RM-Thick. When there are Faraday mixed sources along the line of sight, we preserve the sparsity by using L1 norm in the Faraday depth domain and the wavelet domain simultaneously i.e. CS-RM-Mix.

As shown in the experimental results, the performance of these CS-based methods is markedly superior to the traditional RM synthesis methods (Brentjens & de Bruyn 2005; Heald 2009) in terms of magnitude and angle of the reconstructed Faraday dispersion function. Exemplified by Fig. 2, both Brentjens & de Bruyn’ method and RM-CLEAN do not work well in disentangling two closely spaced sources. In contrast, CS-RM-Thin can separate the sources.


1

Download the file “CS_RM.zip” which includes both CS-RM-Thin and CS-RM-Thick algorithms.

Acknowledgments

We thank Jean-Luc Starck for early discussions on Compressive Sampling. In addition, comments, suggestions and derivations made by Dr. Brentjens during the review process are very much appreciated.

References

  1. Beck, A., & Teboulle, M. 2009, SIAM J. Imaging Sci., 2, 183 [Google Scholar]
  2. Becker, S., Bobin, J., & Candès, E. 2011, SIAM J. Imaging Sci., 4, 1 [Google Scholar]
  3. Bobin, J., & Starck, J.-L. 2009, Wavelets XIII, 7446, 74460I [Google Scholar]
  4. Boyd, S., & Vandenberghe, L. 2004, Convex optimization (Cambridge University Press), 716 [Google Scholar]
  5. Brentjens, M., & de Bruyn, A. 2005, A&A, 441, 1217 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Brown, S., & Rudnick, L. 2009, AJ, 137, 3158 [NASA ADS] [CrossRef] [Google Scholar]
  7. Burn, B. 1966, MNRAS, 133, 67 [NASA ADS] [CrossRef] [Google Scholar]
  8. Candès, E. 2006, Int. Congress of Mathematics, Madrid, Spain, 3, 1433 [Google Scholar]
  9. Candès, E., & Romberg, J. 2007, Inverse Problems, 23, 969 [NASA ADS] [CrossRef] [Google Scholar]
  10. Candès, E., & Wakin, M. 2008, IEEE Signal Proc. Mag., 25, 21 [Google Scholar]
  11. Candès, E., Romberg, J., & Tao, T. 2006a, IEEE Trans. inform. Theory, 52, 489 [CrossRef] [Google Scholar]
  12. Candès, E., Romberg, J., & Tao, T. 2006b, Comm. Pure Appl. Math., 59, 1207 [CrossRef] [Google Scholar]
  13. Daubechies, I. 1992, Ten lectures on wavelets, Society for Industrial and Applied Mathematics Philadelphia, PA, USA [Google Scholar]
  14. de Bruyn, A., Katgert, P., Haverkorn, M., & Schnitzeler, D. 2006, Astron. Nachr., 327, 487 [Google Scholar]
  15. Farnsworth, D., Rudnick, L., & Brown, S. 2011, ApJ, 141, 191 [Google Scholar]
  16. Frick, P., Sokoloff, D., Stepanov, R., & Beck, R. 2010, MNRAS, 401, L24 [NASA ADS] [Google Scholar]
  17. Heald, G. 2009, Cosmic Magnetic Fields: From Planets, 259, 591 [Google Scholar]
  18. Högbom, J. 1974, A&AS, 15, 417 [Google Scholar]
  19. Li, F., Cornwell, T., & de Hoog, F. 2011, A&A, 528, A31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Lustig, M., Donoho, D., & Pauly, J. 2007, Magnetic Resonance in Medicine, 58, 1182 [Google Scholar]
  21. Mishali, M., Eldar, Y. C., Dounaevsky, O., & Shoshan, E. 2009, CCIT Report #751, Oct.-09, EE Pub No. 1708, EE Dept., Technion – Israel Institute of Technology [arXiv:0912.2495] [Google Scholar]
  22. Puy, G., Wiaux, Y., Gruetter, R., Thiran, J., & De, D. V. 2010, IEEE International Symp. on Biomedical Imaging: From Nano to macro [Google Scholar]
  23. Schnitzeler, D., Katgert, P., & de Bruyn, A. 2007, A&A, 471, L21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Schnitzeler, D., Katgert, P., & de Bruyn, A. 2009, A&A, 494, 611 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Wakin, M. 2008, IEEE Signal Processing Magazine, 1623 [Google Scholar]
  26. Wakin, M., Laska, J., Duarte, M., & Baron, D. 2006, International Conference on Image Processing, ICIP 2006, Atlanta George, Oct. 2006, 1437 [Google Scholar]
  27. Wiaux, Y., Jacques, L., Puy, G., Scaife, A. M. M., & Vandergheynst, P. 2009, MNRAS, 395, 1733 [NASA ADS] [CrossRef] [Google Scholar]

All Tables

Table 1

Numerical comparison results (rms error).

All Figures

thumbnail Fig. 1

This diagram shows the relationship between parameters in λ2 domain and φ domain, respectively.

In the text
thumbnail Fig. 2

We have tested our methods on a Faraday dispersion function with four Faraday thin sources. From left to right in the first row are: a) original F(φ); b) dirty curve; c) RM-CLEAN. From left to right in the second row are: d) CS-RM-Thin; e) CS-RM-Thick; f) CS-RM-Mix. The thin solid line shows the real value, the dashed line the imaginary part, and the thick solid line the amplitude. All horizontal axis units are rad m-2, i.e. φ, and all vertical axis units are Jy m2 rad-1.

In the text
thumbnail Fig. 3

Reconstructed results of a Faraday dispersion function with two Faraday thick sources. From left to right in the first row are: a) original F(φ); b) dirty curve; c) RM-CLEAN. From left to right in the second row are: d) CS-RM-Thin; e) CS-RM-Thick; f) CS-RM-Mix.

In the text
thumbnail Fig. 4

Reconstructed results of a Faraday dispersion function with two Faraday thin sources and a thick source. From left to right in the first row are: a) original F(φ); b) RMSF of the observing window with 126 observing channels distributed between 0.036 m to 0.5 m; c) dirty curve. From left to right in the second row are: d) RM-CLEAN; e) CS-RM-Thin; f) CS-RM-Thick. From left to right in the third row are: g) thin components fthin by using CS-RM-Mix; h) thick components fthick by using CS-RM-Mix; i) CS-RM-Mix i.e. fthin + fthick. All horizontal axis units are rad m-2, and all vertical axis units are Jy m2 rad-1.

In the text

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

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

Initial download of the metrics may take a while.