# CHANG-ES. XIX. Galaxy NGC 4013: a diffusion-dominated radio halo with plane-parallel disk and vertical halo magnetic fields

Y. Stein, R.-J. Dettmar, M. Weżgowiec, J. Irwin, R. Beck, T. Wiegert, M.

Krause, J.-T. Li, V. Heesen, A. Miskolczi, et al.

► To cite this version:

Y. Stein, R.-J. Dettmar, M. Weżgowiec, J. Irwin, R. Beck, et al. CHANG-ES. XIX. Galaxy NGC 4013: a diffusion-dominated radio halo with plane-parallel disk and vertical halo magnetic fields. Astronomy and Astrophysics - A&A, 2019, 632, pp.A13. 10.1051/0004-6361/201935558. [hal-03160288](https://hal.science/hal-03160288)

HAL Id: hal-03160288

https://hal.science/hal-03160288

Submitted on 30 Aug 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# CHANG-ES

## XIX. Galaxy NGC 4013: a diffusion-dominated radio halo with plane-parallel disk and vertical halo magnetic fields*

Y. Stein<fnref n="1" />,<fnref n="2" />, R.-J. Dettmar<fnref n="2" />,<fnref n="3" />, M. Weżgowiec<fnref n="4" />, J. Irwin<fnref n="5" />, R. Beck<fnref n="6" />, T. Wiegert<fnref n="5" />, M. Krause<fnref n="6" />, J.-T. Li<fnref n="7" />, V. Heesen<fnref n="8" />, A. Miskolczi<fnref n="2" />, S. MacDonald<fnref n="5" />, and J. English<fnref n="9" />

Received 27 March 2019 / Accepted 17 June 2019

## ABSTRACT

**Context.** The radio continuum halos of edge-on spiral galaxies have diverse morphologies, with different magnetic field properties and cosmic ray (CR) transport processes into the halo.

**Aims.** Using the Continuum HAloes in Nearby Galaxies – an EVLA Survey (CHANG-ES) radio continuum data from the Karl G. Jansky Very Large Array (VLA) in two frequency bands, 6 GHz (C-band) and 1.5 GHz (L-band), we analyzed the radio properties, including polarization and the transport processes of the CR electrons (CREs), in the edge-on spiral galaxy NGC 4013. Supplementary LOw-Frequency ARray (LOFAR) data at 150 MHz are used to study the low-frequency properties of this galaxy and X-ray data are used to investigate the central region.

**Methods.** We determined the total radio flux densities (central source, disk, halo and total) as well as the radio scale heights of the radio continuum emission at both CHANG-ES frequencies and at the LOFAR frequency. We derived the magnetic field orientation from CHANG-ES polarization data and rotation measure synthesis (RM synthesis). Furthermore, we used the revised equipartition formula to calculate the magnetic field strength. Lastly, we modeled the processes of CR transport into the halo with the 1D SPINNAKER model.

**Results.** The central point source dominates the radio continuum emission with a mean of ~35% of the total flux density emerging from the central source in both CHANG-ES bands. Complementary X-ray data from Chandra show one dominant point source in the central part. The XMM-Newton spectrum shows hard X-rays, but no clear AGN classification is possible at this time. The radio continuum halo of NGC 4013 in C-band is rather small, while the low-frequency LOFAR data reveal a large halo. The scale height analysis shows that Gaussian fits, with halo scale heights of 1.2 kpc in C-band, 2.0 kpc in L-band, and 3.1 kpc at 150 MHz, better represent the intensity profiles than do exponential fits. The frequency dependence gives clear preference to diffusive CRE transport. The radio halo of NGC 4013 is relatively faint and contributes only 40% and 56% of the total flux density in C-band and L-band, respectively. This is less than in galaxies with wind-driven halos. While the SPINNAKER models of the radio profiles show that advection with a launching velocity of ~20 km s⁻<fnref n="1" /> (increasing to ~50 km s⁻<fnref n="1" /> at 4 kpc height) fits the data equally well or slightly better, diffusion is the dominating transport process up to heights of 1–2 kpc. The polarization data reveal plane-parallel, regular magnetic fields within the entire disk and vertical halo components indicating the presence of an axisymmetric field having a radial component pointing outwards. The mean magnetic field strength of the disk of NGC 4013 of 6.6 µG is rather small. Large-scale vertical fields are observed in the halo out to heights of about 6 kpc.

**Conclusions.** The interaction and the low star formation rate (SFR) across the disk of NGC 4013 probably influence the appearance of its radio continuum halo and are correlated with the low total magnetic field strength. Several observable quantities give consistent evidence that the CR transport in the halo of NGC 4013 is diffusive: the frequency dependence of the synchrotron scale height, the disk/halo flux density ratio, the vertical profile of the synchrotron spectral index, the small propagation speed measured modeled with SPINNAKER, and the low temperature of the X-ray emitting hot gas.

**Key words.** diffusion – X-rays: galaxies – galaxies: halos – galaxies: individual: NGC 4013 – galaxies: magnetic fields – radio continuum: galaxies

<footnote n="1">¹ Observatoire astronomique de Strasbourg, Université de Strasbourg, CNRS, UMR 7550, 11 rue de l'Université, 67000 Strasbourg, France
e-mail: yelena.stein@astro.unistra.fr</footnote>

<footnote n="2">² Ruhr-Universität Bochum, Fakultät für Physik und Astronomie, Astronomisches Institut (AIRUB), 44780 Bochum, Germany</footnote>

<footnote n="3">³ Ruhr-Universität Bochum, Research Department: Plasmas with Complex Interactions, 44780 Bochum, Germany</footnote>

<footnote n="4">⁴ Obserwatorium Astronomiczne Uniwersytetu Jagiellońskiego, ul. Orla 171, 30-244 Kraków, Poland</footnote>

<footnote n="5">⁵ Department of Physics, Engineering Physics & Astronomy, Queen's University, Kingston, Ontario K7L 3N6, Canada</footnote>

<footnote n="6">⁶ Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany</footnote>

<footnote n="7">⁷ Department of Astronomy, University of Michigan, 311 West Hall, 1085 S. University Ave, Ann Arbor, MI 48109-1107, USA</footnote>

<footnote n="8">⁸ Universität Hamburg, Hamburger Sternwarte, Gojenbergsweg 112, 21029 Hamburg, Germany</footnote>

<footnote n="9">⁹ Department of Physics and Astronomy, University of Manitoba, Winnipeg, Manitoba R3T 2N2, Canada</footnote>

<footnote>* A copy of the reduced images is available at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsarc.u-strasbg.fr/viz-bin/cat/J/A+A/632/A13</footnote>
# 1. Introduction

Studies of star-forming galaxies in the local Universe provide knowledge of evolutionary properties and their driving physical processes. The investigation of cosmic rays (CRs) transport mechanisms is crucial to understand the disk-halo interaction and the formation of radio halos in galaxies (e.g., Parker 1992). Whereas diffusion-dominated CR transport in galaxies leads to small radio halos (e.g., NGC 7462; Heesen et al. 2016), advection-dominated galaxies show large radio halos (e.g., NGC 4666; Stein et al. 2019). Polarization observations reveal large-scale magnetic fields in many spiral galaxies (e.g., NGC 4666; Stein et al. 2019). One proposed mechanism to explain these large-scale fields is the mean-field $\alpha$-$\omega$ dynamo (see e.g., Ruzmaikin et al. 1988; Beck et al. 1996; Chamandy 2016; Henriksen 2017; Henriksen et al. 2018). Turbulent gas motions, for example driven by supernova explosions, and the Coriolis force ($\alpha$-process), combined with shear motions due to differential rotation ($\omega$-process) amplify and order the magnetic field.

In this paper, the radio halo of the edge-on spiral galaxy, NGC 4013, is analyzed using data from the Continuum HAlos in Nearby Galaxies - an Evla Survey (CHANG-ES) program observed with the Karl G. Jansky Very Large Array (VLA). NGC 4013 is a spiral galaxy with a nearly edge-on orientation at a distance of 16 Mpc. The basic galaxy parameters are listed in Table 1. It exhibits a boxy bulge in optical light (Tully et al. 1996). There is strong evidence that the galaxy harbours a bar potential from CO and optical wavelength measurements (García-Burillo et al. 1999). Martínez-Delgado et al. (2009) found a tidal stellar stream extending 3 kpc above the plane of the disk, which they claim to be caused by a minor merger. NGC 4013 has a strongly warped outer HI gaseous disk (Bottema 1996) and a lagging HI halo (Zschaechner & Rand 2015). A lagging halo can be caused by an active disk-halo interface, such as expected for a “galactic fountain”. In spite of clear hints of interaction, there is no obvious partner in the group of 16 galaxies of which NGC 4013 is a member. Wang et al. (2015) suggest that a one-to-three major merger would explain all observable features, meaning the boxy bulge, the huge HI warp and the tidal stellar streams.

Comerón et al. (2011) found a second thick stellar disk (in Spitzer 3.6 µm-band data), which is smoothly distributed. They suggest that this third component was formed during galaxy formation, when an initially thick disk was dynamically heated by a past merger and then subsequently a new thick disk formed within. Follow-up analysis by Comerón et al. (2018) confirmed the existence of the two thick disks in NGC 4013.

Ho et al. (1997) classified the nucleus of NGC 4013 as a “transition object”. This class of objects was described by Ho et al. (1993) as having an emission line nucleus with [OI] strengths intermediate between those of HII nuclei and low-ionization nuclear emission-line regions (LINERs). Ho et al. (1993) proposed that transition objects can be explained as LINERs whose integrated spectra are weakened or contaminated by neighboring HII regions.

In the VLA 1.49 GHz (λ20 cm) image from Condon (1987), the galaxy shows a boxy shape, suggesting the existence of a radio halo. But the radio structure is blurred due to the large beam size of 48", making a distinction between disk and halo emission difficult. Conversely, the λ6 cm VLA+Effelsberg radio map by Hummel et al. (1991) shows a strong central source to be visible as well as a thin radio disk. Baldi et al. (2018) investigated eMERLIN array data of NGC 4013 and found radio emis-

<table><caption>Table 1. Basic galaxy parameters of NGC 4013.</caption><tr><td>RA</td><td>11h58m31s38 <sup>(1)</sup></td></tr><tr><td>Dec</td><td>+43d56m47s7 <sup>(1)</sup></td></tr><tr><td>Distance (Mpc)</td><td>16.0 <sup>(1)</sup></td></tr><tr><td>Inclination (°)</td><td>89 ± 2 <sup>(2)</sup></td></tr><tr><td>PA (°)</td><td>65 ± 3 <sup>(2)</sup></td></tr><tr><td>Major axis (arcmin)</td><td>5.2 <sup>(3)</sup></td></tr><tr><td>Minor axis (arcmin)</td><td>1.0 <sup>(3)</sup></td></tr><tr><td>v<sub>sys</sub> (km s<sup>-1</sup>)</td><td>831 <sup>(3)</sup></td></tr><tr><td>v<sub>rot</sub> (km s<sup>-1</sup>)</td><td>195<sup>(4)</sup></td></tr><tr><td>SFR (M<sub>☉</sub> yr<sup>-1</sup>)</td><td>0.5 <sup>(5)</sup></td></tr><tr><td>SFRD (10<sup>-3</sup> M<sub>☉</sub> yr<sup>-1</sup> kpc<sup>-2</sup>)</td><td>2.4 <sup>(5)</sup></td></tr><tr><td>Classification</td><td>HII, LINER <sup>(1)</sup></td></tr></table>

**References.** (1) Irwin et al. (2012). (2) This work. (3) NASA/IPAC Extragalactic Database (NED, https://ned.ipac.caltech.edu). (4) Bottema (1996). (5) Star formation rate (SFR) and star formation rate surface density (SFRD) from Wiegert et al. (2015).

sion in a circular shape with a diameter of 2.5"–3" around the optical center. They state that this emission is probably associated with a nuclear star-forming ring and thus internal to the galaxy.

In this paper, we investigate the radio continuum emission of NGC 4013 of the halo and its central source. We present flux density measurements at the two frequencies (1.5 and 6 GHz) of CHANG-ES as well as measurements of supplemental radio continuum data at 150 MHz from the LOw-Frequency ARray (LOFAR) Two-metre Sky Survey (LoTSS, Shimwell et al. 2017). With supplementary X-ray data from Chandra and XMM-Newton we have hints about the central source of NGC 4013. Interestingly, the XMM-Newton data show an extended X-ray halo confined to the middle part of the galaxy where most of the polarized radio emission was detected. We further analyze radio scale heights and the linear polarization. The synchrotron map (with thermal emission subtracted) was used to determine the magnetic field strength. The CR transport was analyzed assuming 1D pure diffusion and advection for the CR electrons.

This paper is organized as follows. Section 2 gives an overview of the radio and X-ray data. In Sect. 3, we present our results. This section contains the Stokes I maps, the scale height analysis, the polarized intensity maps (created both with and without rotation measure (RM) synthesis), the map of the magnetic field strength using equipartition, and the 1D cosmic-ray transport model. In Sect. 4 main results are discussed and in Sect. 5 the results are summarized and conclusions are drawn.

# 2. Data

## 2.1. VLA

The radio continuum data are part of the CHANG-ES survey (Irwin et al. 2012) observed with the Karl G. Jansky Very Large Array (VLA). Observations were obtained in B-, C-, and D-configurations at L-band (1.5 GHz, 500 MHz bandwidth, with a gap of 144 MHz width where strong radio frequency interference is located), and in the C- and D-configurations at C-band (6 GHz, 2 GHz bandwidth). We used 2048 spectral channels in 32 spectral windows at 1.5 GHz and 1024 channels in 16 spws at 6 GHz. All polarization products (Stokes I, Q, U, and V) were obtained. The D-configuration data (Wiegert et al. 2015)
<table><caption>**Table 2.** Observation parameters of NGC 4013.</caption><thead><tr><th>Dataset</th><th>Observing date</th><th>Time on source<br>(before flagging)</th></tr></thead><tbody><tr><td>L-band B-configuration</td><td>11. Aug 2012</td><td>2 h</td></tr><tr><td>L-band C-configuration</td><td>31. Mar 2012</td><td>30 min</td></tr><tr><td>L-band D-configuration</td><td>18. Dec 2012</td><td>20 min</td></tr><tr><td>C-band C-configuration</td><td>19. Feb 2012</td><td>3 h</td></tr><tr><td>C-band D-configuration</td><td>27. Dec 2012</td><td>40 min</td></tr></tbody></table>

are public<fnref n="1" /> and the B- and C-configuration data will become public soon (Irwin et al., in prep., Walterbos et al., in prep.).

The data reduction for Stokes *I* (total power) and Stokes *Q* and *U* (linear polarization) was carried out for all five data sets of NGC 4013 (see Table 2) separately, using the Common Astronomy Software Applications (CASA) package (McMullin et al. 2007) and following the calibration procedures as described in the CHANG-ES paper by Irwin et al. (2013). We used J1331+3030 (3C286) as the primary calibrator (also as the bandpass and polarization angle calibrator), J1246–0730 as the secondary calibrator, and J1407+2827 as the zero polarization calibrator. The calibrated data from the different configurations were then combined for C-band and L-band and were used for imaging Stokes *I*, *Q*, and *U*.

Large-scale flux density for extended structures is missed by the VLA at C-band in C- and D-configuration for structures larger than 240". In L-band B-configuration missing flux density is expected for structures larger than 120" and C- and D-configuration for structures larger than 970". As the radio extent of NGC 4013 in the CHANG-ES bands is 205", we expect no significant influence by missing spacings in the configuration-combined images of both bands (potential missing spacings due to the B-configuration in the L-band setup are compensated by the inclusion of the other two configurations).

The Stokes *I* maps were produced by cleaning with a robust weighting (Briggs 1995) parameter of zero, as implemented in CASA and uv-tapering of 18 kλ. The polarization and magnetic field orientation maps were created from the Stokes *Q* and *U* maps, which were cleaned with a robust parameter of two in order to be more sensitive to faint structures.

Table A.1 shows all rms values for the different maps produced for this study. The achieved rms from the combined (C- and D-configuration) C-band data for Stokes *I* (rob 0, 18 kλ taper) is 4.3 µJy beam⁻<fnref n="1" /> with a beam size (i.e. resolution/FWHM) of 5.2" × 5.3". The resulting rms of the combined (B-, C- and D-configuration) L-band data for Stokes *I* (rob 0, 18 kλ taper) is 16.7 µJy beam⁻<fnref n="1" /> with a beam size of 5.9" × 6.0". Smoothed images of Stokes *I* of both bands were produced for the scale height analysis and the thermal and nonthermal separation and the following determination of the magnetic field strength and the CR transport analysis. For RM-synthesis, Stokes *Q* and Stokes *U* were imaged for each spw with a robust weighting parameter of two. The technique was performed following Stein et al. (2019).

## 2.2. LOFAR

The LoTSS is a sensitive, high-resolution, low-frequency (120–168 MHz) survey of the northern sky using LOFAR (Shimwell et al. 2017). LoTSS data of NGC 4013 were obtained

<table><caption>**Table 3.** Parameters of the XMM-Newton observations of NGC 4013.</caption><tr><td>Galaxy</td><td>NGC 4013</td></tr><tr><td>ObsID/obs.date</td><td>0306060201/13.11.2005<br>0306060301/15.11.2005</td></tr><tr><td>n<sub>H</sub><sup>(a)</sup></td><td>1.26 × 10<sup>20</sup> cm<sup>-2</sup></td></tr><tr><td>MOS filter/obs.mode</td><td>Medium/Full Frame</td></tr><tr><td>pn filter/obs.mode</td><td>Medium/Full Frame</td></tr><tr><td>Total/clean MOS time</td><td>94/79 ks</td></tr><tr><td>Total/clean pn time</td><td>71/67 ks</td></tr></table>

Notes. (a) Weighted average value after LAB Survey of Galactic H I Kalberla et al. (2005).

from data release 1 (DR1, Shimwell et al. 2019)<fnref n="2" />. In DR1, images with two resolutions of 6" and 20" are provided. We use the 6" image of NGC 4013 and applied Gaussian smoothing of 10" for the scale height analysis and of 18" to show the extent in comparison to C-band.

## 2.3. XMM-Newton

The X-ray analysis of the emission from NGC 4013 was performed with the use of the XMM-Newton archive data (see Table 3) processed with the SAS 15.0.0 package (Gabriel et al. 2004). Using the standard reduction procedures, event lists for two EPIC-MOS cameras (Turner et al. 2001) and the EPIC-pn camera (Strüder et al. 2001) for each of the two observations were produced and subsequently filtered for periods of high-energy emission.

To increase the sensitivity of the resulting spectrum, information from all three EPIC cameras was used in the spectral analysis. Spectra were extracted from each of the event lists from a circular central region having a radius of 12", that allowed a full sampling of the PSF of the instrument. The position of the region was chosen to match the radio center of the galaxy. The background spectra were obtained using blank sky event lists (see Carter & Read 2007), using the same procedures as in the case of the source event lists. For each spectrum ancillary files (response matrices and effective area files) were used to correct for the instrumental effects. The spectra from all three EPIC cameras of each observation were then merged using the SAS task *epicspeccombine*. Merging of the spectra was accompanied by subtraction of the background spectra. It is important to note, that the subtraction takes place only after both source spectra and background spectra are combined. This leads to higher signal-to-noise ratios in both resulting spectra, assuring best subtraction of the background. XSPEC 12 (Arnaud 1996) was used to fit model to the final spectrum of the core region of NGC 4013.

## 2.4. Chandra

Further supplemental X-ray data were gained from Chandra. The Chandra data were used to detect point sources reliably due to their high angular resolution. The parameters of the Chandra observations are presented in Table 4. The data were taken from (Li & Wang 2013), where the details about the observations and data reduction techniques, including calibration, can be found.

<footnote n="1">¹ CHANG-ES data release I available at www.queensu.ca/changes.</footnote>

<footnote n="2">² LoTSS data are publicly available at: https://lofar-surveys.org.</footnote>
<table><caption>**Table 4.** Parameters of the *Chandra* observations of NGC 4013.</caption><thead><tr><th>Obs ID</th><th>Instr.</th><th>Exp. (s)</th><th>PI</th><th>Data mode</th><th>Exp mode</th><th>Start date</th></tr></thead><tbody><tr><td>4739</td><td>ACIS-S</td><td>79.10</td><td>Strickland</td><td>FAINT</td><td>TE</td><td>2004-02-03 09:31:16</td></tr><tr><td>4013</td><td>ACIS-S</td><td>4.90</td><td>Satyapal</td><td>FAINT</td><td>TE</td><td>2003-03-16 17:09:05</td></tr></tbody></table>

<figure><img src="image_2.png" /><figcaption>**Fig. 1.** NGC 4013 total intensity image of combined (C- and D-configuration) VLA C-band data at 6 GHz. Contours start at the $3\sigma$ level with an rms noise $\sigma$ of 3.6 $\mu$Jy beam$^{-1}$ and increase in powers of 2 (up to 256). The beam size is $5.2'' \times 5.3''$ and is shown in the bottom left corner of the image. A tapering of 18 k$\lambda$ was used, while the robust parameter was set to zero.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** NGC 4013 total intensity image of combined (B-, C- and D-configuration) VLA L-band data at 1.6 GHz. Contours start at the $3\sigma$ level with a $\sigma$ of 16.7 $\mu$Jy beam$^{-1}$ and increase in powers of 2 (up to 128). The beam size is $5.9 \times 6.0''$ and is shown in the bottom left corner of the image. A tapering of 18 k$\lambda$ was used, while the robust parameter was set to zero.</figcaption></figure>

# 3. Data analysis and results

## 3.1. Radio continuum total intensity maps

The total power contours from Stokes I are overlaid on the optical SDSS images made from the ugr-filters using the formulas of Lupton (2005)$^3$. The Stokes I contour map of C-band is displayed in Fig. 1, the one of L-band in Fig. 2 and the total intensity LOFAR image in Fig. 3. The radio emission is dominated by the central point source, which is seen in all frequencies. Especially in C-band, the emission from the disk is weak and only a thin radio disk and a small radio halo are observed. The presented CHANG-ES images are produced with a uv-tapering of 18 k$\lambda$ to

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** NGC 4013 total intensity image of LOFAR data at 150 MHz. Contours start at the $3\sigma$ level with a $\sigma$ of 94.3 $\mu$Jy beam$^{-1}$ and increase in powers of 2 (up to 32). The beam size is $6.0'' \times 6.0''$ and is shown in the bottom left corner of the image.</figcaption></figure>

<table><caption>**Table 5.** Integrated flux densities of NGC 4013.</caption><thead><tr><th>Observation</th><th>Freq (GHz)</th><th>Flux density (Jy)</th></tr></thead><tbody><tr><td>This work <sup>(1)</sup></td><td>0.15</td><td>0.218 ± 0.021</td></tr><tr><td>Condon et al. (1998)<sup>(2)</sup></td><td>1.40</td><td>0.0405 ± 0.0023</td></tr><tr><td>This work <sup>(3)</sup></td><td>1.58</td><td>0.0368 ± 0.0018</td></tr><tr><td>This work <sup>(3)</sup></td><td>6.00</td><td>0.0118 ± 0.0006</td></tr></tbody></table>

**Notes.** (1) The LOFAR Two-metre Sky Survey (LoTSS), measured on the $10''$ resolution images. (2) The NRAO VLA Sky Survey (NVSS). (3) CHANG-ES with effective frequency after flagging measured on the $10''$ resolution image.

enhance the extended flux of the galaxy and to achieve similar resolutions for a quantitative comparison.

The L-band radio halo shows an extent of up to 2 kpc above and below the mid plane and the LOFAR radio halo is extended to 3 kpc, where both values were measured within the $3\sigma$-contours. In L-band the radio halo is observable further out in comparison to the C-band map and at 150 MHz even further than in L-band. This is the expected behaviour for spectral aging, where older cosmic ray electrons (CREs) corresponding to lower energies (emitting at lower frequencies) are transported away from their sources within the disk of the galaxy, while losing energy due to various loss processes. No boxy structure is found, in contrast to what the observations of Condon (1987) suggest, but a slightly dumbbell-shaped thick disk can be seen in the L-band and 150 MHz maps.

With the total flux density of $11.8 \pm 0.6$ mJy in C-band and $36.8 \pm 1.8$ mJy in L-band the radio spectral index is $\alpha = -0.85 \pm 0.05$. The spectral index between 150 MHz and 1.58 GHz is $\alpha = -0.76 \pm 0.05$. There is only one additional flux density measurement in the literature, which is from NVSS (see Table 5). This value of $40.5 \pm 2.3$ mJy at 1.4 GHz is in agreement with the flux densities obtained in this work. The observed radio spectral indices are in good agreement with expected values for

<footnote>$^3$ www.sdss.org/dr12/algorithms/sdssUBVRITransform/#Lupton2005</footnote>
<table><caption>**Table 6.** NGC 4013 flux density measurement of the central source.</caption><thead><tr><th></th><th>Flux<br>(mJy)</th><th>Beam<br>size<br>(")</th><th>Region<br>size<br>(")</th><th>Obs. date</th></tr></thead><tbody><tr><td>BL</td><td>10.09 ± 0.44</td><td>10</td><td>20</td><td>Aug. 11/2012</td></tr><tr><td>DC</td><td>4.52 ± 0.20</td><td>10</td><td>20</td><td>Dec. 27/2012</td></tr><tr><td>CC</td><td>4.52 ± 0.20</td><td>10</td><td>20</td><td>Feb. 19/2012</td></tr><tr><td>LOFAR</td><td>28.3 ± 1.7</td><td>10</td><td>20</td><td>Jan. 21/2015</td></tr></tbody></table>

**Notes.** Measurements of the different array configurations B-configuration *L*-band (BL), D-configuration *C*-band (DC) and C-configuration *C*-band (CC)

synchrotron emission. The low-frequency radio spectral index below 1.5 GHz is slightly steeper (with a difference of $2\sigma$) than the value of $\alpha = -0.57 \pm 0.10^4$ found by Chyży et al. (2018). Also, we see only a small break in the radio continuum spectrum with $\Delta\alpha = \alpha_{low} - \alpha_{high} = 0.09 \pm 0.05$, whereas Chyży et al. (2018) found $\Delta\alpha = 0.18 \pm 0.20$. This suggests different CR electron propagation and energy losses in NGC 4013 compared to other spiral galaxies.

## 3.2. Central point source

The galaxy NGC 4013 was classified as a LINER by Ho et al. (1997) and Irwin et al. (2012). There is an additional X-ray observation published by Dudik et al. (2005). In their study of high-resolution X-ray imaging of nearby LINERS observed by Chandra, a nuclear point source of NGC 4013 was not detected. They concluded either a lack of an energetically significant AGN or a highly obscured AGN with column densities of $11 \times 10^{23} \text{ cm}^{-2}$ to $10.2 \times 10^{24} \text{ cm}^{-2}$, implying that the AGN has a luminosity between $2 \times 10^{38} \text{ erg s}^{-1}$ and $9.5 \times 10^{42} \text{ erg s}^{-1}$. Here, we investigate the “central point source”, meaning a centrally peaked source, although it is not an unresolved source in both dimensions as defined in Irwin et al. (2019). Using the radio continuum data of both CHANG-ES bands and LOFAR data as well as the X-ray data from Chandra and XMM-Newton we further constrain the classification of this galaxy and its central source.

### 3.2.1. Radio continuum

A clear detection of a central point source is visible in both bands of the CHANG-ES data as well as in the LOFAR data. To further investigate the central source, variability in the total intensity of the central source was checked. Therefore, flux densities of the central source were measured in images from the two different VLA array configurations of the *C*-band data, which are observed at different times (see Table 6). As one *C*-band observation was done in the beginning of 2012 and the other in the end of 2012, time variability could be tested. Additionally, the flux density of the central source in B-configuration *L*-band and LOFAR 150 MHz was measured. To do so, all images were smoothed to the same beam size of 10 arcsec in order to overcome beam-dependent errors in the flux density measurement. The flux density was then measured via a Gaussian profile fit to the intensity, using CASA.

<figure><img src="image_5.png" /><figcaption>**Fig. 4.** NGC 4013 Chandra image with combined (C- and D-configuration) VLA *C*-band contours (same data as Fig. 1) at 3- and 30-$\sigma$-levels overlaid. The false-color represents the Chandra color from 0.3 to 7 keV. Gray cycles denote Chandra point-like sources. The image is off-centered to show the location of the two point-like sources in the northern halo of NGC 4013.</figcaption></figure>

We chose a circular region two times larger than the beam size. To measure the error, the chosen region was varied from the beam size to four times the beam size. The difference in the flux density measurements were then used to calculate a mean error. The integrated flux density of the central source of both *C*-band configurations is $4.52 \pm 0.20 \text{ mJy}$. Therefore, no measurable variability is found between our *C*-band observations. The flux density level of the galaxy was $0.14 \pm 0.17 \text{ mJy beam}^{-1}$, respectively. The integrated flux density of the central source in B-configuration *L*-band is $10.09 \pm 0.44 \text{ mJy}$ with a galaxy flux density level of $0.35 \text{ mJy beam}^{-1}$. The flux density of the LOFAR observations is $28.33 \pm 1.73 \text{ mJy}$ with a base level of $2.28 \text{ mJy beam}^{-1}$.

A resulting radio spectral index of $\alpha = -0.66 \pm 0.05$ is calculated between the *C*-band flux density and the B-configuration *L*-band flux density as well as $\alpha = -0.44 \pm 0.05$ between the *L*-band and the LOFAR flux densities. The values suggest a flattening of the spectral index towards lower frequencies due to synchrotron self-absorption, which requires a very compact object and hints at an AGN. As some AGN show circular polarization (Irwin et al. 2018), this was investigated as well. In the Stokes V image, no circular polarized emission could be found. In summary, apart from the fact that the central source dominates the emission of the galaxy and that the spectral index is relatively flat towards lower frequencies, no strong indication that the central object might be an AGN can be found from its radio properties.

### 3.2.2. X-ray

**Chandra.** To further constrain the nature of the central point source, Chandra X-ray data (Li & Wang 2013) on Stokes I contours are shown in Fig. 4. Clearly, there is a Chandra point-like source located at the same position as the central radio source, which is the dominant source within the entire surroundings apart from another less strong point-like source next to it. This finding corroborates the LINER classification since we expect X-ray emission emanating from an AGN but not from a star-burst nucleus.

**XMM-Newton.** The core region of NGC 4013 was fit with a complex model consisting of a gaseous component that accounted for the emission from the hot gas and a power-law model, that accounted for the central source. The mekal model is a model of an emission spectrum from hot diffuse gas based on

<footnote n="4">⁴ Quoted uncertainties for results by Chyży et al. (2018) refer to the standard deviation rather than the smaller standard deviation of the mean.</footnote>
<table><caption>**Table 7.** Physical parameters derived from the model fit to the spectrum of the NGC 4013 core region.</caption><tr><td>kT</td><td>0.27<sup>+0.06</sup><sub>-0.04</sub> keV</td><td>(3.13<sup>+0.70</sup><sub>-0.46</sub> × 10<sup>6</sup> K)</td></tr><tr><td>Mekal flux</td><td>2.00<sup>+0.87</sup><sub>-0.79</sub> × 10<sup>-15</sup> erg cm<sup>-2</sup> s<sup>-1</sup></td><td></td></tr><tr><td>Photon Index</td><td>1.10<sup>+0.22</sup><sub>-0.21</sub></td><td></td></tr><tr><td>Core flux</td><td>2.31<sup>+1.69</sup><sub>-0.99</sub> × 10<sup>-14</sup> erg cm<sup>-2</sup> s<sup>-1</sup></td><td></td></tr><tr><td>Core luminosity</td><td>7.1<sup>+5.2</sup><sub>-3.0</sub> × 10<sup>38</sup> erg s<sup>-1</sup></td><td></td></tr></table>

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** XMM-Newton EPIC spectrum of the core region of NGC 4013.</figcaption></figure>

the model calculations of Mewe and Kaastra (Mewe et al. 1985; Kaastra 1992).

Since we assumed that most of the power-law flux would be attributed to the central source (see also Fig. 4) and given the relatively large uncertainties of its value, we used this flux to calculate the luminosity of the central source. We obtained a value of $7.1_{-3.0}^{+5.2} \times 10^{38}$ erg s$^{-1}$ (Table 7).

The spectrum of the core region of NGC 4013 is presented in Fig. 5. No clear AGN characteristics are visible, but significant hard emission (above 2 keV) and a relatively flat photon index suggest a prominent central source.

The average temperature of the hot gas in the core region of NGC 4013 is $0.27_{-0.04}^{+0.06}$ keV. Assuming a spherical emitting volume of the region, the model fit yields a number density of $5.37_{-1.07}^{+0.85} \times 10^{-3}$ cm$^{-3}$ and the corresponding thermal energy density of $3.50_{-1.11}^{+1.44} \times 10^{-12}$ erg s$^{-1}$.

Figure 6 shows the extended emission from NGC 4013 in the 0.2–1 keV energy band. The spectral region used for the analysis is marked with a red circle. Most of the emission of the hot gas comes from the central parts of the galactic disk with some outflows into the halo visible.

## 3.3. Intensity profiles of the radio disk and halo

### 3.3.1. Scale heights – Strip fitting

The total intensity profiles perpendicular to the major axis of edge-on spiral galaxies are fit by exponential or Gaussian functions to determine scale heights of the thin and thick disk. The thin disk refers to the disk and the thick disk refers to the radio continuum halo of the galaxy, which are the used terms in this work. Following Dumke et al. (1995), the observed intensity profile can be described by the convolution of a Gaussian beam profile with the true intensity profile of the galaxy, which could be either exponential or Gaussian.

<figure><img src="image_7.png" /><figcaption>**Fig. 6.** XMM-Newton EPIC emission of NGC 4013 in the energy band 0.2–1 keV (contours) overlaid on the DSS blue image. The contours start at 3, 3.5, 4.5, $6\sigma$ levels. The image is adaptively smoothed with the largest smoothing scale of 10" and the signal-to-noise ratio of 10. The spectral core region is marked with a red circle.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 7.** NGC 4013 screenshot of the strips used with NOD3 (combined C-band).</figcaption></figure>

We used the “BoxModels” tool within the NOD3 software package (Müller et al. 2017) with the combined configuration maps of C-band and L-band as well as the LOFAR map, which were smoothed to the same beam size of 10" × 10" in order to compare the results. The scale height analysis was carried out in an analogous way to Krause et al. (2018) using five strips perpendicular to the major axis (see Fig. 7). The box height was chosen to be approximately half of the beam size, which was 4" in both bands. The box width was 35" to cover the major axis with the five strips. We note that with the LOFAR data it would have been possible to measure seven strips due to the larger extent of the major axis. A two-component exponential and Gaussian fit was executed for all three frequencies to the five strips of radio intensity profiles of NGC 4013. We achieved a fit with similar $\chi^2$ using Gaussian distributions in comparison with exponentials for the C-band intensity profile, but better results for a Gaussian fit for the L-band and LOFAR intensity profiles.

The parameters of NOD3 are summarized in Table 8. The intensity profiles of the five strips and the two-component Gaussian fits are displayed in Fig. 8 for C-band in Fig. 9 for L-band and in Fig. 10 for the 150 MHz LOFAR data.

The resulting Gaussian scale heights of the strip fitting are presented in Table 9. For the final scale height determination, the central strip was omitted in order to exclude contamination by the central source. The scale height analysis of NOD3 shows an increasing thin disk component from higher to lower frequencies. The ratio between the thin disk scale height and the halo scale height is similar for C-band and L-band and decreases slightly towards the lower frequency range of LOFAR. For the
<table><caption>**Table 8.** NGC 4013 parameters for the scale height fitting.</caption><thead><tr><th>Parameter</th><th>C-band</th><th>L-band</th><th>LOFAR</th></tr></thead><tbody><tr><td>Beam size (")</td><td>10.0</td><td>10.0</td><td>10.0</td></tr><tr><td>Effective beam size (")</td><td>10.2</td><td>10.2</td><td>10.2</td></tr><tr><td>Inclination (°)</td><td>89</td><td>89</td><td>89</td></tr><tr><td>Position Angle (°)</td><td>65</td><td>65</td><td>65</td></tr><tr><td>rms (mJy beam<sup>-1</sup>)</td><td>0.010</td><td>0.025</td><td>0.10</td></tr><tr><td>Galaxy diameter (")</td><td>205</td><td>205</td><td>230</td></tr><tr><td>Box width (")</td><td>35</td><td>35</td><td>35</td></tr><tr><td>Box height (")</td><td>4</td><td>4</td><td>4</td></tr><tr><td>Number of boxes in X</td><td>5</td><td>5</td><td>5</td></tr><tr><td>Number of boxes in Y</td><td>20</td><td>20</td><td>24</td></tr></tbody></table>

thin disk ("disk") scale heights, we find 2.8" ± 0.7" (220 pc) in C-band, 4.7" ± 0.3" (360 pc) in L-band, and 6.1" ± 0.6" (470 pc) at 150 MHz. For the thick disk ("halo") component, the scale heights are 15.3" ± 3.0" (1.2 kpc) in C-band, 25.8" ± 1.5" (2.0 kpc) in L-band, and 40.4" ± 4.0" (3.1 kpc) at 150 MHz. Interestingly, the halo scale heights are similar to the measured 3σ extent measured in Sect. 3.1. This underlines the faintness of the radio halo, which we can observe only to an extent of the scale height.

The radio scale heights are influenced by the contribution from thermal emission. Whereas exponential scale heights seem to be the same within the errors in C-band for NGC 4565, the Gaussian scale heights of the total radio continuum emission in the halo are about 1 kpc smaller compared to nonthermal scale heights (Schmidt et al. 2019). C-band data are most affected by the thermal contribution, so that the halo scale height of 15.3" (1.2 kpc) is a lower limit. (The profiles of the nonthermal emission for C-band (see Sect. 3.6) could only be fit with large uncertainties with a resulting Gaussian halo scale height of 19.8" ± 4.0".) The thermal contribution diminishes at lower frequencies, making the ratio of scale heights at L-band and 150 MHz suited for studying CRE propagation (see below).

Our radio scale heights differ from the ones determined by Krause et al. (2018), the main reason being that we fit double Gaussian functions while they fit exponential ones. In order to make a comparison possible, we show the fitting of two-component exponential functions in Figs. A.1–A.3, for C-band, L-band, and LOFAR, respectively. With the combined configuration C-band data, the two-component exponential fits are comparable to the D-configuration C-band fit from Krause et al. (2018). With the exponential fit in L-band (Fig. A.2) we achieve slightly larger scale heights but they are within the uncertainties of Krause et al. (2018). We conclude that Gaussian scale heights cannot easily be compared to exponential scale heights.

We now compare the exponential scale heights with other values from the literature. Hummel et al. (1991) found a radio scale height of 4" from VLA C-band observations. This is half of our exponential scale height measured for the halo component of 8.1", which, however, equates only to 0.6 kpc and thus shows the halo to be difficult to detect at C-band. With Hi data of NGC 4013, a scale height of 3" (210 pc) was found by Zschaechner & Rand (2015), which is in good agreement with the Gaussian C-band scale height of the thin disk and is also similar to the mean scale height of the thin disk at 5 GHz, for which Heesen et al. (2018) found 0.25 ± 0.13 kpc. Spitzer data at 3.6 µm, analyzed by Comerón et al. (2011), suggest a thin disk of 122.5 pc and a thick disk of 557.5 pc. These values are roughly

half the values found here for the Gaussian C-band radio scale heights.

In summary, the scale height analysis shows Gaussian functions to better represent the vertical radio intensity profiles than exponential functions at L-band and the LOFAR frequency. For C-band, the fit quality is similar. Nevertheless, the Gaussian profile represents the data slightly better visually.

### 3.3.2. Amplitudes and integrated flux densities of Gaussian profiles

With the Gaussian fit to the intensity profile, we gather the scale heights as well as the amplitudes of the Gaussian functions in each of the five strips. In Table 9 we present the mean disk amplitude and halo amplitude as well as the ratio of them together with the mean scale height values. As expected, the disk amplitude dominates in all cases. But even more interesting is the question of which component of the total flux density is the dominating one, i.e. the disk or the halo. With the parameters of the Gaussian fit, we determine the integrated flux density, i.e. calculate the integral of the Gaussian function. To do so, we derive the mean amplitudes and mean scale heights. This is done for four strips omitting the central strip, because we want to find the general Gaussian function describing the galaxy without the central source.

With the Gaussian distribution of the intensity profile for the mean of the four strips:

$$w(z) = \bar{w}_0 \exp(-z^2 / \bar{z}_0^2) \qquad (1)$$

the integrated flux density becomes:

$$W(\text{mJy} \cdot \text{arcsec/beam}) = \sqrt{\pi} \bar{w}_0 \bar{z}_0, \qquad (2)$$

where the mean amplitude $\bar{w}_0$ is measured in mJy per beam area and the mean scale height $\bar{z}_0$ in arcsec. We divide $W$ by the beam area $b$ in arcsec ($b = 113.3 \text{ arcsec}^2$) to obtain the integrated flux density $W$ in mJy per arcsec, i.e. the flux density for a layer of one arcsec thickness perpendicular to the strip, and averaged the values for all strips. To determine the total flux density of the thin disk and the halo, which is to be compared with the measured total integrated flux densities, we multiply with the extent $e$ along the major axis:

$$F(\text{mJy}) = \sqrt{\pi} \bar{w}_0 \bar{z}_0 \cdot e / b. \qquad (3)$$

The extent $e$ is chosen to best represent the mean Gaussian profile of the galaxy. First, we estimate the extent of the galaxy at each frequency (Table 8), then we lower this extent to account for the contributing size of the central source. In the last step we consider the total flux densities of Table 5 as a reference point and increase the value until the calculated and measured total flux densities are similar. We obtain $e = 160'' \pm 15''$ in C-band, $e = 190'' \pm 15''$ in L-band, and $e = 195'' \pm 15''$ at 150 MHz. The total emission of the "disk + halo" is the integrated flux density of Table 5 minus the central source flux density of Table 6, which is $(11.8 \pm 0.6 \text{ mJy}) - (4.5 \pm 0.2 \text{ mJy}) = 7.3 \pm 0.6 \text{ mJy}$ in C-band, $(36.8 \pm 1.8 \text{ mJy}) - (10.1 \pm 0.4 \text{ mJy}) = 26.7 \pm 2.2 \text{ mJy}$ in L-band, and $(218 \pm 21 \text{ mJy}) - (28.3 \pm 1.7 \text{ mJy}) = 190 \pm 21 \text{ mJy}$ for the LOFAR data. These values are presented on the right-hand side of Table 9 along with the flux densities of the central source (from Table 6) for comparison. The results of the integrated flux densities of the disk and halo as well as their ratio are presented. As a consistency check, we added up the flux densities of the disk and halo Gaussian intensity models and confirmed that they agree with the measured total flux densities.
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Strip fitting with NOD3 for five strips of NGC 4013 on the combined (C- and D-configuration) C-band data with a two-component Gaussian fit ($\chi^2_{total} = 0.05$). The first plot is at 75" on the major axis and the last plot is at -75" on the major axis from the center. The red dots represent the mean intensity of each box and the black line is the fit to the data done by NOD3. The black horizontal line shows the size of the beam $\theta$. The parameters of the exponential fit function are given in the boxes on top of each strip. The mean of the four scale heights omitting the central strip as well as the beam size are given in the box to the lower right.</figcaption></figure>

The C-band total flux densities are approximately the same for the disk and the central source, whereas the halo flux density is lower. In L-band, the halo flux density dominates mildly, but still agrees with the contribution from the disk within the uncertainties, while the flux density of the central source is slightly smaller. At 150 MHz, the contribution from the halo to the flux density clearly dominates. The ratio of the integrated flux density of the disk to the integrated flux density of the halo shows a trend with frequency: from 6 GHz to 150 MHz, the ratio decreases from 1.7 to 0.3. The amplitude ratios follow the same trend, decreasing from 9.4 to 1.9. The CREs emitting at lower frequencies have longer lifetimes and thus can be transported into the halo more easily.

### 3.4. Polarization and the large-scale magnetic field structure

The results from the C-band polarization data are shown at high and low angular resolution in Figs. 11 and 12, respectively. In Fig. 11, a symmetric structure in the disk is visible east and west of the central point-like source with the peak of the polarized emission in the east. The orientation of the magnetic field is disk-parallel everywhere and no large-scale polarized structures, besides that of the disk field, are detected. The low-resolution map (beam size of 18" × 18") presented in Fig. 12 reveals new structures. As before, a plane-parallel magnetic field in the disk is observed, but now many polarized features reaching far into the halo with vertical field orientations are seen as well. Because the total radio intensity is not very high (11.8 ± 0.06 mJy in C-band), the polarized intensity is even lower (1.63 ± 0.08 mJy in C-band from Fig. 12). However, strong smoothing (18") in

combination with imaging with a robust weighting parameter of 2 made it possible to reveal the large-scale magnetic field in NGC 4013. For comparison, C-band, L-band and LOFAR Stokes I contours are shown at a $3-\sigma$ level in Fig. 12. The polarized emission extends as far as the LOFAR emission and the L-band outer contour shows an elongation towards the northern polarization “blob”, lending credibility to the detection. The rotation curve of the galaxy indicates that the eastern part is approaching the observer (e.g., Bottema 1996). Thus, within the disk, the approaching side shows more polarized emission than the receding side, similar to NGC 4666 (Stein et al. 2019).

The vertical magnetic field structure in the halo cannot be explained by diffusion alone and indicates an additional outflow. To investigate this further, we follow the strategy described in Damas-Segovia et al. (2016), who used the projected distance of polarized structures together with an estimate of the CRE lifetime to measure the outflow speed. The maximum extent of the north-western polarized emission blob above the disk is about 82" (6.4 kpc). The synchrotron lifetime follows from:

$$t_{\text{syn}} = 1.06 \times 10^9 B^{-3/2} v^{-1/2} \text{yr} \quad (4)$$

with the total field strength $B$ in $\mu$G and $v$ in GHz. We assume a decline of the magnetic field strength of the disk towards the halo, going from the disk value of $B \approx 6.6 \mu$G (see Sect. 3.7) to $B \approx 2 \mu$G at the maximum z-extent of 6.4 kpc. In order to account for inverse Compton losses, we added the average magnetic field strength $B_{\text{halo}} \approx 4.0 \mu$G encountered by CREs travelling into the halo and the equivalent magnetic field strength $B_{\text{CMB}}$ in quadrature, with the latter accounting for the energy
<figure><img src="image_10.png" /><figcaption>Fig. 9. Strip fitting with NOD3 for five strips of NGC 4013 on the combined (B-, C- and D-configuration) L-band data with a two-component Gaussian fit ($\chi^2_{total} = 0.92$). Designations are the same as for Fig. 8.</figcaption></figure>

density of the CMB:

$$ \sqrt{(B_{\text{halo}}^2 + B_{\text{CMB}}^2)} = \sqrt{(4.0^2 + 3.25^2)} \approx 5 \, \mu\text{G.} \quad (5) $$

The synchrotron lifetime then is $t_{syn} \approx 3.8 \times 10^7$ yr ($B \approx 5 \mu\text{G}$, $\nu = 6$ GHz). Thus, the velocity which is needed to transport CREs to 6.4 kpc height in the halo of NGC 4013 is about 160 km s$^{-1}$. This is a low advection velocity in comparison to other galaxies (Heesen et al. 2018).

The X-ray emission, presented in Fig. 6, is mainly located where we observe the strongest polarized emission of NGC 4013 (yellow/orange color in Fig. 12). The southern polarization structure in the halo pointing roughly towards the center (Region I in Fig. 13) seems to be present in the XMM-Newton image as well. The polarized “blob” in the halo to the north (Region II) shows the same orientation of the magnetic field vectors in an area of more than two beams. We find two distinct X-ray emission features within that blob, but these could be associated with two point sources confirmed by Chandra (see Fig. 4), probably background sources. Furthermore, a source is located in the upper part of the blob (RA 11:58:26.87, Dec +43:57:46.14), visible in C-band Stokes I (e.g., Fig. 11) and in the optical (e.g., in SDSS) as a diffuse source, possibly a background galaxy. The possibility that these sources could generate an extended polarized blob seems unlikely. Therefore, we assume this polarization feature to originate from NGC 4013. The polarization structure to the northeast of the galaxy (Region III) reaches far into the halo. This location coincides with an extended HI filament (Bottema 1996) showing that this galaxy seems to have local outflows, which have lower velocities than the escape velocity.

In the L-band data, no polarized emission is detected. Polarized emission in L-band is expected to be affected by depolarization due to Faraday rotation. As this effect is wavelength dependent, the depolarization is stronger in L-band than in C-band. Therefore, images were produced with almost natural weighting (robust = 2), strong uv-tapering, and additional use of RM-Synthesis, but still no polarized flux density within the central disk could be recovered, with the exception of a few patches of low intensity emission. These are at a level of $2\sigma$ (with a $\sigma$ of 29 $\mu$Jy beam$^{-1}$) and consistent with the peaks of the C-band emission.

### 3.5. Rotation measures and the magnetic field configuration in the disk

The polarization angle of an electromagnetic wave experiences a frequency dependent rotation while propagating through a magnetized plasma. This effect is called Faraday rotation with the rotation of the polarization angle $\chi$ being proportional to the wavelength squared and the RM.

We perform RM synthesis following Stein et al. (2019), with the important parameters presented in Table 10. From the data in Xu & Han (2014), the Galactic foreground RM of NGC 4013 is determined. This value ($RM_{foreground} = -8.1 \pm 1.2$) is subtracted from the final RM-cube to correct for the foreground contribution.

The maps of polarized intensity and Faraday RM, both produced with RM-synthesis, are shown in Fig. 14. We overlay the total intensity map with polarization contours, and the intrinsic magnetic field orientation, corrected for Faraday rotation, and the RM map with total intensity contours. As in Fig. 11, the
<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Strip fitting with NOD3 for five strips of NGC 4013 on the LOFAR data with a two-component Gaussian fit ($\chi^2_{\text{total}} = 1.07$). Designations are the same as for Fig. 8.</figcaption></figure>

<table><caption>**Table 9.** Integrated flux densities and scale heights of the thin disk and halo of NGC 4013.</caption><thead><tr><th rowspan="2"></th><th colspan="5">Gaussian fits</th><th colspan="3">From total flux densities</th></tr><tr><th>w&#773;<sub>0</sub><sup>1</sup><br>(mJy beam <sup>(1)</sup>)</th><th>z&#773;<sub>0</sub><sup>(1)</sup><br>(")</th><th>z&#773;<sub>0</sub><br>(kpc)</th><th>flux density<br>(mJy)</th><th>disk + halo<br>(mJy)</th><th>disk + halo<br>(mJy)</th><th>central source<br>(mJy)</th></tr></thead><tbody><tr><td><strong>C-band</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Disk</td><td>0.66 ± 0.08</td><td>2.8 ± 0.7</td><td>0.2 ± 0.1</td><td>4.6 ± 1.7</td><td>7.3 ± 2.6</td><td>7.3 ± 0.6</td><td>4.5 ± 0.2</td></tr><tr><td>Halo</td><td>0.07 ± 0.03</td><td>&gt;15.3 ± 3.0 <sup>(2)</sup></td><td>&gt;1.2 ± 0.2 <sup>2</sup></td><td>2.7 ± 2.0</td><td></td><td></td><td></td></tr><tr><td>Ratio</td><td>9.4</td><td>&lt;0.18 <sup>(2)</sup></td><td></td><td>1.7</td><td></td><td></td><td></td></tr><tr><td><strong>L-band</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Disk</td><td>0.9 ± 0.2</td><td>4.7 ± 0.3</td><td>0.36 ± 0.05</td><td>11.5 ± 1.7</td><td>26.4 ± 4.9</td><td>26.7 ± 2.2</td><td>10.1 ± 0.4</td></tr><tr><td>Halo</td><td>0.2 ± 0.1</td><td>25.8 ± 1.5</td><td>2.0 ± 0.1</td><td>14.9 ± 4.6</td><td></td><td></td><td></td></tr><tr><td>Ratio</td><td>4.5</td><td>0.18</td><td></td><td>0.8</td><td></td><td></td><td></td></tr><tr><td><strong>LOFAR</strong></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Disk</td><td>2.3 ± 1.1</td><td>6.1 ± 0.6</td><td>0.47 ± 0.10</td><td>43 ± 17</td><td>191 ± 20</td><td>190 ± 21</td><td>28.3 ± 1.7</td></tr><tr><td>Halo</td><td>1.2 ± 0.1</td><td>40.4 ± 4.0</td><td>3.1 ± 0.3</td><td>148 ± 11</td><td></td><td></td><td></td></tr><tr><td>Ratio</td><td>1.9</td><td>0.15</td><td></td><td>0.3</td><td></td><td></td><td></td></tr></tbody></table>

**Notes.** <sup>(1)</sup> $\bar{w}_0$ is the mean amplitude omitting the central strip, $\bar{z}_0$ is the mean scale height omitting the central strip. <sup>(2)</sup> Lower limit if only the synchrotron emission is considered.

magnetic field orientation is disk-parallel everywhere and the polarized intensity is confined to the central third of the galaxy, which coincides with the inner part of the hot X-ray halo (Fig. 6).

The RM map shows values between -300 ± 25 rad m⁻² and +300 ± 25 rad m⁻². The sign of RM changes along the disk from east to west, where three overlapping patches are visible. The magnetic field of the patch furthest to the east of the cen-

tral source is pointing away from the observer (negative RM), whereas the magnetic field of the westernmost patch is pointing towards the observer (positive RM). This is a signature of a spiral field that is axisymmetric with respect to the center and symmetric (even parity) with respect to the plane. Assuming trailing spiral arms and considering the rotation of the galaxy with the north-eastern side of the disk being the approaching side
<figure><img src="image_12.png" /><figcaption>**Fig. 11.** NGC 4013 combined (C- and D-configuration) C-band Stokes *I* image (same data shown in Fig. 1) with apparent magnetic field orientations and polarized intensity contours at 3, 6, 9, and 12σ levels with a σ of 2.5 µJy beam⁻¹. The beam size is 5.2” × 5.3” and is shown in the bottom left corner of the image. No tapering was used and the robust parameter was set to two. *Left image:* enlargement of the central region.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. 12.** NGC 4013 combined (C- and D-configuration) C-band polarized intensity from Fig. 11, but smoothed to a beam size of 18” × 18” with a σ of 8.6 µJy beam⁻¹. Apparent magnetic field orientations are shown in light grey, where the polarized intensity is at the 2.5σ level or higher, black apparent magnetic field orientations are shown for polarized intensity of 5σ level or higher. Stokes *I* intensity contours of C-band are displayed in white at 3 and 64σ levels with a σ of 13 µJy beam⁻¹ obtained with an uv-taper of 12 kλ and Gaussian smoothing. The yellow contours of the L-band Stokes *I* intensity obtained with an uv-taper of 18 kλ and Gaussian smoothing are displayed at 3σ level with a σ of 50 µJy beam⁻¹. The orange contours represent the Stokes *I* intensity of LOFAR at 3σ with a σ of 204 µJy beam⁻¹.</figcaption></figure>

(Bottema 1996), the RM values suggest that the radial component of the magnetic field vectors is pointing outwards. Dynamo theory predicts that the radial components of axisymmetric spiral fields should point either inwards or outwards with equal probability, which is indeed the case (with similar numbers of cases found each so far, see Beck 2016).

The field of the central region east of the center shows a field pointing towards the observer in the bottom part (positive RMs) and away from the observer in the upper part (negative RMs). The western side shows hints of the opposite trend, where the upper part points towards the observer and the lower left part points away from the observer. This RM pattern could be caused by an antisymmetric (odd parity) magnetic field structure. A similar pattern, i.e. symmetric in the outer disk and antisymmetric in the inner disk, was proposed for the magnetic field in the Milky Way (Han et al. 1997; Han 2006).

The signal-to-noise ratios of our polarization data are small, so uncertainties remain about the large-scale structure of the

magnetic field. Deeper observations with a wider $\lambda^2$-space distribution are necessary to further investigate the magnetic field structure of this galaxy.

## 3.6. Thermal and nonthermal emission

### 3.6.1. Thermal and nonthermal separation

The separation of thermal and nonthermal (synchrotron) emission in spiral edge-on galaxies was explored by Vargas et al. (2018). We followed their recommended approach of calculating the thermal flux by correcting the $H\alpha$ flux for internal absorption by dust with mid-infrared data, in our case WISE 22 µm data, with the revised correction factor for edge-on galaxies of 0.042 (replacing the 0.031 in Eq. (5) of Calzetti et al. 2007). Nonthermal maps were produced for the CHANG-ES bands with a final resolution of 15” × 15” (see Stein et al. 2019, for details). For the LOFAR data, the expected thermal contribution is small and we disregard it.
<figure><img src="image_14.png" /><figcaption>**Fig. 13.** NGC 4013 combined (C- and D-configuration) C-band polarized intensity contours from Fig. 12 with contours at 2.5, 4, 8 and 16σ levels on SDSS color image (retrieved via Aladin sky atlas, https://aladin.u-strasbg.fr). Regions I, II, and III discussed in the text are marked in yellow.</figcaption></figure>

<table><caption>**Table 10.** RM-synthesis parameters.</caption><thead><tr><th></th><th>C-band</th><th>L-band</th></tr></thead><tbody><tr><td>Bandwidth [GHz]</td><td>2</td><td>0.5</td></tr><tr><td>ν<sub>min</sub> to ν<sub>max</sub> [GHz]</td><td>5 to 7</td><td>1.2 to 1.75 [gap]</td></tr><tr><td>Δλ<sup>2</sup> = λ<sub>1</sub><sup>2</sup> - λ<sub>2</sub><sup>2</sup> [m<sup>2</sup>]</td><td>0.00184</td><td>0.0336</td></tr><tr><td>δφ = <sup>2V<sub>3</sub></sup>&frasl;<sub>Δλ<sup>2</sup></sub> [rad m<sup>-2</sup>]</td><td>1882</td><td>103</td></tr><tr><td>λ<sub>min</sub> [m]</td><td>0.042</td><td>0.17</td></tr><tr><td>max<sub>scale</sub> = π/λ<sub>min</sub><sup>2</sup> [rad m<sup>-2</sup>]</td><td>1781</td><td>109</td></tr><tr><td>spw-width δf [MHz]</td><td>125</td><td>16</td></tr><tr><td>δλ<sup>2</sup> [m<sup>2</sup>]</td><td>0.0001735</td><td>0.001447</td></tr><tr><td>RM<sub>max</sub> = <sup>V<sub>3</sub></sup>&frasl;<sub>δλ<sup>2</sup></sub> [rad m<sup>-2</sup>]</td><td>9983</td><td>1197</td></tr></tbody></table>

**Notes.** Channel width: $\delta\lambda^2$, width of the $\lambda^2$ distribution: $\Delta\lambda^2$, shortest wavelength squared: $\lambda_{min}^2$. FWHM of the resolution (RMSF) in $\Phi$ space: $\delta\phi$ (rad m$^{-2}$). Largest scale in $\Phi$ space to which the observation is sensitive: max-scale (rad m$^{-2}$). Maximum observable RM (resp. $\Phi$): $RM_{max}$ (rad m$^{-2}$) (from Brentjens & de Bruyn 2005).

### 3.6.2. Degree of polarization of the nonthermal emission

To further investigate the appearance of the polarized emission, we plot C-band flux densities of the nonthermal intensity, polarized intensity as well as the fraction of both values within boxes along the major axis. The box width and height were chosen to be equal to the beam size of 15". The result is shown in Fig. 15. The values of the polarized fraction in the disk vary between 4% and 12%, which is a normal range for C-band. Generally, the polarized fraction is highest in the center and decreases with increasing distance along the major axis. The steep increase towards the edges can be explained by the lower rms noise of the polarized intensity map in comparison to the total intensity map, so that we detect polarized emission out to larger distances.

### 3.6.3. Nonthermal fractions

Maps of the nonthermal fraction for both C-band and L-band are shown in Fig. 16. The central point-like source dominates the galaxy and shows synchrotron emission of high intensity as

<figure><img src="image_15.png" /><figcaption>**Fig. 14.** NGC 4013 images from RM-synthesis. Top: total intensity color image with polarization contours at 3, 9, 12σ levels with $\sigma$ of 11.6 µJy. Magnetic field orientations are shown in black. Bottom: RM map of NGC 4013 from C-band on Stokes I contours with a uv-taper of 12 k$\lambda$ and Gaussian smoothing resulting in a beam size of 18" × 18". Contours start at the $3\sigma$ level with a $\sigma$ of 15.0 µJy beam$^{-1}$ and increase in powers of 2 (up to 256). The RM map is cut at the $3\sigma$ level of 35 µJy beam$^{-1}$ of the polarized intensity map. The mean RM error is 25 rad m$^{-2}$.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Fig. 15.** NGC 4013 box integration along the major axis with a box size of 15" × 15", showing the nonthermal intensity (blue data points), the polarized intensity (red data points) as well as the degree of polarization (green) of C-band. The grey box shows the region of the central source, where the data are not considered.</figcaption></figure>

well as a high nonthermal fraction of 93% in C-band. Such a high nonthermal fraction suggests a strong synchrotron source including a high magnetic field strength. The nonthermal fraction of the disk has a mean value of 76%. This implies a thermal fraction of 24% in C-band, which agrees well with the mean thermal fraction of (23 ± 13)% in late-type galaxies at 5 GHz (Tabatabaei et al. 2017).

The nonthermal fraction map in L-band looks similar to the C-band map, but with a more extended distribution towards the halo. The central point-like source is again dominating and
<figure><img src="image_17.png" /><figcaption>**Fig. 16.** NGC 4013 nonthermal fractions of C-band (top) and L-band (bottom). The beam size of 15" × 15" is shown in brown.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>**Fig. 17.** Maps of the nonthermal spectral index. Top: calculated between C- and L-band. Bottom: calculated between L-band and LOFAR. The mean error is 0.05; the error is largest towards the edges with values of up to 0.3. The beam size of 15" × 15" is shown in orange.</figcaption></figure>

exhibits a nonthermal fraction of 97%, which is very high and confirms the C-band results. The entire galaxy disk shows a mean nonthermal fraction of 90% with a corresponding thermal fraction of 10% in L-band. This is in very good agreement with the mean thermal fraction of (10 ± 9)% in late-type galaxies at 1.4 GHz (Tabatabaei et al. 2017).

### 3.6.4. Nonthermal spectral index

As in Stein et al. (2019), we derive the nonthermal spectral index map, here from the nonthermal maps of C-band and L-band as well as from the nonthermal map of L-band and LOFAR, which are presented in Fig. 17. We also created a map of the uncertainties, using error propagation, where we found a mean

<figure><img src="image_19.png" /><figcaption>**Fig. 18.** NGC 4013 magnetic field strength calculated pixel-by-pixel via the equipartition formula using the nonthermal spectral index map and the nonthermal intensity map. The beam size of 15 × 15" is shown in orange.</figcaption></figure>

uncertainty of 0.05. The spectral index map was blanked in places where the uncertainty equals or exceeds 0.3.

The mean spectral index of the disk between C-band and L-band (top of Fig. 17) of $\alpha_{nt} = -0.84 \pm 0.05$ agrees well with the expected theoretical spectral slope of synchrotron emission. The central region shows a slightly flatter spectral index of $\alpha_{nt} = -0.63 \pm 0.03$. At the south-western end of the disk, a region with a flat spectral index of about $\alpha_{nt} = -0.3 \pm 0.2$ is visible. The large error is a result of the weak intensity levels towards the edge of the galaxy. Nevertheless, the flat spectral index indicates a large thermal fraction. This is corroborated by the fact that this region coincides with a star-forming region seen in H$\alpha$ emission. This origin is additionally supported by the lower nonthermal fraction in this region seen in L-band.

The mean spectral index of the disk between L-band and LOFAR (bottom of Fig. 17) of $\alpha_{nt} = -0.69 \pm 0.05$ is flatter in comparison to the higher frequency spectral index. The central source shows a flat spectral index of $\alpha_{nt} = -0.48 \pm 0.03$.

### 3.7. Magnetic field strength via the equipartition assumption

The map of the magnetic field strength was determined using the revised equipartition formula by Beck & Krause (2005). The calculations were done in a pixel-by-pixel fashion by combining the nonthermal intensity map of C-band and the spectral index map between C-band and L-band (top of Fig. 17). As path length, we used 90% of the extent in the radio continuum of 15 kpc. Further, a proton to electron ratio of $K_0 = 100$ is assumed. The uncertainties in the assumptions lead to an error in field strength of about 20% (e.g., Beck 2016).

The resulting map of the magnetic field strength is presented in Fig. 18. The magnetic field strength of the central point source is 12.1 µG. The mean field strength in the disk is weak with 6.6 µG, omitting the central region. In the inner disk, the magnetic field strength is 7.0 µG and the field strength decreases towards the edges of the disk. This map can explain our results from the analysis of the polarized emission. The polarized emission is mostly concentrated in the inner third of the disk along the major axis, where, according to the map shown here, the magnetic field is comparably strong.

### 3.8. Cosmic-ray transport model

Using the 1D CR transport model SPpectral INdex Numerical Analysis of K(c)osmic-ray Electron Radio-emission (SPINNAKER, Heesen et al. 2016, 2018) the transport processes of CRs from the disk into the halo (perpendicular to the galaxy disk) are derived. Analogous to the analysis of NGC 4666
<figure><img src="image_20.png" /><figcaption>**Fig. 19.** Best-fitting CRE advection model with an increasing speed for NGC 4013. Top panel: nonthermal intensities at 150 MHz (LOFAR), middle panel: 1.5 GHz (L-band) data, and bottom panel: nonthermal radio spectral index. In all panels, the solid red line represents the model simulation. The errors are represented by a weighted standard deviation. The best-fitting parameters are listed in Table 11 (SPINNAKER and SPINTERACTIVE are available at: www.github.com/vheesen/Spinnaker).</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Fig. 20.** Best-fitting CRE advection model with an increasing speed for NGC 4013. Top panel: nonthermal intensities at 1.5 GHz (L-band), middle panel: 6 GHz (C-band) data, and bottom panel: nonthermal radio spectral index. In all panels, the solid red line represents the model simulation. The errors are represented by a weighted standard deviation. The best-fitting parameters are listed in Table 11.</figcaption></figure>

(Stein et al. 2019), we carried our an interactive analysis on the intensity profiles of NGC 4013 in order to find the best-fitting model. Because the galaxy is very thin, especially in C-band, only five data points could be measured in the vertical direction (box size of 6'' × 165'') for the L- and C-band data (Fig. 20), and nine for the L-band and LOFAR data (Fig. 19). Since we are using the nonthermal radio continuum maps, we have an angular resolution of 15'' (1.1 kpc).

The models assume either pure diffusion or advection for the CRE, which are injected with a power-law in the disk. The CRE number density is then $N(E) \propto E^{-\gamma}$, where E is the CRE energy. The power-slope can be related to the nonthermal spectral index via $\alpha_{nt} = (1 - \gamma)/2$. Then, a steady state solution is reached between CRE injection and energy losses (synchrotron and inverse Compton losses). The magnetic field strength is described as a double exponential function with the superposition of a thin and thick disk:

$$B(z) = B_1 \exp(-z/h_{B1}) + (B_0 - B_1) \exp(-z/h_{B2}), \quad (6)$$

where $B_0$ is the magnetic field strength in the disk, $B_1$ is the thin disk magnetic field strength, and $h_{B1}$ and $h_{B2}$ are the scale heights of the thin and thick disk, respectively.

To find the best-fitting model, we minimize $\chi_T^2$, which is the quadratic sum of the $\chi^2$-values for the intensity and spectral index profiles. We tested three different models, which are:

(i) Cosmic ray diffusion.

(ii) Cosmic-ray advection with constant speed.

(iii) Cosmic-ray advection with an increasing speed.

The accelerating wind speed is parametrized as follows:

$$V(z) = V_0 \left( 1 + \frac{z}{h_V} \right). \quad (7)$$

For the diffusion model, we assumed no energy dependence for the diffusion coefficient. We fit the three different frequencies

simultaneously pairwise, combining LOFAR and L-band data as well as L-band and C-band data.

In this section, we present the best-fitting profiles for advection with increasing speed [model (iii)] in Figs. 19 and 20, with the other models presented in Appendix A.4. We find that diffusion [model (i)] can describe our L-band and C-band data equally well as advection [models (ii) and (iii)] with each $\chi_T^2$ less than one and thus no possible quantitative statement about which model is the best. While the profile of the radio spectral index between L-band and C-band has a parabolic shape that is typical for diffusion, this can also be explained with advection due to the convolution with the large effective beam ($FWHM = 1.1$ kpc). For advection, we expect the spectral index profile to be linear, but this assumes the vertical profile to be well resolved. For the profile in Fig. 20 the influence of the angular resolution becomes visible, so that, due the Gaussian convolution with an effective beam size of $FWHM = 1.1$ kpc, the resulting spectral index profiles have a shape similar to the one for diffusion.

In contrast, the shape of the spectral index profile between LOFAR and L-band has a linear shape and thus cannot be well fit with diffusion. In Fig. 19 the best-fitting advection model assuming an increasing advection speed is shown, which is preferred over the constant speed since it ensures that energy equipartition is approximately fulfilled in the halo. Nevertheless, the model fits to the LOFAR data are better approximated with diffusion (first row in Fig. A.5) with a smaller $\chi^2$ in comparison to advection. The best-fitting advection speed of $V_0 = 18 - 22$ km s$^{-1}$ (increasing to only ~50 km s$^{-1}$ at 4 kpc height) is very low, not comparable to wind speeds of typically hundreds of kilometers per second. Furthermore, it is lower compared to the outflow velocity estimated in Sect. 3.4, which suggests a basic difference between the global transport mechanisms averaged over the entire galaxy (modeled with SPINNAKER) and local mechanisms that can lead to higher velocities. The best-fitting advection model has $\chi_T^2 = 0.5-1.1$, which shows that the data can be adequately fit with this simple model.
<table><caption>**Table 11.** NGC 4013 SPINNAKER best-fitting parameters.</caption><thead><tr><th>Parameter<br>B<sub>0</sub> (μG)</th><th>6.6 (fixed)</th></tr></thead><tbody><tr><td colspan="2"><strong>Model (i) – Diffusion</strong></td></tr><tr><td>γ</td><td>2.6</td></tr><tr><td>B<sub>1</sub> (μG)</td><td>4.9</td></tr><tr><td>h<sub>B1</sub> (kpc)</td><td>0.1</td></tr><tr><td>h<sub>B2</sub> (kpc)</td><td>9.5</td></tr><tr><td>D (10<sup>28</sup> cm<sup>2</sup> s<sup>-1</sup>)</td><td>0.55–0.65</td></tr><tr><td>χ<sup>2</sup><sub>T</sub></td><td>0.9–2.0</td></tr><tr><td colspan="2"><strong>Model (ii) – Advection (constant speed)</strong></td></tr><tr><td>γ</td><td>2.4</td></tr><tr><td>B<sub>1</sub> (μG)</td><td>4.7</td></tr><tr><td>h<sub>B1</sub> (kpc)</td><td>0.1</td></tr><tr><td>h<sub>B2</sub> (kpc)</td><td>5.5</td></tr><tr><td>V<sub>0</sub> (km s<sup>-1</sup>)</td><td>19–22</td></tr><tr><td>χ<sup>2</sup><sub>T</sub></td><td>0.5–1.2</td></tr><tr><td colspan="2"><strong>Model (iii) – Advection (increasing speed)</strong></td></tr><tr><td>γ</td><td>2.2</td></tr><tr><td>B<sub>1</sub> (μG)</td><td>2.6</td></tr><tr><td>h<sub>B1</sub> (kpc)</td><td>0.4</td></tr><tr><td>h<sub>B2</sub> (kpc)</td><td>15.0</td></tr><tr><td>V<sub>0</sub> (km s<sup>-1</sup>)</td><td>18–22</td></tr><tr><td>h<sub>V</sub> (kpc)</td><td>3</td></tr><tr><td>χ<sup>2</sup><sub>T</sub></td><td>0.5–1.1</td></tr></tbody></table>

# 4. Discussion

## 4.1. Cosmic ray transport

The main cosmic ray transport of this galaxy is analysed via various indicators discussed in this section.

### 4.1.1. Halo scale heights

The frequency dependence of the halo scale height ($\bar{z}_0 \propto v^\gamma$) gives information about the type of CRE propagation, in the case of dominating synchrotron losses (Mulcahy et al. 2018; Krause et al. 2018; Stein et al. 2019). An exponent $x = -0.25$ is expected for diffusive transport and $x = -0.5$ for advective transport. A scenario, where the CREs emit most of their energy within the galaxy halo, is referred to as an *electron calorimeter*. We are working with this assumption (even though it may not be fully valid for the LOFAR 150 MHz data) while keeping in mind that this galaxy has local outflows which are not characterizing the global mechanisms. Using $\bar{z}_0$ for the halo from Table 9 gives $x > -0.38 \pm 0.15$ between 6 GHz and 1.5 GHz. Considering the nonthermal Gaussian scale height in C-band it is $x = -0.20 \pm 0.10$. Thus the smaller C-band/L-band exponent can be attributed to the influence of thermal emission on the scale height at C-band. Together with an exponent of $x = -0.19 \pm 0.05$ between 1.5 GHz and 150 MHz this gives a clear preference for diffusive CRE propagation.

### 4.1.2. Fitting the total intensity distribution

The scale height analysis performed here shows that Gaussian functions fit better to the vertical radio intensity profiles than exponential functions at L-band and the LOFAR frequency. Although the quality of both fits was similar for C-band, the Gaussian profile seems to represent the data slightly better. This suggests that the main transport mechanism of the CREs is diffusion, whereas

<table><caption>**Table 12.** Mean amplitude ratios $\bar{w}(w_0/w_1)$ (Disk/Halo) and mean scale height ratios $\bar{z}(z_0/z_1)$ (Disk/Halo) of the thin disk/halo from Gaussian fits for NGC 4013 and from exponential fits for NGC 4666.</caption><thead><tr><th></th><th>C-band</th><th>L-band</th><th>LOFAR</th></tr></thead><tbody><tr><td colspan="4"><strong>NGC 4013</strong></td></tr><tr><td>\bar{w}(w<sub>0</sub>/w<sub>1</sub>)</td><td>9.40</td><td>4.50</td><td>1.90</td></tr><tr><td>\bar{z}(z<sub>0</sub>/z<sub>1</sub>)</td><td>0.18</td><td>0.18</td><td>0.15</td></tr><tr><td>\bar{w}*\bar{z}</td><td>1.7</td><td>0.81</td><td>0.29</td></tr><tr><td colspan="4"><strong>NGC 4666</strong></td></tr><tr><td>\bar{w}(w<sub>0</sub>/w<sub>1</sub>)</td><td>1.38</td><td>0.97</td><td>—</td></tr><tr><td>\bar{z}(z<sub>0</sub>/z<sub>1</sub>)</td><td>0.26</td><td>0.34</td><td>—</td></tr><tr><td>\bar{w}*\bar{z}</td><td>0.36</td><td>0.33</td><td>—</td></tr></tbody></table>

exponential intensity profiles suggest advection as the main transport process (Heesen et al. 2016). Dominating CRE diffusion implies that no wind or only a slow advection velocity is present in NGC 4013. This results in slow escape of CREs, so that the calorimetric assumption may be used. Therefore, the scale height ratio corroborates the finding of dominating diffusion.

### 4.1.3. Integrated flux densities of disk and halo

The derived ratios of the integrated flux densities between disk and halo (Table 9) decrease from C-band towards lower frequencies with values of 1.7 for C-band, 0.8 for L-band, and 0.3 for the LOFAR frequency. This result shows that the halo becomes more dominant towards lower frequencies. To compare the determined values with another galaxy, we have undertaken the same analysis for the two-component exponential fit of NGC 4666. Stein et al. (2019) derived scale heights for the upper and the lower halo separately for one large strip. The mean amplitudes and scale heights for the two CHANG-ES bands are presented in Table A.2. Interestingly, the ratios of the integrated flux densities for the thin disk and the halo of NGC 4666 are 0.35 in C-band and 0.33 in L-band, much smaller than in NGC 4013 and almost independent of frequency. The halo of NGC 4666 contributes about 74% and 75% of the total flux density in C-band and L-band, respectively, while the corresponding values for the halo of NGC 4013 are about 37% and 56%. The contribution of the halo of NGC 4013 becomes dominant at 150 MHz with 77% of the total flux density. NGC 4666 is transporting the CRs much faster into the halo by an advective wind emerging from the entire galaxy. The halo of NGC 4666 is relatively bright and the disk/halo flux density ratio does not vary with frequency. NGC 4013 is clearly different, with a relatively faint halo and a strong frequency dependence of the disk/halo flux density ratio. A possible reason is diffusive transport of CREs which leads to a more rapid aging due to stronger (energy-dependent) synchrotron losses, explaining why the contribution of the halo decreases with increasing frequency.

As a result we propose to use the multiplication of the mean amplitude ratios $\bar{w}(\bar{w}_0/\bar{w}_1)$ (Disk/Halo) and mean scale height ratios $\bar{z}(\bar{z}_0/\bar{z}_1)$ (Disk/Halo) of different frequencies to get hints on the propagation of cosmic rays in the galaxy. In Table 12 the main values are given. The values of the last row $(\bar{w} * \bar{z})$ show different frequency dependencies for the two galaxies NGC 4013 (diffusion-dominated) and NGC 4666 (wind-dominated).

### 4.1.4. Velocity of cosmic rays

A local transport velocity was determined in Sect. 3.4 via the extent of an emission blob and the synchrotron lifetime,
resulting in a value of 160 km s⁻¹. This velocity is lower than the rotational velocity of the galaxy (see Table 1), whereas other galaxies have advection speeds of (2–3) × vrot (Heesen et al. 2018). Consequently, NGC 4013 hosts only a slow wind that does not reach the escape velocity within a few kpc from the disk. Furthermore, the calculated velocity refers to a local outflow and is therefore different from global advection velocities found via SPINNAKER.

### 4.1.5. Transition from diffusion to advection

The CRE transport is a superposition of diffusion and advection. Our combined C- and L-band data are equally well modeled by diffusion and advection. In contrast, our best-fitting model for the combined L-band and LOFAR data is advection with an accelerating velocity, which fits slightly better than diffusion. Even then, diffusion may play a large role beside advection, as the resulting advection speed of about 20 km s⁻¹ is fairly slow near the disk, but accelerates to about 50 km s⁻¹ at the edge of the observed halo at 4 kpc height. More precisely, diffusion is expected to dominate over advection near the disk, where the CREs are relatively young (observed at L- and C-band), while at larger heights, where the CREs are older (as observed with LOFAR), advection is the dominant transport. The reason for this is that the height that cosmic rays travel scales as z ∝ √Dt for diffusion<fnref n="5" />, whereas for advection the travelling distance is z ∝ Vt. Hence, diffusion dominates over advection for young CRE, whereas for old CRE the opposite is true. Consequently, there is a certain time (CRE age), from which on advection dominates over diffusion.

For a constant wind speed, the transition from the diffusion-dominated regime near the disk to the advection-dominated regime in the halo occurs at a height of (Recchia et al. 2016):

$$s_{\star} \approx (0.3-0.6) \times \frac{D_{28}}{V_{100}} \text{ kpc,} \qquad (8)$$

where D28 is the diffusion coefficient in units of 10²⁸ cm² s⁻¹ and V100 is the advection speed in units of 100 km s⁻¹. With an advection speed of V100 = 0.2 and a diffusion coefficient of D28 = 0.6, we see that for our case s\* ≈ 1–2 kpc. Hence, for the low advection speeds we measure, diffusion would dominate the cosmic-ray transport up to a height of about 2 kpc.

### 4.1.6. X-ray emission

Due to sensitivity of the data and the overall low X-ray brightness of NGC 4013, it was only possible to obtain a reliable spectral model fit to the emission from the galactic core. Nevertheless, since most of the X-ray emission seems to be associated with the inner part of the disk, the derived parameters of the hot gas should be representative for the entire area of the central disk, especially within uncertainties.

The derived temperature of the hot gas in the central region of NGC 4013 of about 3 × 10⁶ K (Table 7) is surprisingly low for the disk of a star-forming galaxy, as such temperatures are more typical for galactic halos (e.g., Strickland et al. 2004; Tüllmann et al. 2006). The low thermal energy density and hence low internal pressure may not be sufficient to drive an efficient outflow. This would confirm the findings from the analysis of the radio data.

### 4.2. Correlation of magnetic field strength with star formation rate

Comparing NGC 4013 to other galaxies, we found NGC 4594 (M104; the Sombrero galaxy) to be quite similar. They both have a comparably low star formation rate (SFR). They both have a dust lane in the optical and have undergone some interaction, with NGC 4013 having a huge warp and tidal stellar streams and NGC 4594 having a dominant nuclear bulge. The radio continuum and polarization intensities of both galaxies are comparable, without signs of an ongoing interaction. The radio continuum observations of NGC 4594 (Bajaja et al. 1988) show a similarly thin disk but with polarization only observable in the inner part of the galaxy rather than along the entire disk, as shown in our RM-synthesis results. In more recent observations by Krause et al. (2006), a plane-parallel disk field, as well as vertical magnetic field components were found, comparable to our findings. In their study the magnetic field strength of NGC 4594 was determined to be 5.0 µG, which is quite small in comparison to other star-forming spiral galaxies (Beck 2016; Krause et al. 2018, for edge-on galaxies). We found here a similarly small magnetic field strength for NGC 4013 of 6.6 µG.

A possible correlation between the SFR and the magnetic field strength is revealed from the low magnetic field strength value in NGC 4013 and NGC 4594 (the Sombrero galaxy) and their low SFR of 0.5 M☉ yr⁻¹ in NGC 4013 (Wiegert et al. 2015) and 0.13 M☉ yr⁻¹ in NGC 4594 (Kharb et al. 2016). Other galaxies following this trend are UGC 10288, with a SFR of 0.4 M☉ yr⁻¹ and a mean magnetic field strength in the disk of 6.3 µG (Krause et al. 2018), and NGC 4565, with a SFR of 0.74 M☉ yr⁻¹ (Wiegert et al. 2015) and a mean magnetic field strength in the disk of 6.0 µG (Schmidt et al. 2019). A correlation between the radio continuum and the SFR is expected, which follows from the well-defined FIR-radio-correlation (e.g., Schleicher & Beck 2013), where the emission is believed to originate in both wavelength regimes from the same population of young stars. To further study this correlation, a larger sample of galaxies should be analyzed and the different components of the magnetic field, namely the turbulent and ordered magnetic field component in galaxies, should be considered (e.g., Drzazga et al. 2011).

The efficiency of the α-effect of the mean-field α-ω dynamo is coupled to the SFR of a galaxy (e.g., Beck et al. 1996). The SFR influences the supernova rate and thus the SN-induced turbulence of the ISM and the magnetic field. The dependence of the magnetic field strength on the SFR is also expected from theoretical simulations (e.g., Bendre et al. 2015). A small-scale dynamo can efficiently amplify the turbulent magnetic field that serves as the seed for the large-scale dynamo (Arshakian et al. 2009). If the SFR is small, the small-scale dynamo is working less effectively.

### 4.3. Central source

The total intensity distribution of NGC 4013 is dominated by the central point source, with 40% of total flux density in C-band, 30% in the L-band, and 13% at LOFAR 150 MHz coming from that region. The radio continuum observations show a central intensity peak at this position. The radio spectral index flattens towards the lower frequencies to α = −0.44 between L-band and the LOFAR frequency. The high resolution L-band B-configuration data of CHANG-ES (including the data of NGC 4013 we use here), are used in Irwin et al. (2019) to search for AGN features in 35 edge-on spiral galaxies, in the radio regime. NGC 4013 does not show any of the characteristic

<footnote n="5">⁵ For 1D anisotropic diffusion, or z ∝ √Dt for 3D isotropic diffusion.</footnote>
features associated with AGNs and was concluded not to be a

smaller towards the galaxy’s center. This could be explained by synchrotron losses as the main energy loss of the CRE particles. The smaller extent of observed radio emission towards the center is due to the larger magnetic fields near the center of the galaxy, which leads to stronger energy losses.

The complementary *XMM-Newton* data show a hard X-ray spectrum and a relatively flat photon index but no clear AGN fingerprint. A dominant X-ray point source coinciding with the location of the radio peak is detected by *Chandra*. These data are suggestive of an AGN.

The fact that we do not see clear signs of an AGN in NGC 4013 could be related to large amounts of dust in the galactic disk (Mosenkov et al. 2018), with the highest concentration in the galactic core. According to Dudik et al. (2005) the central source is possibly a Compton-thick AGN. Due to the high internal absorption in such a case, the obtained X-ray luminosity for the central source of $7.1_{-3.0}^{+5.2} \times 10^{38}$ erg s$^{-1}$ can be considered as a lower limit.

**5. Summary and conclusions**

We present radio continuum CHANG-ES VLA and LOFAR data of the edge-on spiral galaxy NGC 4013 to analyze the total intensity characteristics, the central source and the large scale magnetic field as well as the cosmic ray transport. The total radio continuum intensity distribution of NGC 4013 is dominated by the central point source. The radio and supplemental X-ray data hint at the central source being a Compton-thick AGN as suggested by Dudik et al. (2005). The disk of NGC 4013 shows faint radio emission extending into the halo in C-band with a Gaussian halo scale height of 1.2 kpc. The halo increases in L-band to a scale height of 2.0 kpc and at 150 MHz to 3.1 kpc, so that at low frequencies a large radio halo is observed. The frequency dependence shows a clear preference for diffusive CRE transport. The appearance of a small radio halo in C- and L-band is very different from other spiral galaxies with higher star formation rates, e.g., in NGC 4666 (Stein et al. 2019), where the radio halo at the CHANG-ES frequencies is more extended with a fast advective wind being apparent. It seems that in NGC 4013 fewer or less energetic CREs are reaching the halo, with probably a much weaker wind across the galaxy. This is confirmed by the analysis of the X-ray data, showing a low temperature of the hot gas in the central region of NGC 4013. The resulting lower thermal energy densities of the hot gas could prevent it from leaving the disk plane efficiently.

In this paper, the total flux densities have been derived separately for the thin disk and the halo of a galaxy for the first time. The radio halo of NGC 4013 is relatively faint and contributes only 37% and 56% of the total flux density in C-band and L-band, respectively, compared to 74% and 75% for the bright, wind-driven halo of NGC 4666. The decreasing trend of the ratio of the thin disk to halo with frequency, from C-band to LOFAR 150 MHz, also suggests diffusion to be the main transport process in NGC 4013.

An advective outflow has also been found in the halo of our Milky Way (Everett et al. 2008), which is consistent with the bright radio halo that contributes about 90% of the total flux density at 408 MHz (Beuermann et al. 1985). We propose using either the multiplication of the mean amplitude ratios $\bar{w}(w_0/w_1)$ (Disk/Halo) and mean scale height ratios $\bar{z}(z_0/z_1)$ (Disk/Halo) or the ratio of flux densities between disk and halo at different frequencies as another indicator of the type of CRE transport that is occurring.

Previous observations of Condon (1987) suggested a boxy halo around NGC 4013. This is not found here, but hints of a slightly dumbbell-shaped halo in L-band are seen, where the extent of radio emission perpendicular to the major axis is

Our 1D model suggests that diffusion or a slow and accelerated advective velocity of ~20–50 km s⁻¹ fit equally well as the main CR transport. While a slow advection velocity cannot be ruled out, the cosmic ray transport would be still dominated by diffusion up to 2 kpc height in the halo. Diffusion was found for NGC 7462 by Heesen et al. (2016) with the radio profile being better approximated by a Gaussian than by an exponential fit like in NGC 4013. Scale heights determined by fitting exponential functions are different from Gaussian ones, as a comparison in the three frequencies shows.

Polarized emission in the halo with vertical magnetic field orientations was observed up to about 6 kpc above the plane, indicating an outflow in addition to diffusion. The velocity needed to transport CREs to such height is about 160 km s⁻¹, less than observed in halos driven by a galactic wind. Considering this outflows as a local feature, they are consistent with the result from our transport model for the global transport properties.

In summary, the arguments for a synchrotron-loss dominated CRE halo of NGC 4013 with diffusive propagation are:

- the characteristic frequency dependence of the halo scale heights (Table 9),

- the large disk/halo flux density ratio in C-band, decreasing towards lower frequencies, showing that the halo becomes more dominant over the disk at lower frequencies (Table 9),

- the parabolic shape of the spectral index profile (Fig. 20),

- the low temperature of the X-ray emitting gas (Table 7).

A regular, plane-parallel magnetic field was found to be apparent in the disk of NGC 4013. Additionally, in the low resolution polarized intensity image of C-band (robust 2), halo emission and vertical magnetic field orientations were found. The *XMM-Newton* image shows corresponding emission in the inner part of the galaxy as well as the halo, especially the south and the north-eastern parts show a similar distribution. We propose that the magnetic field has influenced the hot gas dynamics. This is supported by previous findings of other CHANG-ES galaxies like NGC 4631 (Mora-Partiarroyo et al. 2019) or NGC 4666 (Stein et al. 2019), where the magnetic field and X-ray morphologies are also similar.

The mean disk magnetic field strength derived with the equipartition formula of 6.6 µG is rather low. NGC 4013 has low intensity radio continuum and polarization emission and a rather small radio halo. This is probably caused by the low SFR leading to a less effectively working small-scale dynamo. Nevertheless, considering the strong impact of the past interaction(s) on NGC 4013, it is quite remarkable that a large-scale plane-parallel field in the entire galaxy disk as well as a vertical field in the halo is observed in polarization. Our RM data indicate that the large-scale dynamo is still able to generate an axisymmetric spiral field in the disk.

In conclusion, the interaction and the low SFR across the disk probably influence the appearance in the radio continuum of NGC 4013. In the inner third of the galaxy somewhat stronger star formation is apparent, with hints of a star-forming ring, where the polarized emission is strongest and the X-ray emission reaches into the halo.

*Acknowledgements.* We thank the anonymous referee as well as Hans-Rainer Klöckner for constructive and helpful suggestions to improve this work. This research is kindly supported and funded by the Hans-Böckler Foundation. This
research was also supported and funded by the DFG Research Unit 1254 "Magnetisation of Interstellar and Intergalactic Media: The Prospects of Low-Frequency Radio Observations". The LOFAR related research at Ruhr-University Bochum is supported by BMBF Verbundforschung under D-LOFAR IV - FKZ: 05A17PC1. The project has in part also benefit from the exchange programme between Jagiellonian University Kraków and Ruhr-University Bochum. This research has made use of the NASA's Astrophysics Data System Bibliographic Services, and the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration and has made use of "Aladin sky atlas" developed at CDS, Strasbourg Observatory, France.

## References

Arnaud, K. A. 1996, ASP Conf. Ser., 101, 17

Arshakian, T. G., Beck, R., Krause, M., & Sokoloff, D. 2009, A&A, 494, 21

Bajaja, E., Dettmar, R.-J., Hummel, E., & Wielebinski, R. 1988, A&A, 202, 35

Baldi, R. D., Williams, D. R. A., McHardy, I. M., et al. 2018, MNRAS, 476, 3478

Beck, R. 2016, A&ARv, 24, 4

Beck, R., & Krause, M. 2005, Astron. Nachr., 326, 414

Beck, R., Brandenburg, A., Moss, D., Shukurov, A., & Sokoloff, D. 1996, ARA&A, 34, 155

Bendre, A., Gressel, O., & Elstner, D. 2015, Astron. Nachr., 336, 991

Beuermann, K., Kanbach, G., & Berkhuijsen, E. M. 1985, A&A, 153, 17

Bottema, R. 1996, A&A, 306, 345

Brentjens, M. A., & de Bruyn, A. G. 2005, A&A, 441, 1217

Briggs, D. S. 1995, in American Astronomical Society Meeting Abstracts, BAAS, 27, 1444

Calzetti, D., Kennicutt, R. C., Engelbracht, C. W., et al. 2007, ApJ, 666, 870

Carter, J. A., & Read, A. M. 2007, A&A, 464, 1155

Chamandy, L. 2016, MNRAS, 462, 4402

Chyzy, K. T., Jurusik, W., Piotrowska, J., et al. 2018, A&A, 619, A36

Comerón, S., Elmegreen, B. G., Knapen, J. H., et al. 2011, ApJ, 738, L17

Comerón, S., Salo, H., & Knapen, J. H. 2018, A&A, 610, A5

Condon, J. J. 1987, ApJS, 65, 485

Condon, J. J., Cotton, W. D., Greisen, E. W., et al. 1998, AJ, 115, 1693

Damas-Segovia, A., Beck, R., Vollmer, B., et al. 2016, ApJ, 824, 30

Drzazga, R. T., Chyży, K. T., Jurusik, W., & Wiórkiewicz, K. 2011, A&A, 533, A22

Dudik, R. P., Satyapal, S., Gliozzi, M., & Sambruna, R. M. 2005, ApJ, 620, 113

Dumke, M., Krause, M., Wielebinski, R., & Klein, U. 1995, A&A, 302, 691

Everett, J. E., Zweibel, E. G., Benjamin, R. A., et al. 2008, ApJ, 674, 258

Gabriel, C., Denby, M., Fyfe, D. J., et al. 2004, ASPC, 314, 759

García-Burillo, S., Combes, F., & Neri, R. 1999, A&A, 343, 740

Han, J. L. 2006, Chin. J. Astron. Astrophys. Suppl., 6, 211

Han, J. L., Manchester, R. N., Berkhuijsen, E. M., & Beck, R. 1997, A&A, 322, 98

Heesen, V., Dettmar, R.-J., Krause, M., Beck, R., & Stein, Y. 2016, MNRAS, 458, 332

Heesen, V., Krause, M., Beck, R., et al. 2018, MNRAS, 476, 158

Henriksen, R. N. 2017, MNRAS, 469, 4806

Henriksen, R. N., Woodfinden, A., & Irwin, J. A. 2018, MNRAS, 476, 635

Ho, L. C., Shields, J. C., & Filippenko, A. V. 1993, ApJ, 410, 567

Ho, L. C., Filippenko, A. V., & Sargent, W. L. W. 1997, ApJS, 112, 315

Hummel, E., Beck, R., & Dettmar, R.-J. 1991, A&AS, 87, 309

Irwin, J., Beck, R., Benjamin, R. A., et al. 2012, AJ, 144, 43

Irwin, J., Krause, M., English, J., et al. 2013, AJ, 146, 164

Irwin, J. A., Henriksen, R. N., Weżgowiec, M., et al. 2018, MNRAS, 476, 5057

Irwin, J., Wiegert, T., Merritt, A., et al. 2019, AJ, 158, 21

Kaastra, J. S. 1992, An X-Ray Spectral Code for Optically Thin Plasmas, Internal SRON-Leiden Report, updated version 2.0

Kalberla, P. M. W., Burton, W. B., Hartmann, D., et al. 2005, A&A, 440, 775

Kharb, P., Srivastava, S., Singh, V., et al. 2016, MNRAS, 459, 1310

Krause, M., Wielebinski, R., & Dumke, M. 2006, A&A, 448, 133

Krause, M., Irwin, J., Wiegert, T., et al. 2018, A&A, 611, A72

Li, J.-T., & Wang, Q. D. 2013, MNRAS, 428, 2085

Martínez-Delgado, D., Pohlen, M., Gabany, R. J., et al. 2009, ApJ, 692, 955

McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, ASP Conf. Ser., 376, 127

Mewe, R., Gronenschild, E. H. B. M., & van den Oord, G. H. J. 1985, A&AS, 62, 197

Mora-Partiarroyo, S. C., Krause, M., Basu, A., et al. 2019, A&A, 632, A11

Mosenkov, A. V., Allaert, F., Baes, M., et al. 2018, A&A, 616, A120

Mulcahy, D. D., Horneffer, A., Beck, R., et al. 2018, A&A, 615, A98

Müller, P., Krause, M., Beck, R., & Schmidt, P. 2017, A&A, 606, A41

Parker, E. N. 1992, ApJ, 401, 137

Recchia, S., Blasi, P., & Morlino, G. 2016, MNRAS, 462, 4227

Ruzmaikin, A. A., Sokolov, D. D., & Shukurov, A. M. 1988, in Magnetic fields of galaxies, Astrophys. Space Sci. Libr., 133

Schmidt, P., Krause, M., Heesen, V., et al. 2019, A&A, 632, A12

Schleicher, D. R. G., & Beck, R. 2013, A&A, 556, A142

Shimwell, T. W., Röttgering, H. J. A., Best, P. N., et al. 2017, A&A, 598, A104

Shimwell, T. W., Tasse, C., Hardcastle, M. J., et al. 2019, A&A, 622, A1

Stein, Y., Dettmar, R.-J., Irwin, J., et al. 2019, A&A, 623, A33

Strickland, D. K., Heckman, T. M., Colbert, E. J. M., Hoopes, C. G., & Weaver, K. A. 2004, ApJS, 151, 193

Strüder, L., Briel, U., Dennerl, K., et al. 2001, A&A, 365, 18

Tabatabaei, F. S., Schinnerer, E., Krause, M., et al. 2017, ApJ, 836, 185

Tüllmann, R., Pietsch, W., Rossa, J., Breitschwerdt, D., & Dettmar, R.-J. 2006, A&A, 448, 43

Tully, R. B., Shaya, E. J., & Pierce, M. J. 1996, VizieR Online Data Catalog: J/A+A/208/479

Turner, M. J. L., Abbey, A., Arnaud, M., et al. 2001, A&A, 365, 27

Vargas, C. J., Mora-Partiarroyo, S. C., Schmidt, P., et al. 2018, ApJ, 853, 128

Wang, J., Hammer, F., Puech, M., Yang, Y., & Flores, H. 2015, MNRAS, 452, 3551

Wiegert, T., Irwin, J., Miskolczi, A., et al. 2015, AJ, 150, 81

Xu, J., & Han, J.-L. 2014, Res. Astron. Astrophys., 14, 942

Zschaechner, L. K., & Rand, R. J. 2015, ApJ, 808, 153
# Appendix A: Additional tables and figures

## A.1. Imaging values

In Table A.1 we list all achieved rms values for the different images used in this study.

## A.2. Exponential scale height fits of NGC 4013

In Fig. A.1–A.3 we show the exponential fits of the C-band, L-band and LOFAR data.

## A.3. Derived Gaussian flux densities of NGC 4666

With the same approach as for NGC 4013, we derive exponential flux densities of the thin disk and halo for NGC 4666. With the exponential distribution of an intensity profile of a galaxy:

$$w(z) = w_0 \exp(-z/z_0). \qquad (A.1)$$

The integral between zero and infinity of this exponential is.

$$W(\text{mJy} \cdot \text{arcsec / beam}) = w_0 \cdot z_0. \qquad (A.2)$$

<table><caption>**Table A.1.** rms values of NGC 4013.</caption><thead><tr><th>Observation</th><th>Freq (GHz)</th><th>Weighting</th><th>Resolution (")</th><th>rms (μJy beam<sup>-1</sup>)</th></tr></thead><tbody><tr><td rowspan="3">LOFAR</td><td rowspan="3">0.15</td><td>rob 0</td><td>6 × 6</td><td>94.3</td></tr><tr><td>rob 0</td><td>10 × 10</td><td>110</td></tr><tr><td>rob 0</td><td>18 × 18</td><td>204</td></tr><tr><td rowspan="3">CHANG-ES L-band</td><td rowspan="3">1.58</td><td>rob 0, 18 kλ taper</td><td>5.9 × 6.0</td><td>16.7</td></tr><tr><td>rob 0, 18 kλ taper</td><td>10 × 10</td><td>21.9</td></tr><tr><td>rob 0, 18 kλ taper</td><td>18 × 18</td><td>50.0</td></tr><tr><td rowspan="3">CHANG-ES C-band</td><td rowspan="3">6.00</td><td>rob 0, 18 kλ taper</td><td>5.2 × 5.3</td><td>3.6</td></tr><tr><td>rob 0, 18 kλ taper</td><td>10 × 10</td><td>6.8</td></tr><tr><td>rob 0, 18 kλ taper</td><td>18 × 18</td><td>13.0</td></tr><tr><td rowspan="3">Polarization</td><td rowspan="3">6.00</td><td>rob 2</td><td>5.2 × 5.3</td><td>2.5</td></tr><tr><td>rob 2</td><td>15 × 15</td><td>8.0</td></tr><tr><td>rob 2</td><td>18 × 18</td><td>8.6</td></tr></tbody></table>

**Notes.** Images with 10" resolution were used for the scale height analysis, the σ of the polarization was measured with the mean of Stokes Q and Stokes U.

<figure><img src="image_22.png" /><figcaption>**Fig. A.1.** Strip fitting with NOD3 for five strips of NGC 4013 on the combined (C- and D-configuration) C-band data with a two-component exponential fit ($\chi^2_{\text{total}} = 0.03$) for comparison. Designations are the same as for Fig. 8.</figcaption></figure>

The full integral under the exponential curve is twice larger. $w_0$ is measured in mJy beam$^{-1}$ and $z_0$ in arcsec. As the integrated flux density W is in units per beam area, first we divide by the beam area in arcsec ($b = 200.4$ arcsec$^2$). To determine the total flux density F of the thin disk and the halo, we multiply with the extent of the major axis, $e = 145''$ in C-band, $e = 155''$ in L-band, and use again the mean values of the amplitude and scale height:

$$F(\text{mJy}) = 2 \cdot \bar{w}_0 \cdot \bar{z}_0 \cdot e / b. \qquad (A.3)$$

The results are presented in Table A.2 where the entries to the right are taken from Stein et al. (2019) for comparison.

## A.4. CR transport models

In this section, we show the best-fitting diffusion and advective models with a constant wind speed. The diffusion models are shown in Figs. A.4 and A.5 and the advective models in Figs. A.6 and A.7.
<figure><img src="image_23.png" /><figcaption>**Fig. A.2.** Strip fitting with NOD3 for five strips of NGC 4013 on the combined (B-, C- and D-configuration) L-band data with a two-component exponential fit ($\chi^2_{\text{total}} = 1.11$) for comparison. Designations are the same as for Fig. 8.</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>**Fig. A.3.** Strip fitting with NOD3 for five strips of NGC 4013 on the LOFAR data with a two-component exponential fit ($\chi^2_{\text{total}} = 1.86$) for comparison. Designations are the same as for Fig. 8.</figcaption></figure>
<table><caption>**Table A.2.** Integrated flux densities of the thin disk and the halo of NGC 4666.</caption><thead><tr><th rowspan="2"></th><th colspan="5">Exponential fits</th><th colspan="2">From total flux densities</th></tr><tr><th>⟨w̄₀⟩<br>(mJy beam<sup>-1</sup>)</th><th>⟨z̄₀⟩<br>(′′)</th><th>⟨z̄₀⟩<br>(kpc)</th><th>flux density<br>(mJy)</th><th>disk + halo<br>(mJy)</th><th>disk + halo<br>(mJy)</th><th>central source<br>(mJy)</th></tr></thead><tbody><tr><th colspan="8" style="text-align: left;">C-band</th></tr><tr><td>Disk</td><td>6.7 ± 1.5</td><td>3.0 ± 1.2</td><td>0.4 ± 0.1</td><td>29 ± 13</td><td>111 ± 15</td><td>109 ± 1</td><td>1.8 ± 0.2</td></tr><tr><td>Halo</td><td>4.9 ± 0.4</td><td>11.7 ± 0.5</td><td>1.6 ± 0.2</td><td>82 ± 8</td><td></td><td></td><td></td></tr><tr><td>Ratio</td><td>1.38</td><td>0.26</td><td></td><td>0.35</td><td></td><td></td><td></td></tr><tr><th colspan="8" style="text-align: left;">L-band</th></tr><tr><td>Disk</td><td>11.6 ± 0.8</td><td>5.5 ± 0.9</td><td>0.74 ± 0.1</td><td>99 ± 18</td><td>397 ± 44</td><td>398 ± 2</td><td>4.4 ± 0.4</td></tr><tr><td>Halo</td><td>11.9 ± 1.4</td><td>16.2 ± 0.9</td><td>2.1 ± 0.1</td><td>298 ± 40</td><td></td><td></td><td></td></tr><tr><td>Ratio</td><td>0.97</td><td>0.34</td><td></td><td>0.33</td><td></td><td></td><td></td></tr></tbody></table>

Notes. $\bar{w}_0$ being the mean amplitude, $\bar{z}_0$ being the mean scale height.

<figure><img src="image_25.png" /><figcaption>**Fig. A.4.** Best-fitting CRE diffusion model for NGC 4013. Top panel: 150 MHz (LOFAR) data, middle panel: 1.5 GHz (L-band) data, and bottom panel: nonthermal radio spectral index. In all panels, the solid red line represents the model simulation. The errors are represented by a weighted standard deviation. The best-fitting parameters are listed in Table 11.</figcaption></figure>

<figure><img src="image_26.png" /><figcaption>**Fig. A.5.** Best-fitting CRE diffusion model for NGC 4013. Top panel: 1.5 GHz (L-band) data, middle panel: 6 GHz (C-band) data, and bottom panel: nonthermal radio spectral index. In all panels, the solid red line represents the model simulation. The errors are represented by a weighted standard deviation. The best-fitting parameters are listed in Table 11.</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>**Fig. A.6.** Best-fitting CRE advection model with a constant wind speed for NGC 4013. Top panel: 150 MHz (LOFAR) data, middle panel: 1.5 GHz (L-band) data, and bottom panel: nonthermal radio spectral index. In all panels, the solid red line represents the model simulation. The errors are represented by a weighted standard deviation. The best-fitting parameters are listed in Table 11.</figcaption></figure>

<figure><img src="image_28.png" /><figcaption>**Fig. A.7.** Best-fitting CRE advection model with a constant wind speed for NGC 4013. The top panel shows the 1.5 GHz (L-band) data, middle panel: 6 GHz (C-band) data, and bottom panel: nonthermal radio spectral index. In all panels, the solid red line represents the model simulation. The errors are represented by a weighted standard deviation. The best-fitting parameters are listed in Table 11.</figcaption></figure>