Decoding the thermal history of the merging cluster Cygnus A

We report on a detailed spatial and spectral analysis of the large-scale X-ray emission from the merging cluster Cygnus A. We use 2.2 Ms Chandra and 40 ks XMM-Newton archival datasets to determine the thermodynamic properties of the intracluster gas in the merger region between the two sub-clusters in the system. These profiles exhibit temperature enhancements that imply significant heating along the merger axis. Possible sources for this heating include the shock from the ongoing merger, past activity of the powerful AGN in the core, or a combination of both. To distinguish between these scenarios, we compare the observed X-ray properties of Cygnus A with simple, spherical cluster models. These models are constructed using azimuthally averaged density and temperature profiles determined from the undisturbed regions of the cluster and folded through MARX to produce simulated Chandra observations. The thermodynamic properties in the merger region from these simulated X-ray observations were used as a baseline for comparison with the actual observations. We identify two distinct components in the temperature structure along the merger axis, a smooth, large-scale temperature excess we attribute to the ongoing merger, and a series of peaks where the temperatures are enhanced by 0.5-2.5 keV. If these peaks are attributable to the central AGN, the location and strength of these features imply that Cygnus A has been active for the past 300 Myr injecting a total of $\sim$10$^{62}$ erg into the merger region. This corresponds to $\sim$10% of the energy deposited by the merger shock.


INTRODUCTION
Galaxy clusters represent the fossil record of billions of years of accumulated physics driving the formation and evolution of largescale structure in the Universe.Shocks from mergers (Markevitch & Vikhlinin 2007) as well as the effects from any central active galactic nuclei (AGN; Dunn & Fabian 2006;Bîrzan et al. 2012) can heat the intracluster medium (ICM) up to temperatures of 10 keV.Such heating mechanisms are crucial to prevent catastrophic radiative cooling of the ICM and reduce the amount of gas flowing into the core (Bîrzan et al. 2004;Rafferty et al. 2006;Dunn & Fabian 2006).Furthermore, gas from the core can be transported out to larger radii by powerful jets from AGN forming cavities in the process (for a detailed review, see McNamara & Nulsen 2007).Such cavities subsequently rise buoyantly through the ICM transporting energy and metals (McNamara & Nulsen 2007;Sanders et al. 2005;Roediger et al. 2007).Thus, better constraints on the energy input into the ICM from both mergers and AGN feedback is essential to understand the evolution of clusters over cosmic time.
Although the key physics involved in these two mechanisms is broadly known, the details are poorly understood.For example, the fraction of energy released by mergers that is deposited into the ICM ★ E-mail: a.majumder@uva.nl(AM) is not well constrained.Furthermore, any additional heating of the ICM by a central AGN complicates our understanding of the relative contribution of these two processes to the overall energy budget.While some past studies have tried to tackle the problem of energy transport due to AGN from the core to the ICM (Ruszkowski & Begelman 2002;Ruszkowski et al. 2004;Sanders & Fabian 2007;Fabian et al. 2017), others have focused on energy flows from outskirts to the core due to mergers (e.g., Churazov et al. 2003;Brüggen et al. 2012).There is, however, a lack of suitable targets where these two important processes can be studied simultaneously to understand their relative importance.
The nearby powerful FRII-class radio galaxy Cygnus A (Fanaroff & Riley 1974) is an excellent object to study these two physical processes in detail.At a redshift of  = 0.0561, it is 10 3 times brighter than any other AGN at similar distances (Stockton & Ridgway 1996;Carilli & Harris 1996) and embedded in a large reservoir of X-ray emitting cluster gas extending over a Mpc from the center (Arnaud et al. 1984;Reynolds & Fabian 1996).The subcluster containing the radio source Cygnus A is furthermore undergoing a major merger at the largest scales that is driving a shock in the outer cluster atmosphere (Owen et al. 1997;Markevitch et al. 1999;Smith et al. 2002;Ledlow et al. 2005).Due to the proximity and extreme brightness of Cygnus A, it is possible resolve X-ray features from kpc scales near the core out to Mpc scales.Hence, it is possible to investigate sig-natures of AGN feedback, structures associated with the merger, as well as any possible interaction between the two.This makes Cygnus A an ideal laboratory to investigate both AGN feedback and merger effects on the ICM in unprecedented detail.
We use 2.2 Ms of archival Chandra ACIS data to create thermodynamic profiles of this merging system.The high spatial resolution of Chandra allows us to create detailed temperature and density profiles along the merger direction that reveal a very hot and complex ICM.In order to separate signatures of heating from the underlying temperature structure of the undisturbed ICM, we construct a model for the system assuming hydrostatic equilibrium.Any deviations from this quiescent model can be attributed to other physical mechanisms such as feedback from the AGN or the merger.Furthermore, we also make use of a 40 ks archival XMM-Newton exposure of Cygnus A to independently verify the detected heating signatures.By modeling the various contributions to the detected level of heating, we can provide new constraints on the effectiveness of both processes in heating the ICM.
We begin in §2 with a description of the basic data reduction used to create the event files and mosaic images.Our paper carefully takes into account any Galactic foreground emission and this is subsequently discussed.This discussion is followed by a general description of the large-scale X-ray morphology in §3.Thermodynamic profiles of the system are presented in §4 for both the merger and non-merger regions of the system.The construction of our model for the undisturbed cluster system, the process to create the simulated datasets, and the subsequent derived model thermodynamic profiles are described in §5.An analysis of the excess observed temperature and energy is shown in §6.The implications of our analysis are summarised in §7.
In this paper we have assumed a ΛCDM cosmology with  0 = 70 km s −1 Mpc −1 , Ω  = 0.3, Ω Λ = 0.7, yielding a linear scale of 65.3 kpc per arcminute at Cygnus A's redshift.All errors reported in this paper are of 1 significance.

Chandra
We have reprocessed and analysed all existing observations of Cygnus A in the Chandra archive.These include deep observations from a multiwavelength observational campaign by Wise (2014).We have, however, excluded ObsID 359 and 1707 from our analysis as they were taken to study photon pileup from the nucleus with short CCD readout times.The rest of the 71 ObsIDs have a combined exposure time of 2.199 Ms.
We use standard reprocessing of each individual observation using CIAO 4.13 (Fruscione et al. 2006) and CALDB 4.9.4.All ObsIDs were examined for contamination due to strong background flares by extracting light curves from a chip not containing the bulk of the cluster emission.Only two ObsIDs (5830 and 6226) had a mild signature of flares and the affected time intervals (∼ 1.5 ks) were removed.The combined exposure after all the filtering is 2.198 Ms.For VFAINT mode observations, we filter background particles in the event files during reprocessing by setting check_vf_pha=yes.
In our analysis, we treat the non-X-ray background (NXB) and X-ray background (XRB) for Cygnus A separately.We discuss modeling the X-ray background in §2.3.To determine the instrumental background, we use the "stowed" ACIS event files which are available in CALDB.This background was measured when the instrument was not exposed to the sky.We created background event files for each individual ObsID by applying the correct gains to the "stowed" files followed by reprojecting them to the same tangent plane.For VFAINT mode observations, we filter the background files further by only selecting status=0 events.

XMM-Newton
Cygnus A has been observed by XMM-Newton for a total of 40 ks (ObsIDs: 0302800101 and 0302800201) in full-frame mode with medium filter.Both ObsIDs were reduced using the XMM-Newton Science Analysis System (SAS) v18.0.0.We use standard processing to obtain MOS and PN event files from the observation data files.Both the ObsIDs were observed in Full-Frame mode and the PN out-of-time events were taken into account as per standard analysis 1 .We remove contamination due to soft-proton flares by building good time intervals.We fit gaussians to the count rate histogram and reject all time intervals with count rate more than  + 2, where  is the mean and  is the standard deviation.We find the rejected time interval to be negligible.
For the XMM-Newton non-X-ray background, we choose the filter wheel closed event files of XMM-Newton2 .These event files were then reprojected to match sky coordinates of individual ObsIDs.The event files were further filtered to remove flares in an identical manner as the observed data.

Modeling the X-ray background
The Cygnus A event files include a significant contribution at soft energies from Galactic emission because of its proximity to the Galactic plane (l=76.19b=+05.76).This Galactic soft emission was modeled using ROSAT all-sky survey (RASS) spectra generated by the X-ray background tool3 (Sabol & Snowden 2019).The RASS spectra were selected from a region with minimal cluster emission to maximize the detected counts from the actual X-ray background.We select an annular region centered around Cygnus A with an inner radius of 1.5R 200 and outer radius of 1.5R 200 + 1 • .The spectra was extracted from this region and analyzed in XSPEC to estimate the contribution from Galactic soft X-ray emission.This procedure yields 16869 counts in the energy band 0.1−2.4keV.
We fit a tbabs(apec+apec+powerlaw)+apec model to this spectrum to model the Galactic halo, low Galactic latitude stars, Galactic non-thermal emission and the local hot bubble, respectively.The first three components are absorbed by the Galactic neutral hydrogen component while the local hot bubble is unabsorbed.We set the abundance value for each apec component to proto-solar values while the redshift was set to zero.The powerlaw photon index of Galactic non-thermal emission was set to 1.46 and the powerlaw norm was set to 8.88 × 10 −7 photons/s/keV/cm 2 at 1 keV (Snowden et al. 2008).As can be seen in Fig. 3, the XRB is much fainter than the ICM signal at all radii considered here, therefore any uncertainties in the background model do not impact our results significantly.The hydrogen column density of Cygnus A was calculated using the method of Willingale et al. (2013) 4 which takes both atomic and molecular hydrogen into account.The weighted effective column density is   = 4.14 × 10 21 atoms cm −2 .We used the latest available abundance model asplund (Asplund et al. 2009) in XSPEC v12.10.1f.
We chose this model as the He/H ratio5 is close to the primordial value from the Big Bang (Steigman 2007).The fit parameters of the Galactic soft X-ray spectrum are reported in Table 1.The quality of the fit is  2 /dof = 0.46/2.In all subsequent spectral analysis, we scale the normalizations of the soft energy model components to the area of the corresponding extraction regions.The derived parameters for the soft X-ray background model are then held fixed while fitting the cluster emission.

Image processing
To create mosaiced images, all Chandra event files were reprojected to a common tangent plane.Count images and exposure maps for each reprojected event file were created in the 0.7 − 7.0 keV band and combined to produce total image and exposure maps.The combined image was inspected for point sources using the script wavdetect.The detected point sources were visually confirmed before being removed from the combined image and from each event file.
Background images from individual background files were created and scaled to match the 9.0−12.0keV emission in the observed event files.This energy band is suitable as it is dominated by NXB due to the negligible ACIS effective area.The resulting background-subtracted, exposure-corrected mosaic for the full field of view around Cygnus A in the energy range 0.7 − 7.0 keV is shown in Figure 1.These combined counts images, background images, and exposure maps were also used to determine the surface brightness profiles in §4.
We use the XMM-Newton data only for spectral analysis.Hence, no images were constructed for both MOS and PN datasets.

LARGE SCALE PROPERTIES
Figure 1 shows a smoothed, exposure-corrected, background subtracted X-ray mosaic of the Chandra data in the 0.7 − 7.0 keV band.It is clear from the image that the large scale emission exhibits complex morphology and the diffuse emission extends to well over a Mpc from the center of Cygnus A. This extended emission has been observed previously (Arnaud et al. 1984;Reynolds & Fabian 1996;Smith et al. 2002), but here we can see it clearly consists of two distinct subclusters in the Cygnus A merging system.In the proceeding analysis, we shall designate the subcluster containing the Cygnus A radio source as CygA and the other subcluster to the NW as Cyg NW.
The X-ray surface brightness of the CygA subcluster is highly peaked around the radio galaxy itself, so we assume the centroid of the subcluster to be co-located with the central AGN (RA=19:59:28.356, Dec=+40:44:02.097;Fey et al. 2004).On the other hand, the surface brightness of the CygNW subcluster is significantly fainter with an extended central concentration.To estimate the centroid of the emission from CygNW, we take the median position of all events within a 1 arcmin radius centered on the apparent peak of the surface brightness.This procedure gives us a centroid of (RA=19:58:47.313,Dec=40:53:17.280)for CygNW.Based on the choice of centroids, the two subclusters are separated by ∼12 arcmin or 784 kpc, oriented along an axis with a position angle of 140 • as measured East from North.The merger axis has an offset of ∼30 • from the major axis of the Cygnus A radio galaxy as defined by the line connecting the radio hotspots.
The large scale morphology suggests that the two subclusters are undergoing a merger.This assertion is supported by previous studies (Owen et al. 1997;Markevitch et al. 1998Markevitch et al. , 1999;;Ledlow et al. 2005) which argue that the system is in the early stages of a merger ∼0.5 Gyr prior to the initial core passage.The region in between the two subclusters shows enhanced surface brightness and contains filamentary structure in the merger region at radii between 200 − 400 kpc from the center of CygA.The position of these filaments is consistent with the location of a region of enhanced termperature seen in lower angular resolution ASCA data (Markevitch et al. 1999;Sarazin et al. 2013).In the next sections, we will investigate these features in more detail.

Region selection
In order to clearly isolate the physical processes occurring in the merger region, we have separated the cluster into two annular sectors, one centered on the merger axis and the other containing the remainder of the undisturbed cluster emission.We define the merger region by choosing a 90 • wedge centered on the merger axis between the two subclusters (see Figure 2).This region encloses the interaction region and extends from the center of the CygA subcluster to the outskirts of the CygNW subcluster.This region was adaptively divided into annular bins centered on the CygA subcluster such that each bin has at least 50,000 counts (SNR ∼ 223).This count rate was chosen to ensure the error on the fitted temperature is less than 5% in the hot interstitial region.This choice gives us 70 distinct annular bins with a bin width that ranges from 2 arcsec near the core to 5 arcsec across the merger region.
The non-merger region for the CygA subcluster was defined by excluding the merger region from the event files.This region is shown in Figure 2. Spectra for the non-merger portion of the cluster were extracted using the same radial binning as the merger region to allow direct comparison between the thermodynamic profiles for both regions.The non-merger region for the CygNW subcluster was chosen in a way to minimise the contribution of emission from CygA.This region is also shown in Figure 2. We then adaptively divided this region into annular bins centered on the CygNW subcluster such that each bin has at least 10,000 counts (SNR ∼ 100).
Due to the low exposure of the XMM-Newton data, we divided the merger region into annular bins such that each bin has at least ∼1600 counts (SNR ∼ 40).This gives us an error on the fitted temperature of ∼10% and results in a bin width of ∼10 arcsec at the core and ∼20-30 arcsec in the merger region.We only analyzed the XMM-Newton data in the merger region as the exposure is too low to be used in the fainter non-merger region.

Spectral analysis
For each Chandra OBSID, individual source, background spectra, and response files within a given extraction region were obtained using the script specextract.The BACKSCAL header of each of the instrumental background spectrum was then scaled so that the source and background count rate match in the 9 − 12 keV band.For XMM-Newton, the XMM-SAS task evselect was used to extract both source and background spectra with FLAG==0.We only select single to quadruple events in MOS (pattern <= 12) and single events in PN (pattern == 0).The redistribution matrix file (RMF) and ancillary response file (ARF) for each spectra were generated with the tasks rmfgen and arfgen.A user defined detector map was provided for each of the tasks to better model spatial variation of the source.The background spectra were then scaled, as discussed in §2.1, and used in all subsequent analysis.
For both Chandra and XMM-Newton, we fit spectra from all ObsIDs for a given extraction region simultaneously using a tbabs(apec) model in XSPEC along with the Galactic soft emission model described in §2.3.We use asplund abundance model in our fitting and use cstat in all our analysis.This setup gives us good spectral fits with the value of C-statistic comparable to the number of degrees of freedom.For Chandra data, the spectral fits had cstat/dof = 15688/13789 in the outermost bin of the merger region and cstat/dof = 24727/21547 in the outermost bin of the non-merger region.For XMM-Newton data on the other hand, the spectral fits had cstat/dof = 245/253 in the outermost bin of the merger region.

Surface brightness
The radial surface brightness profiles in the merger and non-merger directions of Cyg A were extracted using the same radial grid used for the spectral extractions.A comparison of these profiles is shown in Figure 3.Total contributions from the NXB and XRB are shown

3D 2D
! !"# Sphere: Net volume: Figure 4. Schematic diagram for the calculation of deprojected volumes.The emission from the 3D spherical shell (shown in red) excluding the cylinder and spherical caps is detected on the 2D annulus (also shown in red).We assume that the emission from the outer spherical shell (shown with dashed line) does not contribute to the same annulus.This approximation is based on the fact that the ICM is less dense in the outer shell.
along with the resulting background subtracted profile.The Galactic foreground is calculated in the 0.7 − 7.0 keV band by using our Galactic emission model normalized assuming Chandra on-axis effective area and response.The surface brightness of this component is assumed to be constant across the whole field of view.At ∼50 arcsec, we see a sharp jump in temperature as well as a sharp discontinuity in the surface brightness profile.These features are the result of the well-studied cocoon shock (Begelman & Cioffi 1989;Carilli & Barthel 1996;Clarke et al. 1997;Wilson et al. 2006;Snios et al. 2018) being driven into the surrounding ICM by the central AGN (see Figure 1).We also see see a bump in all the profiles at ∼ 800 arcsec because of the presence of CygNW.There are also some fluctuations in the profile along the merger region which correspond to the location of the surface brightness features seen in Figure 1.In contrast, the surface brightness profile along the non-merger direction (Figure 3) is relatively smooth.This smoothness suggests that this region is consistent with being in hydrostatic equilibrium and representative of the undisturbed cluster atmosphere on larger scales.

Temperature
The temperature profiles along the merger and non-merger directions of CygA are shown in Figure 3.It is immediately obvious that the gas on the merger side is hotter everywhere by 0.5 − 2.5 keV compared to the non-merger side.As mentioned above, we also notice the presence of a clear jump in temperature at 60 − 70 arcsec in the merger profile that is consistent with the cocoon shock (see Snios et al. 2018).
It is interesting to note that there are a number of features present in the merger profile between ∼ 100 − 400 arcsec.The temperature is fluctuating in this region by 0.5 − 1.0 keV on top of the underlying profile.We estimate and analyze these fluctuations in detail in §6.The non-merger side however does not show any of these features and has a smooth profile.These fluctuations form a series of ripples.The origin of such temperature features is unclear but they correspond to the filamentary structures discussed in §3 and seen in Figure 1.The implications of such temperature variation is explored in §6.

Density
Using the normalization of the APEC model, we can estimate the gas density by using the following equation: where   is the angular diameter distance to the source in cm,  is the redshift, V is the volume of the region along the line of sight, and  APEC is the normalization of the APEC model obtained from the spectral fit.The electron to proton density ratio is 1.18 for the assumed ASPLUND abundance profile (  ≈ 1.18  ).We estimate the volume, , as the intersection between a spherical shell and a cyllinder, with radii corresponding to the projected 2D annulus in the plane of the sky.Each of the projected annular bins along the 90 • wedge has a line of sight volume of (see Figure 4): where   is inner radius of the annular bin and  +1 is the outer radius as measured from the centre of Cygnus A. This volume is an approximation assuming all the emission in the annulus comes from the spherical shell between   and  +1 .We neglect emission from larger radii that lie along the line of sight behind and in front of this spherical shell.
The pseudo-deprojected electron density profiles along both merger and non-merger directions of Cyg A are shown in Figure 3.The presence of the CygNW subcluster is clearly visible at ∼800 arcsec in the merger profile.The non-merger profile, shown in Figure 3, has a slight bump at the position of CygNW.While we have taken care to exclude CygNW, some minor contribution could remain due to the very extended nature of CygNW.We also note that the actual geometry near CygNW likely differs from the spherical geometry assumed in this section.The exact volume along the line of sight in this region is likely an intersection between a sphere centered on CygA and another sphere centered on CygNW.The volume expression used in Eq. 2, however, does not take into account such a complicated geometry.As such, systematic uncertainties can be associated with density calculations near CygNW.

Abundance
We show the iron abundance profile in the merger and non-merger direction of Cyg A again in Figure 3.There is a sharp drop in abundance at ∼60 arcsec, consistent with the location of the cocoon shock.From previous work, it is known that there is significant non-thermal emission around the cocoon shock (Snios et al. 2018) and it is likely that fitting spectra from regions around the shock can lead to lower abundance measurements, due to imprecision in estimating the thermal continuum.Alternatively, the temperature of the ICM across the shock boundary is different and fitting such multi-temperature gas spectra with a single temperature APEC model can lead to a drop in abundance around the shock.Near the core, the abundance profile is peaked which is expected for a cool-core cluster like Cygnus A (De Grandi & Molendi 2001;Mernier et al. 2017).Although the errors on the profile become large at the outskirts, it is clear that beyond 200 arcsec, the abundance on both the merger and non-merger sides is similar and flat.This supports the early-enrichment scenario where metal enrichment occurs during the proto-cluster phase (Werner et al. 2013).Overall, the abundance profile of Cygnus A is similar in shape to previous statistical studies of abundance profiles (Mernier et al. 2017).

Comparison of temperature profile with XMM-Newton data
The temperature profile of Cygnus A reveals features consistent with the location of filamentary structure on the merger side as well as an overall enhancement due to the merger.As a further check, we compare our Chandra merger temperature profile with XMM-MOS and PN data as it provides an independent verification that the observed features are physical and not due to issues with the Chandra data extraction or analysis.For this comparison, we must take into account the systematic differences in fitted temperatures determined with MOS, PN and ACIS spectra.Such differences occur because of cross-calibration uncertainties of the effective area between different detectors.Hence, we need to scale the Chandra and PN temperature profile before comparing it with the XMM-MOS profile.The scaling relationship is obtained from Schellenberger et al. (2015): where the values of  and  depend on the type of detector and energy band.For the hard band and ACIS-EMOS detectors, the best-fit values are:  = 1.028,  = −0.048.The values for EMOS-EPN detectors in the same band, on the other hand, are:  = 0.940,  = 0.038.We used the hard band values for the Cygnus A merger region as it is very hot (≥ 6 keV).Using these values, the scaled Chandra profile is calculated and shown in Fig. 5.The XMM-MOS profile is also shown in the figure for comparison.The overall large scale temperature profile from both Chandra and MOS are in agreement with each other.All the temperature features discussed in §4.2 are also visible in the MOS profile.Although there are significant uncertainties, it is evident that the ripples are situated in a comparable location to those in Chandra.We further show the XMM-MOS profile overlaid on XMM-PN profile in Fig. 5.The PN temperature profile is consistent with MOS as well.This gives us an independent confirmation that the features are real and not due to data analysis techniques or artifacts of one particular instrument.

SIMULATION SETUP
The source of heating of the ICM along the merger direction is not obvious.As the system is undergoing a merger and there is AGN activity at the center of Cygnus A, it is possible that the ICM has been heated by both of these two processes.The central objective of our work is to study how much energy is added to the ICM by such processes and what is their relative contribution.Thus, it is imperative to compare the CygA-CygNW merging system with a model where the ICM is neither being heated due to interactions of the two subclusters nor by the central AGN.We achieve this by building a model of both CygA and CygNW subclusters based on their locations and physical properties as determined outside the merger region.We use this model to create simulated Chandra datasets which we then analyze in the same way as the actual data.The temperature difference between the actual data and the model will be the excess heating due to the combined effects of the merger and the AGN.The creation of the simulated Chandra data and subsequent analysis is described in more detail in the following sections.

Source model
We construct a simple spherical model of a cluster using pyXSIM 6 , an implementation of the PHOX algorithm, designed to create mock X-ray observations (Biffi et al. 2012(Biffi et al. , 2013)).We set up a spatial grid of varying resolution, with finer grids at the center of the subclusters and coarser grids at outskirts.The resolution is chosen such that the pixels are no larger than 1/3rd of radial size of the spectral regions (as discussed in §3.1).This is achieved using an adaptive mesh refinement (AMR) grid (Berger & Colella 1989) that is available inside the yt module 7 (Turk et al. 2011).The density, temperature and metallicity in each cell is set according to the following profiles (Vikhlinin et al. 2006;Mernier et al. 2017): (5) where  = (/  )   and r 500 is the distance at which mean density of the cluster is 500 times the critical density of the Universe.We set density, temperature and abundance to zero beyond  200 which we set as our simulation boundary for each cluster.Using values from Halbesma et al. (2019), we set the boundary of CygA and CygNW to 1831 and 1609 kpc respectively.Next, we create photons from each cell in this sphere according to an APEC model (Smith et al. 2001) based on AtomDB (Foster et al. 2012) in the energy range 0.7 − 7.0 keV.These photons are then projected onto the sky plane and absorbed with a hydrogen column density of 4.14 × 10 21 atoms cm −2 .We use the XSPEC abundance model asplund and hydrogen column model wabs for our simulation.Finally, the projected photons are written to a SIMPUT 8 file as a photon list.This SIMPUT file is then used as input to the MARX instrument simulator for creating realistic Chandra observations.

Instrument simulator
We use the instrument simulator MARX (Davis et al. 2012)   aspect solution, exposure, detector type and detector offsets.The start time of the simulated ObsIDs is matched with that of the actual ObsIDs to correctly set the effective area at the time of observation.
The MARX tool marx2fits is then used to convert simulation files to a FITS event file.In order to reproject all the event files to a common tangent point, we use the CIAO tool reproject_obs that takes care of the offsets between various ObsIDs and modifies the pixel coordinates accordingly.Finally, we produce counts images and exposure maps of each reprojected ObsID, and add them up using the tool flux_obs.We set the detsubsys parameter of mkinstmap to UNIFORM in order to model the lack of spatial variation of quantum efficiency in MARX simulations.These datasets were then used for all subsequent analysis.

Input thermodynamic profiles
To define the input properties of the models, we fit Equations 4, 5 and 6 to the thermodynamic properties in the non-merger regions of both the CygA and CygNW subclusters.However, these input profiles are in 3-dimensional space and hence can not be directly fitted to the projected profiles shown in section 4. As such, we first convert the APEC normalization      from spectral fits to      by dividing out the extraction area.Equation 4 can then be integrated along the line of sight and fitted to this new profile.Similarly, the 3D temperature and metallicity profiles in equations 5 and 6 can be converted to projected 2D profiles by: where  () can be a 3D profile like temperature or metallicity,     () their 2D analog and dV the volume along line of sight.
The weighting function  is assumed to be proportional to emissivity of each gas element along the line of sight  ∝     √  (Sarazin 1986).These profiles can now be easily fitted to the projected temperature and metallicity profiles along the non-merger direction of CygA and CygNW.These fits are shown in Figure 6.As there is evidence of contamination of CygNW in the density profile beyond 600 arcsec (see §4.3), we limit the CygA fits below this distance.The CygNW profiles are calculated up to the edge of field of view.The fitted parameter values obtained from such fits are reported in §A.These parameter values completely define the input profiles and can be used to start the simulation.

Simulation results
The exposure-corrected image from the simulation is shown in Figure 7.The two subclusters overlap in the region in between which suggests that the surface brightness and temperature profile in this region should be different from individual subcluster profiles.This effect is shown in Figure 8 where a comparison of SB and temperature profiles is shown.These profiles have been calculated by extracting simulated counts along the overlapping region using the same annular bins discussed in section 3.1.The temperature profiles are calculated by fitting a single temperature absorbed APEC model.We use the wabs model for absorption and keep the column density fixed at   = 4.14 × 10 21 atoms cm −2 .We obtain fits with cstat/dof = 38543/30167 in the outermost bin.
We note that the surface brightness in the 'Two subcluster model' is simply a sum of the surface brightness of two individual subclusters.The surface brightness from the 'Two subcluster model' is lower than that derived from the Chandra data by a factor of a few (See Figure 8).This may mean that there is additional gas in the merger region than predicted by the input thermodynamic profiles that were obtained by fitting density profiles on the non-merger side of the two subclusters (Figure 6).Furthermore, there maybe additional heating due to the merger.These two effects explain why the model surface brightness does not match the surface brightness from the data in the merger region.
On the other hand, the temperature profile for the two subcluster scenario is not the average of the temperatures from two clusters as one may trivially expect.This is because the net temperature profile for the two subcluster model is determined by the underlying 3D temperature and density profile.We find that the two subcluster temperature profile agrees with the single subcluster profiles near the core of the subclusters.However, it is lower in the overlapping region than what one may expect when only the CygA subcluster is present.This is a non-trivial effect due to the presence of two non-interacting gaseous spheres.

ENERGY EXCESS IN THE MERGER REGION
From §5.4, it is clear that the interaction in the merger region creates a substantial change in the properties of the gas.We investigate this interaction by subtracting the two subcluster model temperature profile from the Chandra temperature profile along the merger direction.This excess temperature is shown in Figure 9.We further show the excess density in the merger region by subtracting the density in the two subcluster model from the Chandra density profile.With both excess temperature and excess density, we can derive the excess energy in the merger region using the following equation: Here (Δ) Data is obtained by multiplying density and temperature from the Chandra data.(Δ) Model is similarly obtained by multiplying density and temperature from the model. corresponds to effective volume along the line of sight (See Eq. 2).We note that assuming a spherical geometry is a major simplification and the true geometry may differ significantly.Such information about the geometry can only be obtained using numerical simulations assuming different geometries and comparing to observation.This is beyond the scope of this paper.Assuming such a simplification however, the excess energy from the merger is shown in Figure 9.It is hard to visualise any feature from this figure however.Hence, we divide this excess energy by the thermal energy predicted by our two subcluster model and show it in the same figure.We also show a gaussian smoothed Δ/ Model curve for better visualisation of the fluctuations.
In Figure 9, we see several interesting features.Firstly, there is an overall temperature, density and energy excess over the entire merger region.Secondly, there are multiple fluctuations in both temperature and in the Δ/ Model plot.These fluctuations are associated with the prominent filamentary structures seen in the Chandra image.The positive excess everywhere is expected from detailed numerical simulation of the CygA system (Halbesma et al. 2019), but smaller scale temperature fluctuations observed in Figure 9 are not easily explained by the merger.We provide a hypothesis for the cause of these fluctuations and model the Δ/ Model curve in the next section.

Merger and AGN contribution to excess energy
The excess energy can be due to two different possible sources.A numerical simulation of the Cygnus A merging system (Halbesma et al. 2019) shows that a smooth increase of temperature in the merger region is expected.This implies that the merger can partly explain the energy excess we see in the merger region.We model such a smooth increase with a simple 2nd-order polynomial expression.We assume that the effect of the merger is negligible at the cocoon shock as any energy increase in this region is likely to be dominated by the current epoch of AGN activity.The relative energy increase due to the merger can thus be modeled as: The observed energy fluctuations however are not easily explained by the merger.We hypothesize that these fluctuations are caused by previous outbursts from the central AGN.There is evidence that the power injected by the central AGN is significantly higher than the core cooling luminosity (See Fig. 8 of Ubertosi et al. 2023).Excess power from outbursts may thus leave their imprint in the ICM in the form of energy fluctuations.We model these outbursts with a set of gaussians on top of the smooth temperature increase due to the merger.
In Figure 10, we show such toy models with different numbers of gaussian-shaped outbursts.Multiple outburst components are necessary to reduce the fit residuals, as shown in the figure.However, we take caution to not overfit the data with more components than what is necessary.To do this, we use a statistical technique called Bayesian Information Criterion (BIC).The Bayesian information criterion can be written as: Here  is number of parameters,  is number of data points and  2 =  (  −   ) 2 / 2  , where   is the th data point,   is the th model value and   is the error on th data point.The BIC penalises  2. We also derive properties of all outburst components for completeness sake.Including more outburst components leads to an increase in the BIC value and thus they are not favored.
We note that the excess energy bump at the CygNW position is likely due to residual emission from CygNW rather than any effect due to the AGN.CygNW is an extremely diffuse cluster undergoing tidal disruption and hence the true geometry of CygNW is probably different from the spherical geometry assumed in the simulation.
We multiply the relative excess energy with the energy predicted by our two subcluster model to get the absolute excess energy due to the merger and outbursts.We calculate total excess energy of each component by adding the excess energy value at every radial bin.We find that the merger is adding energy up to (3.2 ± 0.3) × 10 62 erg into the ICM.The total energy injected by all the outbursts is (3.1 ± 0.7) × 10 61 erg, which is ∼10% of the merger energy.Furthermore, we calculate the full-width at half maximum (FWHM) for each gaussian component and derive outburst timescale as Δ = 2 × FWHM/  , where   is local sound speed.The age of each outburst is calculated by noting the distance of the fitted gaussian peak from the central AGN and integrating the time taken by a sound wave to reach that distance assuming that the density and temperature varies according to the profiles discussed in §4.Finally, we calculate the power of each outburst by dividing its energy output by the corresponding outburst timescale, Δ.All the derived properties for each outburst are reported in Table 2.The mean values and their errorbars are calculated by assuming that fit parameter values are gaussian distributed around their best-fit values and calculating 1000 realisations of the fit parameters.outburst timescale − older and younger events have similar outburst periods within 2 errors.It does seem though that the power of the outbursts may have increased with time.Such a scenario is consistent with the result found by Chon et al. (2012).However, we have only 3 possible outburst candidates so any evidence is marginal at best.Also, if the 2nd and 3rd Gaussian in Table 2 are real, then we can say that CygA has been consistently active for the past ∼ 300 Myr.

The results in
With the exception of the most recent episode associated with the cocoon shock, the power for other two outbursts does not show any significant evolution over the past ∼ 300 Myr.The powers quoted in Table 3 are lower limits on the outburst energies since they assume that the underlying large-scale temperature increase is entirely due to the merger.Although well motivated by numerical simulations, the shape and amplitude of the merger parametric model as a function of radius is somewhat arbitrary.For instance, if we assume that there is little merger contribution beyond a few hundred arcsec, the derived outburst energy values will be higher.Numerical simulations will be needed to explore such possibilities and are beyond the scope of this paper.
In this study, we used the thermodynamic profiles from the nonmerger side of CygA and CygNW to create our simulated subclusters.In doing so, we assumed the non-merger side of CygA does not contain AGN heating.We make this assumption as we do not see signatures of outbursts on the non-merger side with the present data quality.Nevertheless, the AGN's outbursts may also be heating the ICM on the non-merger side.In that case, the derived outburst energy in this section is a lower estimate than the actual total outburst energy.

CONCLUSIONS
In this work, we have used 2.2 Ms Chandra and 40 ks XMM-Newton data to study the Cygnus A merging system.Our primary objective was to study the hot ICM in between the CygA and CygNW subclusters and investigate the mechanisms responsible for heating the ICM.We resolved the X-ray features in the ICM from kpc scale to Mpc scale and our results can be summarised as follows:  • The Cygnus A merger system shows diffuse emission extending over a Mpc.The merger is between two distinct subclusters that we call CygA and CygNW.The merger region shows enhanced surface brightness and filamentary structures.In previous low resolution studies (Markevitch et al. 1999;Sarazin et al. 2013), this region has been associated with higher temperature due to the merger.
• We confirm the results seen in previous analysis finding a significant temperature enhancement in the merger region relative to the surrounding cluster gas outside this region.With a much deeper exposure, however, we also notice temperature fluctuations of the order of 0.5 − 2.5 keV on top of the underlying profile.These temperature fluctuations correspond spatially with the filamentary structures seen in the surface brightness maps.
• In order to separate the contribution of the merger from any additional heating mechanisms, we have constructed a model of CygA and CygNW subclusters based on their locations and physical properties as determined outside the merger region.The subclusters in our model do not interact with each other, providing a baseline scenario with respect to which any heating in the merger region of actual data can be compared.The temperature profile in our two subcluster model differs substantially from the case of a single subcluster due to the presence of multi-phase gas.We also conclude that the merger region has more gas than that predicted by smooth hydrodynamic profiles.
• The temperature difference between the actual data and our model, Δ, shows a positive excess at all radii from the center of CygA.The temperature excess varies between 0.5 − 2.5 keV with significant temperature fluctuations.The derived energy excess also shows such fluctuations and when smoothed, the energy excess suggests the presence of ripples on top of a smooth underlying profile.We hypothesise that this excess energy can be best explained by heating due to a merger shock and a series of AGN outbursts that have propagated a distance of few hundred kpc from the center of Cygnus A.
• The energy excess is best modeled by a set of two gaussians representing AGN outbursts and a smooth polynomial representing the ongoing merger.There is weaker evidence that there may have been two more outbursts in the history of Cygnus A. If these additional outbursts are indeed present, Cygnus A has been active over the past 300 Myr with a duty cycle of 30 − 200 Myr.The smooth merger component, on the other hand, injects up to ∼10 62 erg into the ICM.Furthermore, the total outburst energy can be ∼ 10% of the merger energy or higher.Assuming our model of merger heating + AGN outburst is correct, this suggests that AGN outbursts can be a significant driver of ICM heating over the lifetime of a cluster.

Figure 1 .
Figure 1.Exposure-corrected, background-subtracted 0.7-7.0keV surface brightness mosaic for the full field of view covered by the existing 2.2 Ms Chandra data.The field measures 0.8 × 0.9 deg 2 and has been smoothed with a  = 10 ′′ Gaussian to enhance the appearance of the diffuse emission.The merging subclusters, CygA and CygNW, are visible along with a filamentary structure of material in the interstitial region between them.The bright X-ray emission associated with the Cygnus A radio galaxy is clearly visible in the centre of the CygA subcluster.These features are discussed in §3 and §4.

Figure 2 .
Figure 2. The regions used for analysing the merger side of CygA subcluster are shown.Middle: The regions used for analysing the non-merger side of CygA subcluster are shown.Bottom: The regions used for analysing the nonmerger side of CygNW subcluster are shown.All three images have been exposure-corrected and background-subtracted in the 0.7-7.0keV band.They have been smoothed with a  = 10 ′′ Gaussian to enhance the appearance of the diffuse emission.

Figure 3 .
Figure 3.The left panels correspond to profiles along the merger region while the right panels correspond to profiles along the non-merger region of CygA.First row: Surface brightness profile as a function of projected radius.Second row: Projected gas temperature as a function of projected radius.Third row: Pseudo-deprojected electron density as a function of projected radius.Fourth row: Metallicity profile as a function of projected radius.The dotted and dashed lines show the location of cocoon shock in CygA and the peak of the CygNW subcluster, respectively.

Figure 5 .
Figure 5. Left: Projected temperature profile across the CygA-CygNW merger region based on spectral fits to 40 ks observations from XMM-MOS.The scaled Chandra profile has been shown for comparison.Right: Projected temperature profile across the CygA-CygNW merger region based on spectral fits to 40 ks observations from XMM-PN.The MOS temperature has been scaled as well to compare with the XMM-PN profile.The dotted line shows the approximate location of the cocoon shock while the dashed line shows the location of CygNW subcluster.

Figure 6 .
Figure 6.The left panels show profiles along the non-merger direction of CygA while the right panels show profiles along the non-merger direction of CygNW.First row: Projected density profiles as a function of projected radius.Second row: Projected gas temperature as a function of projected radius.Third row: Projected metallicity profile as a function of projected radius.The fitted, projected profiles discussed in §5.3 are shown in black.The corresponding 3D profiles are shown in green.The residuals for each fit are shown immediately below each figure.

Figure 7 .
Figure 7. Exposure-corrected 0.7-7.0keV surface brightness mosaic for the 2.2 Ms MARX simulation.The field of view is the same as that covered by the observed Chandra mosaic (measuring 0.8 × 0.9 deg 2 ).The image has been smoothed with a  = 10 ′′ Gaussian to enhance the appearance of the diffuse emission.The subclusters, CygA and CygNW, have been modeled as spherically symmmetric.

Figure 8 .
Figure 8. Left: Surface brightness along the merger region for our two subcluster model, CygA only model and CygNW only model.Right: Temperature profile along the merger region for our two subcluster model, CygA only model and CygNW only model.The presence of two clusters of different temperature alters the net temperature in our two subcluster model.The Chandra SB and temperature profiles along the merger direction have been shown for comparison.

Figure 9 .
Figure 9. Upper Left: Excess temperature along the merger direction.Upper Right: Excess total density along the merger region.Lower Left: Excess energy along the merger direction.Lower Right: Excess energy along the merger direction divided by the energy predicted by the two subcluster model.The blue line shows the gaussian smoothed excess energy.The dotted line shows the location of the cocoon shock while the dashed line shows the location of the CygNW subcluster.

Figure 10 .
Figure 10.Top: Residual temperature fitted with only a 2nd order polynomial.This represents the merger heating.Middle: Fit with polynomial + 1 gaussian.AGN outbursts are approximated as gaussians.Bottom: Fit with polynomial + 2 gaussians.The dotted and dashed lines show the location of the cocoon shock and CygNW cluster.

Table 1 .
X-ray foreground components constrained by the RASS spectrum.The temperature of 'low Galactic latitude stars' component was held fixed during the fit.The normalizations are scaled to a 1 arcmin 2 area.
to simulate the Chandra observations.We simulate each ObsID with matching

Table 2 .
Properties of merger and outburst signatures.Δ/ Model is the amplitude of the outburst.Radius is distance of the peak of the outburst from the center of CygA subcluster.FWHM is the full-width at half maximum of the outburst.isthe age of the outburst.Δ is the outburst timescale.isenergy output from the outburst.isthe power generated by the outburst.The merger excess energy is modeled by a polynomial and we report the output energy,   .The 1st gaussian peak is residual from CygNW (see §6). with the  log() term and optimizes the number of parameters necessary to reduce the residuals.The number of components for which the BIC is minimum is considered the best model.We note the BIC values of our models in Table3.The mean values and their errorbars are calculated from BIC values of 1000 realisations that assumes the Δ/ values are gaussian distributed around their mean.It is clear that four gaussians along with the merger component leads to the minimum BIC value.But the BIC values for 'Polynomial + 3 gaussians' and 'Polynomial + 4 gaussians' models are within 1 error of the BIC value for 'Polynomial + 2 gaussians' model.Thus, we note these two components as possible outbursts in Table Table 2 have interesting implications regarding the history of Cygnus A. The derived ages of the outbursts imply that Cygnus A has been active for the past ∼ 300 Myr with duty cycle of ∼ 30 − 200 Myr.These values of duty cycles are similar to the values determined for AGN-induced cavity systems in other cluster cores (for reviews, see McNamara & Nulsen 2007, 2012).The measured ages in Table 2 suggest that the duty cycle between outbursts could have reduced over time.We do not see any temporal evolution of

Table 3 .
Bayesian Information Criterion (BIC) values for different models.