A CHANDRA STUDY OF TEMPERATURE SUBSTRUCTURES IN INTERMEDIATE-REDSHIFT GALAXY CLUSTERS

In the spectral analysis that follows, we utilized the Chandra blank-sky templates for the S3 CCD as the background. The templates were tailored to match the actual pointings, and the background spectra were extracted and processed identically to the source spectra (Section 3.2.2). After rescaling each background spectrum by normalizing its high-energy end to the corresponding observed spectrum, we further enhance the accuracy of the background by adding an unabsorbed APEC model with fixed temperature (0.2 keV) and metal abundance (1 Z_☉; Snowden et al. 1998) to account for variations in the Galactic X-ray background (GXB; Vikhlinin et al. 2005). We note that the GXB component is significant only in A2204, where it has a flux of 2.0 ×10⁻¹⁴ erg s⁻¹cm⁻²arcmin⁻² in 0.3–1.0 keV and is nearly uniform over the S3 chip. This was also reported by Sanders et al. (2009), who indicated that it is probably due to the excess Galactic foreground.

Here we calculate the azimuthally averaged temperature profiles, which will be used to create the reference temperature maps applied in wavelet detection (Section 3.3.2) and calculation of excess energies in the temperature substructures (Section 4). We divided the inner 400 h⁻¹₇₀ kpc of each cluster into concentric annuli, masked all the detected point sources (≃15 point sources per cluster), and extracted the spectrum from each annulus. The corresponding spectral redistribution matrix files (RMFs) and auxiliary response files (ARFs) were generated using the CIAO tools mkwarf and mkacisrmf. The spectra were grouped with >20 counts per bin to allow the use of χ² statistics. We used XSPEC 12.4.0 package to fit the spectra, by modeling the gas emission with an absorbed APEC component and employing the PROJCT model to correct the projection effect. The absorption column density N_H was fixed to the Galactic value (Dickey & Lockman 1990) in all cases, except for A478, for which N_H was set to be about twice the Galactic value (Sanderson et al. 2005; Table 3); since the absorption excess in A478 exhibits a spatial extension beyond the cluster's core region, Pointecouteau et al. (2004) argued that it is likely to be associated with a local absorption substructure in our Galaxy. In order to calculate the error ranges of the model parameters, we executed the XSPEC command steppar and carried out multiple iterations at every step, to ensure that the actual minimum χ² is found. The best-fit gas temperature and abundance profiles are illustrated in Figure 3, in which multiple sets of annuli are used for each cluster to cross-check the results. We note that the inclusion of an additional gas component could only slightly improve the fits for the central regions of A478, A1650, A2204, A2244, and A3112, which is insignificant as indicated by the F-test.

Zoom In Zoom Out Reset image size

A significant temperature decline by a factor of >2 is observed in the central regions of A478, A1068, A2204, A2556, and A3112, four of which (except A2556) possess radio halos and/or lobes in the center (Section 2). In the remaining four clusters without apparent radio emission, the inward temperature drop is somewhat ambiguous, as the difference between the core temperature and mean cluster temperature, which refer to the temperatures measured in the central 100 h⁻¹₇₀ kpc (∼0.1r₅₀₀) and 100–400 h⁻¹₇₀ kpc regions (Figure 3), respectively, is found within 3σ error ranges. If we apply the classification criteria of Sanderson et al. (2006), the latter four clusters are non-cooling core clusters, whose central regions have not yet fully condensed. These results reflect the trend for the cooling core clusters to host central radio sources (Burns 1990). As shown in Figure 3, we also adopt Equations (4)–(6) of Vikhlinin et al. (2006) to describe the obtained radial temperature profiles in a smooth form, which will be used in Section 3.3.2 to calculate the reference temperature maps in order to examine the significances of detected temperature substructures, and in Section 4 to calculate the excess energies in these substructures.

In A478, A1068, A1650, A2204, A2556, and A3112, the abundance increases inward significantly (68% confidence level) from ≈0.2 Z_☉ in the outermost annuli to >0.6 Z_☉ in the innermost regions (≲ 30 h⁻¹₇₀ kpc). Emission measure-weighted gas temperature and abundance of each cluster are also calculated (Table 3), which are in good agreement with previous results obtained with Chandra and XMM-Newton (Sun et al. 2003; Takahashi & Yamashita 2003; Takizawa et al. 2003; Pointecouteau et al. 2004; Wise et al. 2004; Donahue et al. 2005; Sanders et al. 2005; Sanderson et al. 2005).

Following the methods of, e.g., O'Sullivan et al. (2005), Maughan et al. (2006), and Kanov et al. (2006), we define a set of knots r_i (>5000 knots per cluster) in the central 300 h⁻¹₇₀ kpc region, which are randomly distributed with a separation of Δ_ij < 5'' between any two adjacent knots i and j. To each knot, we assign an adaptively sized circular cell, which is centered on the knot and contains more than 800 counts in 0.7–8 keV after the background is subtracted. The typical cell radius ranges from ≲5'' at the center to about 10'' at r ⩾150 h⁻¹₇₀ kpc, which is always larger than the associated Δ_ij. We calculate the projected temperature of each cell T_c(r_i), which is assigned to the corresponding knot, by fitting the spectrum extracted from the cell with an absorbed APEC model, using the absorption column given in Table 3; adding another thermal component can neither improve the fit significantly, nor make the parameters constrained more tightly. In the fittings, the source spectrum is grouped to have >20 photons in each energy bin. In addition to the standard χ² statistics, we also have applied the maximum-likelihood statistics, and find that the best-fit parameters obtained with the two statistics agree nicely with each other within the 1σ error ranges. The calculated 1σ errors of gas temperature are typically <10% in the central 100 h⁻¹₇₀ kpc, which increase progressively to ≲15% (A478, A1650, A2244, and A3112) or ≲20% (the remaining five clusters) in 100–300 h⁻¹₇₀ kpc.

For any given position r, we define a scale s(r), so that there are at least 800 counts contained in a circular region centered at r, whose radius is s(r). We calculate the projected temperature at r by integrating over all the knots r_i located in the circular region,

$$\begin{eqnarray} T({\bf r})& = &\sum _{{\bf r}_i} (G_{{\bf r}_i}(R_{{\bf r},{{\bf r}_i}})T_{\rm c}({{\bf r}_i}))\big/\sum _{{\bf r}_i} G_{{\bf r}_i}(R_{{\bf r},{{\bf r}_i}}),\nonumber\\ &&\quad {\rm when} \quad \it R_{{\bf r},{{\bf r}_i}} < s({\bf r}), \end{eqnarray} \tag{ 1 }$$

where $R_{{\bf r},{{\bf r}_i}}$ is the distance from r to r_i, and $G_{{\bf r}_i}$ is the Gaussian kernel whose scale parameter σ is fixed at s(r_i). Since s(r) is essentially proportional to the square root of the local counts, the obtained temperature map T(r) (Figure 4) is less affected by the statistical uncertainties caused by surface brightness fluctuations. In addition, angular resolutions of ∼10'' are guaranteed by the use of compact Gaussian kernel. Maps for the upper and lower limits of T(r) are calculated in a similar way.

Zoom In Zoom Out Reset image size

Wavelet Transform. A visual inspection of the obtained temperature maps (Figure 4) indicates that there exist significant substructures, which typically have linear sizes of ∼100 h⁻¹₇₀ kpc, in all of the sample clusters. In order to describe these substructures in a quantitative way, we analyze the statistical properties of the temperature maps by applying the wavelet transform, a powerful tool in imaging process that identifies and disentangles small embedded features as a function of scales (e.g., Vikhlinin et al. 1997; Starck & Pierre 1998). In an analogous study, Bourdin & Mazzotta (2008) conducted the wavelet transform to achieve an optimized spatial binning of temperature maps. We employ a nonorthogonal version of discrete wavelet transform, the à Trous wavelet (Shensa 1992; Starck et al. 1995; Vikhlinin et al. 1998), which defines a straightforward decomposition rule. To be specific, at position r = (x, y) and scale level i (i = 1, 2, ...), the à Trous wavelet w_i(x, y) equals the signal difference between two smoothed maps c_{i − 1}(x, y) − c_i(x, y), where c₀(x, y) is the obtained temperature map T(x, y), and c_i(x, y) connects with c_{i − 1}(x, y) via c_i(x, y) = ∑_m,nh(m, n)c_{i − 1}(x + 2^{i − 1}m, y + 2^{i − 1}n), in which h(m, n) is the convolution kernel and (m, n) is its grid. Here h(m, n) is chosen to be a compact Gaussian $h(m,n) = e^{-(m^2+n^2)/2/2\pi }$. Each temperature map is thus decomposed into a set of subimages w_i(x, y) on nine scales (2ⁱ pixels, i = 1, 2, ..9), respectively.

Based on the obtained w_i(x, y) maps, we calculate the significant wavelet coefficient W_i(x, y) to detect the temperature substructures on different scales. For each cluster, we create a set of random knots r_i as used in Section 3.3.1, and assign a randomized temperature T'_c(r_i) to each knot, which is constrained by the measured temperature profile T(r) and its error ranges (Figure 3). By applying Equation (1) to the T'_c(r_i) map, we obtain a reference temperature map T_ref(r), in which the effects of the radial temperature gradients and knot tessellation are both involved. A series of such reference temperature maps are created by repeating this random process, and each reference map is decomposed into subimages w_i,ref(x, y) on nine scales by performing the same à Trous wavelet transform as described above. On each scale i, we determine the significant coefficient as W_i(x, y) ≡ w_i(x, y) when w_i(x, y) ⩾ 3σ_i,ref and W_i(x, y) ≡ 0 when w_i(x, y) < 3σ_i,ref (e.g., Slezak et al. 1994; Starck et al. 2003), where σ_i,ref is the standard deviation of the simulated w_i,ref(x, y) maps so that 3σ_i,ref is equivalent to a confidence interval of 1.3 × 10⁻³ for Gaussian noise. A similar method was used to determine the actual significances of diffuse sources on ROSAT images (e.g., Rosati et al. 1995; Grebenev et al. 1995; Damiani et al. 1997).

Since the à Trous transform is essentially a nonorthogonal algorithm, the obtained coefficient w_i(x, y), and thus W_i(x, y), may be systematically biased from the exact solution (Starck & Pierre 1998). This bias can be corrected by using an iterative method (Van Cittert 1931; Biviano et al. 1996; Starck et al. 2003; Bourdin et al. 2004), the iteration routine of which is defined as

$$\begin{eqnarray} &&W_i^{n}(x,y) \nonumber\\ &&\; =\left\lbrace \begin{array} {@{}l@{}} W_i^0(x,y) + W_i^{n-1}(x,y) - PW_i^{n-1}(x,y) \quad W_i^{0}(x,y)\ne 0, \\ W_i^{n-1}(x,y) \qquad\qquad\qquad\qquad\quad\quad \; W_i^{0}(x,y)=0, \end{array} \right.\nonumber\\ \end{eqnarray} \tag{ 2 }$$

where Wⁿ_i(x, y) is the corrected detection on scale i after n iteration(s) (W⁰_i(x, y) ≡ W_i(x, y) is the initial run), and P is a nonlinear synthesized operator of inverse wavelet transform, wavelet transform, and filtering operation using 3σ_i,ref as the threshold. Typically the calculation converges within five iterations, which yields a corrected W_i(x, y) map as the unbiased reconstruction of the temperature substructures on scale i (Figure 4). We define the wavelet spectrum as 〈 ∣ W_i(x, y) ∣ ²〉 (Torrence & Compo 1998), and plot it in Figure 5 as a function of scale for each cluster.

Zoom In Zoom Out Reset image size

As a cross-check, we also employ the two-dimensional B-spline function (e.g., Bourdin & Mazzotta 2008) as the scaling function ϕ(x, y) to perform a compact à Trous transform in the Fourier domain (u, v) so that the Fourier transform of the wavelet coefficient w_i(x, y) is directly calculated as $\hat{w}_{i}(u,v) = (1-\hat{h}(2^{i-1}u,2^{i-1}v))\hat{c}_{i-1}(u,v)$, where the filter $\hat{h}(u,v)=\hat{\phi }(2u,2v)/\hat{\phi }(u,v)$ when |u|, |v| ⩽ 1/2 (in this paper the Fourier transform of any function f is written as $\hat{f}$). By carrying out the inverse Fourier transform to $\hat{w}_{i}(u,v)$, w_i(x, y) maps are obtained, and then used to calculate the corrected W_i(x, y) maps following the procedure described above. The corresponding wavelet spectra are illustrated in Figure 5 for comparison, which agree nicely with those obtained with the Gaussian wavelet. We find that the most popular scales of the projected temperature substructures are in the range of about 50–200 h⁻¹₇₀ kpc, and about 70% of these substructures are located at 100–200 h⁻¹₇₀ kpc from the cluster center, exhibiting irregular shapes and spatial distributions (Figure 4).

Power Spectrum. As a straightforward approach, we have also studied the power spectra of the two-dimensional gas temperature distributions. Adopting the method of Schuecker et al. (2004), we determine the flat fields $\bar{T}({\bf r})$ by applying a low-pass wavelet filter with a scale of 2⁸ pixels on the original temperature maps, and calculate the normalized temperature fluctuation distributions as $\delta T ({\bf r}) = T({\bf r})/\bar{T}({\bf r})-1$. By carrying out the Fourier transform of δT(r)

$$\begin{equation} \delta ({\bf k})={1\over S}\int _{S} \delta T({\bf r})e^{i{\bf r}\cdot {\bf k}}d{\bf r}, \end{equation} \tag{ 3 }$$

where S is the projected area, we obtain the power spectrum P(k),

$$\begin{equation} P(k) = \langle\mid {\delta ({\bf k})}\mid ^2 \rangle. \end{equation} \tag{ 4 }$$

Following the method of Schuecker et al. (2004), we carry out Monte Carlo simulations for each cluster to estimate the temperature fluctuation at any given knot r_i based on the temperature error range measured at the same knot (Section 3.3.1). After a randomized value T_err(r_i) is assigned to each knot, a random temperature error map T_err(r) is created by applying Equation (1) to the T_err(r_i) map, in which the uncertainties introduced by data statistics, spectral model fitting, and cell tessellation are all involved. Then, by combining T(r) and T_err(r) maps we create a random temperature map T_ran(r). The 1σ error bars shown in Figure 6 are calculated from the variances of power spectra obtained from 10 such T_ran(r) maps for each cluster. In general, the power spectra P(k) show a significant peak in 50–200 h⁻¹₇₀ kpc, and drop quickly beyond 200–250 h⁻¹₇₀ kpc. This strongly supports our conclusion derived with the wavelet analysis (Figure 5).

Zoom In Zoom Out Reset image size