© P. K. Humire et al. 2022

1. Introduction

Methanol (CH₃OH) is a molecule prone to population inversion under specific excitation conditions in the interstellar medium (ISM) (e.g., Cragg et al. 1992), causing maser emission. In particular, methanol masers are unique tools for studying physical properties of dense gas associated with young stellar objects (YSOs). Given their brightness and compactness (Menten 1991), their positions can be determined with high precision astrometry (i.e., at milli-arcsecond accuracy with very long baseline interferometry) and over vast distances.

Thousands of such methanol masers have been detected in the Milky Way (Cotton & Yusef-Zadeh 2016; Green et al. 2017; Yang et al. 2019). However, in nearby galaxies we can only account for a handful of successful detections (e.g., McCarthy et al. 2020). In particular, the brightest Galactic methanol maser transition at 6.7 GHz (Breen et al. 2015), remains elusive in extragalactic objects outside the Local Group (Ellingsen et al. 1994; Darling et al. 2003), with the only exception being NGC 3079 (Impellizzeri et al. 2008), and possibly Arp220 (Salter et al. 2008), where this line is detected in absorption.

The early discovery that methanol masers can be divided into two classes, a collisionally pumped Class I and a radiatively pumped Class II (Batrla et al. 1987; Menten 1991), allows us to trace either stellar-induced outflows (Class I) or ultra-compact H II regions (Class II). Class I methanol masers have been observed toward high and low-mass stars (Kalenskii et al. 2006, 2010; Rodríguez-Garza et al. 2017), while Class II masers have been observed only toward high-mass YSOs (Breen et al. 2013). Unlike H₂O and OH masers, Class II methanol masers seem to be exclusively correlated with star-forming regions (Walsh et al. 2001; Breen et al. 2013).

Because Class II masers are usually brighter than Class I masers in our Galaxy, the former type of masers have been studied in great detail, leading to surveys targeting exclusively their relation with the surrounding conditions (Yang et al. 2017; Billington et al. 2020). However, outside our Galaxy, Class II masers were only detected in the Magellanic Clouds (Sinclair et al. 1992; Green et al. 2008; Ellingsen et al. 2010) and the Andromeda galaxy (Sjouwerman et al. 2010), with luminosities not surpassing those in our Galaxy. On the other hand, extragalactic Class I masers can be more luminous than those of Class II, and have been successfully observed beyond the Local Group, particularly in nearby barred spiral galaxies such as NGC 253, IC 342, or NGC 4945 (e.g., Ellingsen et al. 2014; McCarthy et al. 2017; Gorski et al. 2018).

There are two types of methanol. For E-type methanol, one of the protons in the hydrogen atoms of the methyl (CH₃) group has an antiparallel nuclear spin with respect to the others, analogous to the case of para-NH₃. In A-type methanol, the nuclear spins of the three protons in the methyl group are parallel, as in the case of ortho-NH₃. As the two methanol types have different transition frequencies and may arise in different physical environments, we decided to analyze them separately.

Hereafter we use the conventional notation for A⁺ and A⁻ introduced by Lees & Baker (1968), related to a combination between the A–CH₃OH overallparity and Mulliken symbols¹A₁ and A₂. This is done to discriminate between splitted +K and −K levels (doublets), with K being the projection of the angular momentum along the molecular symmetry axis. Furthermore, +K and −K levels are torsionally degenerate for the case of A-type methanol, contrary to the case of E-type methanol (see also Cragg et al. 1993).

As shown by Lees (1973), among E-type methanol transitions (E–CH₃OH) Class I population inversion is favored in the K = −1 relative to the K = 0 or K = −2 ladders. This leads to the prominence of the J₋₁ → (J − 1)₀ − E series (see, e.g., Leurini et al. 2016, their Fig. 3). Indeed, the 4₋₁ → 3₀ − E (36.2 GHz) and 5₋₁ → 4₀ − E (84.5 GHz) lines have been recently discovered to be masing in one extragalactic object, NGC 253 (Ellingsen et al. 2014; McCarthy et al. 2018). For A-type methanol, population inversion in the K = 0 relative to the K = 1 ladder is favored, playing out in the J₀ → (J−1)₁ − A⁺ series, of which the 44.1 GHz line, at J = 7, has been also detected in NGC 253 (Ellingsen et al. 2017). The emission of Galactic Class II methanol masers are more compact than those of their Galactic Class I cousins (Moscadelli et al. 2003; Matsumoto et al. 2014), and that could make Class II methanol masers more difficult to detect at extragalactic distances. In addition, Class I masers require lower densities and temperatures than Class II masers (Menten 2012), making them more numerous. Specifically for the J₋₁ → (J − 1)₀ − E series with J = 4 and 5, their intensities were predicted to be on the order of 50 mJy in the case of NGC 253 (Sobolev 1993).

Class I methanol masers may be associated with a variety of phenomena, such as supernova remnants (Plambeck & Menten 1990; Pihlström et al. 2014), massive protostellar induced outflows (Cyganowski et al. 2018), and interactions of expanding H II regions with surrounding molecular gas (Voronkov et al. 2010), namely regions where shocks compress and heat the gas. In the central molecular zone (CMZ) of our Galaxy (i.e., the inner ∼200 pc in radius; Morris & Serabyn 1996), cosmic ray interactions with molecular clouds have been claimed to be an additional source of methanol production. While a high methanol abundance alone is certainly not sufficient to trigger maser emission, there seems to be indeed a clear enhancement of Class I masers in this region (Yusef-Zadeh et al. 2013; Cotton & Yusef-Zadeh 2016; Ladeyschikov et al. 2019). In particular, extended strong emission showing characteristics of maser action had been found by Haschick & Baan (1993) and Salii et al. (2002) in the G1.6−0.025 region at the periphery of the CMZ and by Szczepanski et al. (1989) and Liechti & Wilson (1996) in a region in which the supernova remnant Sgr A East interacts with a Giant Molecular Cloud (GMC). Noteworthy, CMZ conditions in the Milky Way should provide some guidance to the central regions of starburst galaxies (Belloche et al. 2013).

Class I methanol masers at 36, 44, and 84 GHz have already been reported in NGC 253 by Ellingsen et al. (2014, 2017), and McCarthy et al. (2018), respectively. The former, at 36 GHz, has been also observed in NGC 4945, IC 342, NGC 6946, and Maffei 2 (McCarthy et al. 2017; Gorski et al. 2018; Humire et al. 2020). A 84 GHz mega-maser (≥10⁶ times more luminous than typical Galactic masers; Lo 2005) has been reported in NGC 1068 (Wang et al. 2014), but a confirmation would be needed to put this onto a firm basis. Although a tentative detection of methanol mega-maser emission was claimed in Arp 220 (Chen et al. 2015), later studies ruled it out (Humire et al. 2020; McCarthy et al. 2020). Making use of the unprecedented spectral coverage available by the Atacama Large Millimeter/sub-millimeter Array (ALMA) Large Program ALCHEMI (Martín et al. 2021), we have been able to make a comprehensive study toward one of the best candidates to search for extragalactic maser emission, the Sculptor galaxy NGC 253.

NGC 253 is a nearby (D ∼ 3.5 Mpc, Rekola et al. 2005) highly inclined (i ∼ 70°–79°, Pence 1980; Iodice et al. 2014) SAB(s)c galaxy (de Vaucouleurs et al. 1991) with a systemic heliocentric velocity (v_sys) of ∼258.8 km s⁻¹ (Meyer et al. 2004). Its large-scale bar feeds the nuclear region producing stars at an approximate rate of ∼1.7 M_⊙ yr⁻¹ (Bendo et al. 2015) within the central 20″ × 10″, where 1″ corresponds to ∼17 pc. The position angle (PA) of the large-scale bar is 51° (Pence 1980). The PA continues till the central ∼170 parsecs. Further into the core, the isovelocity contours of the gas change their orientation by about 90° due to the existence of a nuclear bar (Cohen et al. 2020). Like our Galaxy, NGC 253 is characterized by a particularly strong and diverse molecular emission in its central 500 pc (Sakamoto et al. 2006), which we therefore identify as its CMZ. Previous studies do not suggest that an Active Galactic Nucleus (AGN) is important for the properties of the molecular gas in the CMZ of NGC 253 (Müller-Sánchez et al. 2010).

This paper is organized as follows: our observations are described in Sect. 2. In Sect. 3 we show the selected regions to be studied in the CMZ of NGC 253. In Sect. 4 we introduce the two methods used to identify methanol masers. In Sect. 5 we discuss the results and provide a number of conditions to explain maser emission. We finally draw the conclusions of our findings in Sect. 6.

2. Observations

We use ALMA observations of NGC 253 taken in Cycles 5 and 6 as part of the large program ALCHEMI (ALMA Comprehensive High-resolution Extragalactic Molecular Inventory, 2017.1.00161.L, followed up by program 2018.1.00162.S, see Martín et al. 2021). The central region of the galaxy (50″ × 20″ = 850 × 340 pc, with a position angle of 65°) was covered in bands 3–7 (84.2–373.2 GHz), with both the 12 m and 7 m antenna arrays, achieving a final homogeneous angular resolution of $1\stackrel{″}{.}6$. The mentioned region was covered with a single pointing in band 3, with a primary beam ranging between 57″ and 68″, and Nyquist-sampled mosaic patterns of 5 up to 19 pointings for the remaining bands, ensuring a homogeneous sensitivity across the selected region². The phase center of the observations is α = 00^h47^m33^s.26, $\delta =-25°{17}^{\prime }17\stackrel{″}{.}7$ (ICRS³).

The absolute flux calibration accuracy is on the order of 15% for all ALMA bands, since most percentages are lower across individual bands. The flux density RMS noise ranges from 0.18 to 5.0 mJy beam⁻¹, and the averaged sensitivity is 14.8 mK. For more details about the data reduction procedures and image processing, see Martín et al. (2021).

3. Selected positions

In Fig. 1, we show the integrated intensity map of the 6₋₁ → 5₀ − E methanol emission (132.9 GHz rest frequency) using a 3σ clip threshold and applying the formalism suggested by Mangum & Shirley (2015, their Appendix C) to obtain the uncertainties. The selection of this methanol line is motivated by its high signal-to-noise ratio (S/N), the absence of notable blending candidates, and maser emission in some regions (see below). In NGC 253, the inner Lindblad resonance (ILR) was found to be cospatial with its circumnuclear ring (CNR) by Iodice et al. (2014). The total ILR extension measured by them is 0.3 ± 0.1 − 0.4 ± 0.1 kpc (deprojected size). We draw concentric ellipses at 0.2 and 0.5 kpc to denote the ILR limits in Fig. 1 (dash-dotted gray ellipses) as well as its center at 0.35 kpc (solid black ellipse). Along the text, we use the ILR position as a proxy to the x₂ orbit in NGC 253 (see Sect. 5).

Fig. 1. Integrated intensity map for the 6−1 → 50 − E methanol line at 132.9 GHz (within a velocity interval of vsys ≲ ±200 km s−1). Selected positions for spectral extraction (Sect. 3) are encircled in beam-sized apertures and numbered from 1 to 10. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. Superbubbles identified by Sakamoto et al. (2006) are indicated as orange circles. The ALMA beam ($1\stackrel{″}{.}6$) is shown in blue in the bottom left corner.

The presence of methanol masers depends on the physical conditions prevailing across the CMZ of NGC 253. We decided to select ten distinct positions within the CMZ (see Fig. 1 and Table 1). These positions were established based on intensity peaks of the CS J = 2 − 1 and H¹³CN J = 1 − 0 transitions, applying a GMC identification approach based on Leroy et al. (2015), that is to say, using the CPROPS software (Rosolowsky & Leroy 2006). Our region coordinates are not exactly the same as those used by Leroy et al. (2015), even though our numbering is close to theirs. We preferred to use those GMC locations instead of the peak intensities of our methanol lines at a given frequency, because the peak location varies depending on the chosen transition (see also, Zinchenko et al. 2017; McCarthy et al. 2018). We list the coordinates and velocities of all these positions in Table 1. The velocities were obtained from preliminary LTE radiative transfer models, similar to the ones described in Sect. 4.2, and averaged over each methanol symmetric type.

For each region, we extracted the full spectrum (84–373 GHz, see Sect. 4.2) from the data convolved to a common circular beam of $1\stackrel{″}{.}6$ diameter, equivalent to the maximum angular resolution of the observations (see Sect. 2). At higher resolution, however, GMCs further divide into molecular clumps with a size range of $0\stackrel{″}{.}07-0\stackrel{″}{.}25$ (1.2–4.3 pc; Leroy et al. 2018). Contrary to the spectrum of the other nine studied regions, we have found a very complex spectrum in region 5 located next to the dynamical center of the galaxy (Müller-Sánchez et al. 2010). Radio recombination line observations revealed an S-shaped pattern with complex kinematics including a counter-rotating core in the inner 2″ suggestive of a secondary bar (Anantharamaiah & Goss 1996) and evidencing a black hole mass of ∼10⁷ M_⊙ (Cohen et al. 2020). This highly-perturbed environment is likely causing the crowded spectrum observed in region 5. In this position, the broad line emission, with a full width at half maximum (FWHM) above 100 km s⁻¹, prevents an accurate Gaussian line fitting for the rotation diagram analysis in Sect. 4.1, but can be used for the radiative transfer modeling in Sect. 4.2. In addition, toward region 5 and mostly in the low frequency range (∼84–163 GHz, ALMA bands 3 and 4) we observe absorption components.

Methanol lines in absorption or self-absorption (in the case of the line at 358.6 GHz) along the entire ALCHEMI spectral coverage are observed in at least one region for the following four transitions: the 3₁ → 4₀ − A⁺ line at 107.0 GHz, the 4₋₂ → 4₋₀ − E line at 190.1 GHz, the 5₋₂ → 5₋₀ − E line at 190.3 GHz (in absorption in region 10 only), and the 4₁ → 3₀ − E line at 358.6 GHz. They are shown in Fig. 2.

Fig. 2. Methanol lines observed in absorption or self-absorption in one or more of the selected regions (Table 1). Red and blue lines are synthetic spectra for methanol A- and E- type, respectively (see Sect. 4.2). Line frequencies are labeled in the top of the figure, while the velocity range (in km s−1, after applying the radio convention (https://web-archives.iram.fr/IRAMFR/ARN/may95/node4.html) and subtracting the vLSR velocities from Table 1) is indicated in the bottom-left corner. Regions 1–10 are ordered from top to bottom as indicated in the leftmost panel (y-axis: R from “Region” plus the corresponding number). For the location of the regions, see Fig. 1 and Table 1.

4. Methanol maser emission identification

In some⁴ cases, it can be difficult to assess whether a given line is excited under thermal or maser conditions. Methanol maser lines in Galactic star-forming regions are usually very bright (in the sense of brightness temperature) and narrow. However, in the case of extragalactic objects, beam dilution makes them apparently weaker and less distinguishable from thermal emission due to line blending from surrounding gas. In the spectral dimension this also complicates a proper line identification since the contribution of other species might affect the true brightness of the putative maser lines (see the last column of Table A.1).

In this section we explore two methods to identify methanol emission lines outside local thermodynamic equilibrium (LTE) as an indicator of potential maser emission, namely rotation diagram analysis and comparisons with synthetic spectra from radiative transfer modeling. For both methods, we constrained the upper energy above the ground level (E_up/k) of the lines to ≤150 K, since the entirety of lines above that limit are blended with other methanol transitions below that limit, and are likely not contributing significantly to the overall emission (see the analysis in Sect. 4.2).

In this study, we consider as LTE conditions those reached when a single excitation temperature, T_ex, characterizing the energy level population according to the Boltzmann distribution, is sufficient to explain the line emission of methanol along all its observed transitions. Namely, an equivalency between kinetic (T_kin) and excitation temperatures is not strictly required.

4.1. Rotation diagram

As a first approach, we use the rotation diagram method to compare the relative intensities of methanol lines in NGC 253. The details of this method can be found in Goldsmith & Langer (1999) and are also summarized in our Appendix A.

Since this method assumes that the gas is under LTE conditions and that the lines are optically thin, any transition deviating from a straight-line fit in the rotation diagram indicates that, for values below the fitted line, opacities are not low and that, for values above the fitted line, maser emission may be a suitable explanation.

If a given line is affected by blending, we could overestimate its integrated intensity, which would lead to a wrongly classified maser line. Therefore, after identifying all the methanol transitions in our observations and before performing the rotation diagrams, we have initially searched both in the CDMS (Müller et al. 2005) and JPL (Pickett et al. 1998) databases for potentially contaminating lines from other molecular species (last column in Table A.1). Special care was taken on previously reported molecules in NGC 253 (e.g., Martín et al. 2006; Meier et al. 2015; Ando et al. 2017), as well as preliminary line identification performed on the $1\stackrel{″}{.}6$ resolution ALCHEMI data (Martín et al., in prep.). Blended nonmethanol line candidates are required to fall inside the FWHM of the given methanol line, which varies for each symmetric type and region.

To produce the rotation diagrams we have used the CASSIS⁵ software. In particular, we have used the spectroscopic VASTEL database, which comes from the JPL catalog⁶, since it distinguishes between A and E methanol forms. Covering transitions with E_up/k < 150 K, we did not consider methanol lines separated from each other by less than their FWHM. Because of that, we did not include in our analysis the J_K → (J − 1)_K line series at ∼96.7 GHz (J = 2), ∼145.1 GHz (J = 3), ∼193.5 GHz (J = 4), ∼241.8 GHz (J = 5), and ∼290.1 GHz (J = 6).

A summary of the considered transitions for each region is shown in Table A.1, where line parameters were taken from the CDMS database (the same holds for the other tables throughout this work). For each of the selected lines, we have obtained the integrated line intensity ∫T_mbdv by fitting a single Gaussian profile. The resulting rotation diagrams are presented in Fig. 3, for A– and E–CH₃OH symmetry species separately and for each studied region (Table 1), where we conservatively assumed a flux calibration uncertainty of 15%, as recommended by Martín et al. (2021). This uncertainty plus Gaussian fitting uncertainties were added in quadrature for each line.

Fig. 3. Rotation diagrams based on our selected methanol transitions (Table A.1) and separated by methanol symmetric type (E–CH3OH and A–CH3OH). Straight black lines represent our best fits to the data considering only methanol lines not blended with other methanol lines (see Sect. 4.1) and following LTE conditions. Dotted gray lines are the slopes of the complementary methanol symmetric type in the same region. Maser candidates (see Sect. 4.1) are labeled with symbols. Points are color-coded depending on their critical densities (in units of cm−3), as indicated by the color bar to their right (transitions with unavailable critical densities are in black). In the bottom left corner of each panel, the column density (NE − CH3OH and NA − CH3OH, in cm−2) and the excitation temperature (Tex, in Kelvin) derived from the fit are indicated. Non-blended lines in LTE are marked by circles, while methanol lines blended with other methanol lines (at Eup/k > 150 K) are assigned by stars. Shaded areas represent 1σ uncertainties (see Sect. 4.1).

The rotation diagram method assumes the Rayleigh–Jeans (RJ) approximation (see Appendix A). In our case, given the low T_ex measured by the rotation diagrams in certain regions (see Table A.2), the RJ approximation is not a good assumption for high-frequency transitions, leading to an overestimation of ∼40% for the upper level population of some transitions in the worst case scenario. As explained in Appendix A, applying a correction factor to all our transitions, as an attempt to go beyond the RJ approximation, we obtain slight changes in T_ex but total column densities decrease by up to ∼40% in extreme cases. The value that would be derived from the Planck function for the corresponding transition assuming the same T_ex originally derived with the RJ approximation for each region and methanol symmetric type, is shown inside the error bars in Fig. 3. The application of the Planck formula increases the lower-limit uncertainty as a function of the T_ex (see Eq. (A.3)) found for each region through RJ and the individual frequency of each transition and therefore mostly affects the regions where T_ex is lower than 12 K, i.e., regions 1, 2, and 8–10 in E–CH₃OH. The parameters derived after applying this correction are also included in the uncertainties presented in Table A.2.

In Fig. 3 we also indicate the critical densities of the different transitions shown in the rotation diagrams. To this end, we have used the most recently available collisional rates (C_ij), interpolating their values, given in steps of 10 K by Rabli & Flower (2010), by the corresponding T_ex obtained through the rotation diagrams in each region and for each methanol symmetric type. The critical densities thus obtained, in units of cm⁻³ are then indicated by colors in the rotation diagrams (Fig. 3; transitions with unavailable C_ij are in black).

We started including all the considered transitions (Table A.1) below a certain E_up/k (see below) to perform our rotation diagrams and obtain the best linear fit. However, we soon realized that some transitions did not fall on the linear fit line derived from the ensemble of methanol integrated intensities. We excluded these lines that do not follow LTE conditions but are not likely to correspond to masers, see Appendix B, to perform, in a second iteration, a linear fit to the remaining data in order to estimate the column density and T_ex. Our fit results are shown at the bottom of each rotation diagram in Fig. 3 with outliers indicated in the legends when they correspond to maser candidates (see below). Column densities and T_ex with their uncertainties are listed in Table A.2. In general, the distribution of the data points can be well described by a single T_ex.

Without considering the outliers, the standard deviations (1σ) of the residuals of the data points from the fitted lines are of a difference (ln(N_{up,  data points}) – ln(N_{up,  bestfit})) between 0.2 (E–CH₃OH in region 10) and 0.5 (A–CH₃OH in region 4). These values translate into a factor (N_{up,  data points}/N_{up,  bestfit}) between exp(0.2) = 1.2 and exp(0.5) = 1.7.

We conservatively rounded up the average standard deviation (σ = 0.35–0.4) and obtained a factor of 3.3 scatter, in the exponential scale of our rotation diagrams (exp(3σ)). We then classify as outliers all methanol transitions surpassing by more than a factor of 3.3 the expected value from LTE conditions.

These outliers are listed in Table 2, where their transitional quantum numbers, frequencies, symmetry types, maser classes, E_up/k, A_ij, and optical depths are indicated in Cols. (1)–(7), respectively. In this table we also provide information about possible line blending with other methanol transitions falling within the FWHM of the main transition; they possess higher E_up/k and lower A_ij, implying a lower contribution to the observed spectrum due to the need for more extreme conditions to be emitted. In Table 2, the information given for these lines is the transition, frequency, symmetric type, E_up/k, and A_ij along Cols. (8)–(12), respectively.

Deviations from LTE are expected in the ISM, in particular for complex level diagrams and if radiative and collisional excitation and de-excitation compete. But in general, given the success in the fitting, it can be inferred that potential blending lines in Table A.1 are, in most cases, not significantly contaminating.

We wish to emphasize that a sophisticated model, including hundreds or thousands of transitions from many molecules, would be needed to determine methanol line blending with transitions of other species in a thorough way, providing a percentage of contamination to methanol made from other species. Since this is clearly beyond the scope of this paper (but we sometimes quote initial results from region 5, obtained from Martin et al., in prep.), we refer in the following to the case of blending of methanol with other methanol lines, unless contamination with other species is explicitly mentioned.

Based on the line width, CASSIS is able to discriminate between unblended methanol lines and those which are contaminated by adjacent lines of methanol. In Fig. 3 we represent unblended lines with circles and blended lines with stars. Since line widths change per region, the same transition may be denoted by circles or stars depending on the region. Besides of that, we have found that some methanol lines presented in Fig. 3 are blended with other highly energetic methanol lines (E_up/k ≥ 1000 K, see Table 2) and therefore significant blending is highly unlikely. However, we decided to show them also as blended lines (stars) in order to be impartial along all the regions and transitions. In addition, we note that lines with E_up/k ≤ 150 K, being blended with methanol transitions with E_up/k > 150 K, are plotted with the lower E_up/k value in Fig. 3. Finally, the LTE slope is adjusted exclusively considering unblended and nonmasing methanol lines (for a definition of the latter, see below).

4.1.1. Rotation diagram results

To distinguish maser candidates among outliers, that is, lines located beyond a 3σ scatter; the respective line must belong to known methanol maser line series already discovered in NGC 253 at lower frequencies. This is because maser action is more prominent at lower frequencies than the ones covered in this work: as A_ij is proportional to the cube of the frequency, lower A_ij leads to a longer time lapse to accumulate inverted populations. In addition, we also require that the candidate line shows maser behavior in more than one region.

In the following we measure the departure from LTE of the maser lines as the ratio between their nominal upper level column densities in the rotation diagrams (N_up, maser) over the expected upper level column density in LTE (N_up, LTE). In the computation of N_up one generally assumes a proportional relation with the integrated intensity (see Eq. (A.2)). For maser lines, however, intensities and abundances are not related, as their negative opacity amplifies radiation from the background. Therefore this N_up, maser/N_up, LTE factor must be taken as an intensity difference instead of abundance difference.

4.1.2. E-type methanol masers

Among the E–CH₃OH maser line candidates belonging to the J₋₁ → (J − 1)₀ − E series (at 84.5, 132.9, and 229.8 GHz), none of them were detected out of LTE in regions 3–6 (see Table 2), which are therefore excluded from the following analysis.

We present the spectra and velocity integrated intensity maps of the two transitions with the highest S/Ns (5₋₁ → 4₀ − E at 84.5 GHz and 6₋₁ → 5₀ − E at 132.9 GHz) in Figs. 4 and 5.

Fig. 4. Central panel: velocity integrated intensity of the 5−1 → 40 − E methanol line at 84.5 GHz obtained by integrating the channels inside the colored areas shown in the side panels. Regions where we detect Class I maser emission in the J−1 → (J − 1)0 − E series (see Sect. 4.1.2) are labeled in red and with slightly larger numbers. Side panels: methanol spectra of the different regions as labeled in the top right corner of each panel (a vsys of 258.8 km s−1 was subtracted). Vertical lines indicate the peak velocity of the 13CO J = 1 − 0 line (110.2 GHz) extracted in each region, as a reference for the systemic velocity at that position. Coordinate system, beam size and other parameters are the same as in Fig. 1.

Fig. 5. Same as for Fig. 4, but for the 6−1 → 50 − E methanol line at 132.9 GHz.

The 5₋₁ → 4₀ − E line at 84.5 GHz departs from LTE by factors (N_up, maser/N_up, LTE) ranging from 4.2 ± 1.3 to 13.2 ± 1.2 in regions 10 and 1, respectively (see Fig. 3). This is the first detection of the 84.5 GHz maser line in more than one position, after its first discovery by McCarthy et al. (2018). The higher J masers in the J₋₁ → (J − 1)₀ − E series are usually observed to originate in the same regions, suggesting similar distributions and pumping mechanism.

The 6₋₁ → 5₀ − E line at 132.9 GHz departs from its expected integrated intensity in LTE by factors ranging from 4.6 ± 1.2 to 10.6 ± 1.2, in regions 2 and 1, respectively. Contrary to the previous transition at 84.5 GHz, it is not blended with other methanol lines.

The 8₋₁ → 7₀ − E line at 229.8 GHz is observed to mase in the same regions as the ones at 84.5 and 132.9 GHz except in regions 2 and 10, which display the lowest LTE fitted T_ex in the E-type methanol species and also the lowest S/Ns. In these couple of regions, the emission is too weak to be classified as a maser. The 8₋₁ → 7₀ − E line is the highest frequency transition detected as a maser by the rotation diagram method; its intensity departs from LTE by 9.8 ± 1.2, in region 7 and up to 32.0 ± 1.2 in region 1.

4.1.3. A-type methanol masers

With respect to A–CH₃OH, maser line candidates are part of the J₀ → (J − 1)₁ − A⁺ series (9₀ → 8₁ − A⁺ and 10₀ → 9₁ − A⁺ transitions). We note that a lower J (and lower frequency) transition in this series is the 7₀ → 6₁ − A⁺ line at 44.1 GHz, which is the strongest Galactic Class I maser.

The 9₀ → 8₁ − A⁺ line at 146.618 GHz departs from LTE by factors in the range of 8.4 ± 1.2 (in region 8) to 29.8 ± 1.2 (in region 1), excluding region 4 because of line blending (see Appendix B). This maser transition has been suggested to be part of the Class I family of methanol masers (e.g., Yang et al. 2020, and references therein).

The 10₀ → 9₁ − A⁺ line at 198.4 GHz departs from LTE in regions 3, 6, and 7, being the only masing line candidate for regions 3 and 6. In those regions, this line departs from LTE by factors of 4.3 ± 1.3, 3.6 ± 1.2, and 8.5 ± 1.2, respectively.

We have found that all our maser candidates, except the 10₀ → 9₁ − A⁺ line at 198.4 GHz, are masers at the outskirts of the CMZ of NGC 253 (see Table 2). The physical conditions in the regions giving rise to those masers are characterized (see Fig. 6) by total column densities for A–CH₃OH lower than 6 × 10¹⁴ cm⁻², an excitation temperature lower than 15 K in E–CH₃OH, and differences between A– and E–CH₃OH excitation temperatures larger than 5.0 K. Conversely, we do not see a correlation between the presence of masers and total column densities in E-type methanol, the temperature described by A–CH₃OH, or differences between column densities of the two methanol types.

Fig. 6. Physical conditions for the methanol maser candidates (see Table 2) derived from our rotation diagrams. The x-axis indicates the number of the region. The y-axis represents, from top to bottom, (1) the ratio between the nominal upper level column densities of the lines Nup, MASER and the expected values from our best fit to LTE conditions Nup,LTE, (2) the total column density, (3) the difference between total E– and A–type column densities, (4) the LTE excitation temperatures and (5) the difference between the A– and E–type excitation temperatures. Shaded areas highlight regions were maser emission in the J−1 → (J − 1)0 − E line series (at 84, 132, and 229 GHz), and the J0 → (J − 1)1 − A+ line series (at 146 and 198.4 GHz), are observed, according to the first panel.

With the exception of the line at 198.4 GHz, it is worth noting that all the regions where we detect Class I maser emission lie inside the Lindblad resonances (see, e.g., Fig. 4), with the only exception being region 7. We further explore this point in Sect. 5.2.

Based on LTE modeling, we have also found other transitions that are possibly experiencing maser activity. They will be described below.

4.2. Radiative transfer modeling

In comparison with the rotation diagram method, synthetic spectra offer a number of improvements, such as the possibility to reproduce line profiles both in LTE and out of LTE, which may also subtly indicate the number of gas components present in the observed source (see, e.g., Fig. C.1). Also, as the lines are observed on a linear scale, the effect of slight changes in the assumed gas conditions are noticed in more detail compared with the logarithmic scale used in the rotation diagram method. Apart from those advantages, radiative transfer modeling allows us to fit not only single methanol lines, but also blended ones, as their spectra are summed up when trying to obtain the desired intensities. Alternatively, this provides an estimate for the degree of contamination caused by other species. This allows us to include many more lines, increasing the sample size. In fact, as we see below, there are 57 lines used for the fitting in our radiative transfer models, in contrast to the 39 finally used in the rotation diagrams (see Sect. 4.1.1).

The possibility to perform non-LTE modeling helps us to discard maser line candidates in case negative optical depths are not required to reproduce their emission. Finally, the nature of blended methanol lines with other methanol lines at known maser frequencies, such as the transitions at 95.169 and 278.305 GHz, could not be unveiled through the rotation diagram method. As mentioned in Sect. 4.1, we use the term blending to refer to contamination of methanol lines with other methanol transitions, unless contamination by other species is explicitly mentioned.

We have computed synthetic spectra for each of the 10 selected regions in NGC 253, covering the entire ALCHEMI frequency range (∼84–373 GHz, see Sect. 2), by using the CASSIS software, capable of producing LTE and non-LTE spectral modeling. For the non-LTE case, CASSIS is used as a wrapper of RADEX (van der Tak et al. 2007), a one-dimensional non-LTE radiative transfer code based on the escape probability formulation. In CASSIS it is possible to create a physical model defined by 6 parameters: column density of the species (N(Sp)), excitation (T_ex, LTE case) or kinetic (T_kin, non-LTE case) temperature, full width at half maximum (FWHM) of the lines, velocity of the source in the Local Standard of Rest (V_LSR) system, size of the source in arcseconds, and the H₂ volumetric density (n_H2) in the case of non-LTE modeling. For simplicity, in this section we have considered that only a single physical component for each methanol type is responsible for the observed emission; for a more detailed analysis, see Appendix C. CASSIS makes use of the Monte Carlo Markov chain (MCMC) method (Hastings 1970) to explore a user-predefined range of values for each of the parameters previously mentioned. By means of the χ² minimization method, CASSIS is able to find the best ensemble of solutions. Computing the synthetic spectra we assumed a beam filling factor of one, by selecting a source size of $1\stackrel{″}{.}6$ (∼27 pc), and a slab geometry (appropriate for shocks, e.g., Leurini et al. 2016).

Model solutions could be strongly influenced by the initial conditions. Therefore we start by selecting unblended methanol lines to be fitted and then add blended lines whose total line profile is successfully reproduced with the starting models. This was achieved by exploring the success of CASSIS in reproducing the lines throughout our selected regions.

After our initial attempts to fit methanol lines along the entire ALCHEMI spectral coverage, assuming either LTE or non-LTE conditions, it became clear from our 10 regions that it is impossible to properly fit all the lines with a single physical set of parameters, even after separating between A– and E–CH₃OH flavors. This is especially important at frequencies below ∼156 GHz, and it is possibly due to a couple of factors: the presence of a series of lines out of LTE, and the higher number of maser candidates.

Using our LTE model, we found that the J_K → (J − 1)_K transitions (for both methanol symmetric types), which were avoided in Sect. 4.1 due to line blending, are not following LTE conditions. Within the ALCHEMI frequency coverage, these line series have the following frequencies: ∼96.7 GHz (J = 2), ∼145.1 GHz (J = 3), ∼193.5 GHz (J = 4), ∼241.8 GHz (J = 5), and ∼290.1 GHz (J = 6). Performing a non-LTE model in region 8, they are satisfactorily reproduced, as can be seen in Fig. C.2 without invoking the presence of masers.

At frequencies below ∼156 GHz we also cover four maser candidates, three of them reported in Sect. 4.1.2 plus an additional one that we see below (Sect. 4.2.1). Therefore, at lower frequencies our models fail to reproduce an important quantity of available, not blended, methanol lines.

We discard possible software issues by doing a sanity check with another radiative-transfer code capable of producing LTE models, MADCUBA (Martín et al. 2019). It shows similar results including convergence primarily toward lines above 156 GHz.

Having said the above, we do not include the mentioned J_K → (J − 1)_K transition series nor the maser candidates in our LTE modeling. Even with this restriction, we observe that these series are well reproduced (> 50%) in region 4 by merely fitting the other transitions, and this can be improved with a two-component LTE model (see Appendix C).

When inspecting and comparing the results for all regions between the LTE and non-LTE models, considering a single component for each methanol symmetric type, we observe that they are in agreement within an uncertainty of about 15% (although this agreement is not observed when we consider a two-component non-LTE model, see Appendix C). This is not surprising since high H₂ densities of > 10⁷ cm⁻³ are needed in the non-LTE models to reproduce most line profiles. At such densities collisions play a major role, tending to constrain the spread in excitation temperatures between different lines. Therefore, we conclude that the simpler LTE conditions are sufficient to represent the observed spectra in NGC 253. From now on, we therefore mainly refer to our LTE modeling, except in a few exceptional cases where this is explicitly mentioned.

As described in Sect. 4, we set an upper E_up/k threshold of 150 K for the synthetic spectra. This was determined through comparing observations to model fits with and without the higher energy levels (E_up/k > 150 K). A total of 600 LTE models for each individual spectrum, one per selected region (see Table 3), were computed, reaching a convergence after 300–400 iterations. We find that we only require models that include E_up/k < 150 K to fit the observed spectra. Thermal line emission with E_up/k > 150 K may be there, but is too faint to affect model fits or to be separated from line blends. In the Galaxy, methanol masers with levels around 150 K above the ground state are scarce. Within the ALCHEMI frequency coverage we can mention the case of the 11₋₁ → 10₋₂ − E line at 104.3 GHz (Leurini & Menten 2018). However, this Class I maser is rarely seen (Voronkov et al. 2012). We do not obtain any strong emission at this frequency in NGC 253.

Best fit parameters determined from our LTE models for each region and methanol species are listed in Table 3. Differences between A– and E–CH₃OH symmetry species are present in terms of temperature and density. In velocity they differ by a few km s⁻¹, although this discrepancy never exceeds our spectral resolution of ∼8–9 km s⁻¹. Thus, the averaged V_LSR of these two species is adopted as the velocity of the region, in the same way as it has been previously established in Table 1.

A comparison between synthetic and observed spectra allows us to perform a deep scan of the methanol lines, highlighting some lines that slightly deviate from LTE conditions in the rotation diagrams or that are blended with other methanol lines and that were previously not discussed. The result of this inspection is henceforth described.

4.2.1. LTE modeling results

As can be seen in our Fig. 7, our synthetic spectra fit reasonably well the observed spectra in most of the regions. Region 5 is the most difficult to reproduce due to the large FWHM (∼140 km s⁻¹) of the lines, with its spectrum almost reaching the confusion limit. Fortunately, based on the remaining regions, we selected a large number of lines that are reproduced (see Tables C.1 and C.2, and Fig. 7). We attempted to fit the same lines in region 5, preventing in this way false line identifications and subsequent erroneous fitting attempts.

Fig. 7. Zoom into successfully fitted lines for each region, which can be considered to follow LTE conditions. Red and blue colors indicate methanol A- and E-type, respectively, while the superposition of the two types, the final fit, is indicated in green. Line frequencies (in GHz units) are labeled at the top of the figure and are also indicated as a green dashed vertical line inside panels, while the velocity range (in km s−1 with respect to the systemic velocity of the individual regions, see Fig. 4) is indicated in the bottom-left corner. Regions 1–10 are ordered from top to bottom as indicated in the leftmost panel (R from “Region” plus the corresponding number).

The best solution for Region 5 was reached by fitting both methanol symmetric types simultaneously. When this is done in CASSIS, the whole set of parameters between A- and E-type methanol is forced to be equal, allowing to change only the ratio between their column densities. As can be seen in Table 3, the column density of E–CH₃OH is 2.1 times higher than that of A–CH₃OH. A lower column density E/A ratio leads us always to a worse fitting and is therefore not attempted to reach. We have also discarded from the initial fit the A–CH₃OH transitions between 303.3 and 309.2 GHz, as their inclusion always leads to an overestimation of a number of lines (e.g., transitions at 239.7, 241.9, 338.6, and 350.9 GHz, see Tables C.1 and C.2).

All the proposed maser candidates listed in Table 2 that were initially unveiled through the rotation diagrams were confirmed by our models.

In addition to the outliers previously detected (see Table 2), we found several transitions with intensities not reproduced by our synthetic spectra. Among those outliers, listed in Table 4, there are a couple of maser candidates. We describe them below.

4.2.2. Maser line candidates

Contrary to the remaining outliers found through the radiative transfer modeling (shown in Fig. 8 and described at some level in Appendix D), our maser candidates have negative optical depths, depart significantly from LTE (see below), and belong to the same transition series as the maser candidates detected previously through the rotation diagram method (Sect. 4.1.3), being Class I methanol masers.

Fig. 8. Outliers to our LTE modeling (see Appendix D) with positive optical depths. The only exceptions are the J1 → J0A−+ line series, recently proposed to present methanol maser Class II activity by Zinchenko et al. (2017). Line frequencies are labeled at the top of the Figure. Most of the lines are unclassified so far. Labels and colors are the same as for Fig. 7.

The 8₀ → 7₁ − A⁺ transition line at 95.2 GHz shows intensities > 12.8 times stronger than predicted by the LTE modeling in all the regions but region 4, where it is three times stronger. This large departure in the inner regions of the CMZ remind us of the case of the 10₀ → 9₁ − A⁺ methanol line at 198.4 GHz, which shows emission 3.3 times larger than expected in regions 3 and 6.

The 9₋₁ → 8₀ − E line at 278.305 GHz is the last methanol transition in the J₋₁ → (J − 1)₀ − E line series covered in this work. As previously mentioned in Sect. 4.1.2, it is strongly contaminated by the 2₋₂ → 3₋₁ − E line at 278.342 GHz and was therefore excluded from the analysis with rotation diagrams. In regions 1, 7, 8, and 9, however, its line profile is clearly distinguished from the contaminating line when being checked by the synthetic models. According to our radiative transfer modeling, the companion line at 278.342 GHz is the only one that should be observed (while the masing line intensity should be negligible under LTE, see Fig. 9), with its peak velocity about 40 km s⁻¹ lower than the maser line in region 9, a difference four times larger than our spectral resolution.

Fig. 9. Our proposed methanol masers belonging to previously known maser transitions. Labels and colors are the same that for Fig. 7.

In summary, all the outliers in the LTE models have intensities above the expected one under LTE conditions, but most of them are likely able to be reproduced under non-LTE conditions if the results obtained for region 8 (Appendix C) are maintained for the other regions. The only clear maser candidates are the lines at 95.2 and 278.3 GHz.

Given the spatial resolution and software limitations (see van der Tak et al. 2007, their Sect. 3.6), we cannot account for a proper model to our maser candidates. Instead, we reproduce their line profiles and check whether negative optical depths are required or not. When no other results are plausible, we can assure that the transitions are effectively experiencing a population inversion. Unfortunately, as Class I maser emission arises from spot sizes on the order of 10⁻⁵–10⁻³ pc (Voronkov et al. 2014; Matsumoto et al. 2014), any attempt to model our observations will be devoid of a real physical meaning. Negative optical depths are determined in the lines belonging to the J₋₁ → (J − 1)₀ − E and J₀ → (J − 1)₁ − A⁺ series (see Fig. 9) and they constitute our maser candidates. All of them are Class I methanol masers and, for the case of those belonging to the J₋₁ → (J − 1)₀ − E line series, depart from LTE at the outskirts of the CMZ of NGC 253.

We summarize the observed methanol transitions along the entire ALCHEMI coverage in Fig. 10, where LTE lines are labeled in straight lines and maser candidates are in dashed lines. We indicate with colors the ALMA band for each transition. From Fig. 10 it becomes clear that the 7₋₁ → 6₀ − E line at 181.295 GHz should be pumped by the same conditions than the other maser candidates in the J₋₁ → (J − 1)₀ − E series. Our non-LTE models performed in region 8 yield a negative optical depth (−2.4) for this transition suggesting maser emission. Unfortunately, our spectral and angular resolution is not sufficient to discriminate between the 181.295 GHz line and the J = 2 − 1 HNC line at 181.324 GHz.

Fig. 10. Diagram of the energy levels for A-type (left) and E-type (right) methanol transitions covered by ALCHEMI. The x-axis gives the quantum number K, while the number for each level gives J. The transitions used in the analysis are indicated with a color code corresponding to ALMA bands, as given in the top right corner. The maser candidate lines are presented as dashed lines, with their frequency in GHz. We note the energy difference of 7.9 K between the ground state of A and E-type.

5. Discussion

5.1. Conditions for maser emission

Among all the lines proposed to be masers, the ones detected initially by our rotation diagrams are the most plausible ones. Here we focus on the Class I methanol masers in the J₋₁ → (J − 1)₀ − E series (at 84, 132, 229, and 278 GHz), as they describe a very clear difference in maser occurrence: a complete LTE behavior in the central regions (3–6) and strong maser activity at the outskirts of the CMZ of NGC 253, where important T_ex differences between A- and E-type methanol take place (see Fig. 6), and also where Lindblad resonances are located (with the only exception of region 7).

Unfortunately, based on the rotation diagram method we do not have enough evidence to account for the presence of Class II masers. The only exception might be the 3₁ → 4₀ − A⁺ transition line at 107 GHz (E_up/k = 28 K), that follows LTE conditions in regions 2, 5, and 10 (∼1.5, ∼1.5, and 2σ, respectively). It is slightly above the fit in regions 3, 4, 6, and 7, and is present in absorption in regions 1 and 9 (3 and 2σ, respectively), and maybe also in region 8 (2σ); although in this region the line seems to exhibit emission in the middle of the absorption feature.

Maser activity by the line at 107 GHz was first discovered by Val’tts et al. (1995). They found this line either in absorption or showing quasi-thermal emission. In the first case, this line is spatially correlated with the 5₁ → 6₀ − A⁺ Class II maser line at 6.7 GHz, that belongs to the same family of lines. Other lines in the J₁ → (J + 1)₀ − A⁺ series are the transitions at 57 GHz (J = 4), 156.6 GHz (J = 2), and 205.8 GHz (J = 1). By means of the rotation diagram method, we note that the transition at 156.6 GHz is above the LTE trend in the same regions where the line at 107 GHz surpasses the LTE modeling (regions 3, 4, 6, and 7). However, as indicated in Appendix B, the line at 156.6 GHz may also follow the LTE conditions of A–CH₃OH.

In NGC 253, absorption against the continuum in region 5 is noticeably seen in the rotational ground state transitions of dense gas tracers such as the formyl cation (H¹³CO⁺; Harada et al. 2021). Absorption features are also very prominent in the J = 2 − 1 (e.g., Meier et al. 2015) and J = 3 − 2 transitions of silicon monoxide (SiO). All those features are uniquely observed in region 5, contrary to the case of the 3₁ → 4₀ − A⁺ transition line at 107 GHz, that presents absorption mostly in regions 1 and 9.

Therefore, the edges of the CMZ, namely regions 1, 2, 7, 8, 9, and 10, appear to provide suitable conditions for Class I maser emission in the J₋₁ → (J − 1)₀ − E line series. We are aware that these conclusions only refer to the average conditions in our selected regions, as our linear resolution is on the order of 27 pc, larger than typical clump sizes by a factor of 6–23 (see, e.g., Leroy et al. 2018, for estimations of clump sizes in the CMZ of NGC 253) but usually small enough to resolve GMCs.

Although Lindblad resonances (e.g., ellipses in Fig. 11) have been claimed as a good candidate for shocks (García-Burillo et al. 2000; Ellingsen 2018), for galaxies with a strong bar there is theoretical support against the existence of Lindblad resonances (Regan & Teuben 2003) or their relationship with the circumnuclear ring (CNR) position, which is more accurately defined by x₂ orbits (Kim et al. 2012; Li et al. 2015; Schmidt et al. 2019). The interplay between the nuclear dust/gas lanes and the CNR looks as a more likely explanation for the production of bi-symmetric shock/active regions: at the outskirts of the CMZ, where we find higher levels of HNCO as compared to SiO (regions 1, 2, and 7–10; see Fig. 11), in a similar way as in M 83 (Harada et al. 2019), leading to the appearance of methanol masers in the J₋₁ → (J − 1)₀ − E line series. We should see maser action encircling the CMZ if they were caused by the ILR.

Fig. 11. Strong over weak shock tracers as accounted by SiO J = 2 − 1/HNCO (40, 4–30, 3) (left) and SiO J = 5 − 4/HNCO (100, 10–90, 9) (right) line ratios obtained from ALCHEMI data (Martín et al. 2021). A sigma clip of 3.0 was applied for both SiO and HNCO before obtaining the line ratios. Regions where we detect Class I maser emission in the J−1 → (J − 1)0 − E series (see Sect. 4.1.2) are labeled in red and with slightly larger numbers. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. A square root stretch has been applied for an easy visualization.

As an attempt to account for gas disturbances in the CMZ of NGC 253, we obtained line intensity ratios from two well-known shock tracers, SiO and HNCO. Their moment 0 maps, along with data reduction and a thorough analysis will be described in a forthcoming paper by Huang et al. (in prep.). Here our only interest is to relate the presence of shocks to methanol maser emission in the J₋₁ → (J − 1)₀ − E line series, which is mostly observed toward intersecting areas between the circumnuclear ring and leading-edge dust/gas lanes of the bar.

We selected transitions in the same ALMA bands (3 and 6) in order to reduce instrumental uncertainties. Based on previous studies, here we consider SiO as a tracer of strong shocks and HNCO as a tracer of weak shocks (e.g., Meier et al. 2015; Yu et al. 2018), and therefore the SiO/HNCO ratio may be used as an indicator of the shock strength (Kelly et al. 2017).

We note that HNCO (4_0, 4–3_0, 3), at 87.9 GHz, is enhanced in practically the same regions where we detect Class I maser activity in the J₋₁ → (J − 1)₀ − E line series (at 84, 132, 229, and 278 GHz); this is also true for the previously detected methanol transition in this series, at 36 GHz, if we dismiss the weak maser candidate at the very center shown in Gorski et al. (2018). The resulting SiO J = 2 − 1/HNCO (4_0, 4–3_0, 3) ratios (see Fig. 11, left panel) are enhanced in the central regions, especially in regions 4 and 6, where no maser activity, in the J₋₁ → (J − 1)₀ − E line series, is detected.

Considering the SiO and HNCO transitions at frequencies ∼220 GHz (ALMA band 6), the largest SiO J = 5 − 4/HNCO (10_0, 10–9_0, 9) line ratios come from region 1 (see Fig. 11, right panel), the region where we see the largest departure between maser Class I emission in the J₋₁ → (J − 1)₀ − E line series (and also in the A-type methanol line at 146.6 GHz) and LTE conditions (see Fig. 6). On the other hand, the lowest SiO J = 5 − 4/HNCO (10_0, 10–9_0, 9) ratios are observed in region 3 and this coincides with the lowest LTE departures of potential Class I masers in E–CH₃OH lines (see Fig. 6).

One, however, also needs to take into consideration the varying gas properties likely causing the highly varying trend of these SiO/HNCO ratios. A more complete investigation and further discussion upon the gas properties traced by HNCO and SiO is covered in a forthcoming paper (Huang et al., in prep.). In addition, there is something missing in describing the intermediate regions. For example, regions 4 and 6, where we do not see Class I maser emission in E–CH₃OH, present SiO J = 5 − 4/HNCO (10_0, 10–9_0, 9) ratios similar to regions 1, 2, 7, 8, and 9, where we do observe Class I maser emission in E–CH₃OH.

The missing element to be considered for the presence or absence of methanol masers can be the occurrence of photodissociation regions (PDRs), since photodissociation is the main mechanism for methanol destruction (see, e.g., Hartquist et al. 1995). We can use PDR tracers such as CN, whose abundance is found to be high in PDR regions (e.g., Fuente et al. 1993; Kim et al. 2020), to investigate whether PDRs may be dominant. As mentioned by Meier et al. (2015), the high rate of star production in the central regions of NGC 253 can be identified, for instance, by means of CN/C¹⁷O line ratios, where C¹⁷O traces dense molecular gas (e.g., Thomas & Fuller 2008).

We obtained the CN/C¹⁷O line ratios by taking the ratio of moment 0 maps of the chosen molecules. A 3σ clipping was applied in the creation of those moment 0 and the resulting line ratio maps (see Fig. 12).

Fig. 12. CN (N = 1 − 0)/C17O (J = 1 − 0) (left) and CN (N = 2 − 1)/C17O (J = 2 − 1) (right) line ratios obtained from ALCHEMI data (Martín et al. 2021), indicating the presence of photodissociation regions. A 3σ clipping was applied. Red contours indicate the Hα emission in logarithmic scale from MUSE archival data (ID:0102.B-0078(A), PI:Laura Zschaechner). Regions where we detect Class I maser emission in the J−1 → (J − 1)0 − E series (see Sect. 4.1.2) are labeled in red and with slightly larger numbers. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. Higher ratios indicate a higher photodissociation level. A square root stretch has been applied for an easy visualization.

In Fig. 12 we also plot the Hα emission from recent MUSE⁷ observations obtained from archival data (ID:0102.B-0078(A), PI:Laura Zschaechner) after continuum subtraction by means of the STATCONT software (Sánchez-Monge et al. 2018). The stronger Hα emission is related to the starburst induced outflow emerging from the nuclear bar, its morphology resembles the inner parts of previously observed outflows in X-rays (Strickland et al. 2000) as well as the OH plume to the northwest (Turner 1985). We created six red contours, on a linear scale, covering emission in the range (1–7) × 10⁻¹⁴ erg s⁻¹ cm⁻². The CN/C¹⁷O line ratios were obtained from ALMA band 3 and 6, namely CN (N = 1 − 0)/C¹⁷O (J = 1 − 0) and CN (N = 2 − 1)/C¹⁷O (J = 2 − 1) ratios, respectively. These ratios are plotted in Fig. 12, where high values correspond to a high rate of photodissociation.

It can be noted from Fig. 12 that high values of CN/C¹⁷O line intensity ratios are closely following the root of the large-scale outflowing gas. Overall, we note the highest level of photodissociation at the central regions of the CMZ (regions 4, 5, and 6), coinciding with regions where we observe a lack of methanol masers in the J₋₁ → (J − 1)₀ − E line series, and also the A-type transition at 146.6 GHz, belonging to the J₁ → (J + 1)₀ − A⁺ line series. The only exception involves the comparison between regions 3 and 7. Masers are observed in region 7, even though the CN/C¹⁷O ratio is slightly higher than in region 3, where no maser emission is encountered.

Although region 7 is located away from the nuclear ring, it is associated with a hot core cluster indicating an active star formation activity, as deduced from an increment of Complex Organic Molecules (COMs) such as CH₃COOH (Ando et al. 2017). As tracers of YSO outflows, Class I methanol masers are expected to be detected in environments like our region 7. Additionally, region 7 is not showing a high level of photodissociation as compared to regions 4 or 5, where the star formation activity is further confirmed by the presence of the H40α line (Bendo et al. 2015).

In conclusion, the PDRs at some level can destroy methanol molecules in the central regions of NGC 253. This is one of the scenarios suggested by Ellingsen et al. (2017). Either starburst induced outflows and/or shocks where the leading-edge dust/gas lanes on the bar are connected to the circumnuclear ring, are causing the prevalence of Class I methanol masers in the J₋₁ → (J − 1)₀ − E line series (at 84, 132, 229, and 278 GHz) to regions farther away from the core (regions 1, 2, 7, 8, 9, and 10). Region 7, although away from the resonances, presents a strong star formation but a rather moderate photodissociation level, giving way to the formation of methanol masers both in the J₋₁ → (J − 1)₀ − E and the J₀ → (J − 1)₁ − A⁺ line series.

5.2. Maser emission distribution

Based on our rotation diagrams, we selected a couple of E–CH₃OH lines that follow LTE conditions in all the studied regions and compare their averaged integrated intensities with the ones observed in the maser line at 84.5 GHz (E_up/k = 32.5 K) in each spaxel. The selected lines in LTE are the 4₂ → 3₁ − E and 3₁ → 2₀ − E transitions at 218.4 (E_up/k = 37.6 K) and 310.1 GHz (E_up/k = 27.1 K), respectively (see Table A.1). They are not blended with other methanol transitions and do not appear to be significantly contaminated by blending lines from other species (Table A.1), given their good fit in the rotation diagrams. Additionally, these two lines have E_up/k values around the one of the maser line at 84.5 GHz. The result of dividing the integrated intensity of the 84.5 GHz maser line by the averaged intensity of the two lines proposed to be in LTE (at 218.4 and 310.2 GHz) is presented in the left panel of Fig. 13. In this figure the spaxels with higher intensity ratios are in line with regions where we previously identified maser emission through the rotation diagram method. The distribution of maser emission we found is also in good agreement with latest detections of methanol maser emission at 36.2 GHz in NGC 253 (see McCarthy et al. 2020, their Fig. 2). We have further checked an equivalency between the intensity ratios and upper level column density ratios in our Appendix F. We find that a line intensity ratio of > 0.1244 corresponds to regions with maser emission in the J₁ → (J − 1)₀ − E line series (upper level column density ratios > 3.3 with respect to LTE (Fig. 6, first panel)). This intensity ratio may appear low, but it is only because the spatially widespread thermal lines have higher intensities than the likely much less widespread nonthermal methanol transition at 84.5 GHz. Furthermore, it is important to emphasize that the maser’s column densities are nominal and do not reflect column densities in a linear way due to amplification effects.

Fig. 13. Left: maser emission distribution per spaxel as observed by dividing integrated intensities of the maser methanol line at 84.5 GHz by the mean intensity value of the surrounding (in terms of Eup/k) methanol lines in LTE at 218.4 and 310.1 GHz (see Sect. 5.2). A 3σ clipping was applied in all the transitions involved to produce the figure. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. Regions where we observe methanol maser emission in the J−1 → (J − 1)0 − E line series are labeled in red. Right: same spaxels (with their values as colors) as in the figure to the left, but this time distributed according to their ratios in the SiO/HNCO (y-axis) and CN/C17O line ratio maps in ALMA band 3. A logarithmic stretch has been applied in both panels for an easy visualization. Colors are in common for both panels. Our threshold of 3.3 above the rotation diagram fit given in Fig. 6 (first panel) corresponds to an intensity ratio of > 0.1244 (see Appendix F).

Taking into account the previous Sect. 5.1, we noticed that, in order to explain the lack of methanol maser emission in the J₋₁ → (J − 1)₀ − E line series (Class I) in the central regions of NGC 253, it is necessary to invoke a combination of strong photodissociation in region 5 as well as strong rather than weak shocks in regions 3, 4, and 6. In other words, suitable conditions of Class I maser emission in these line series are weak shocks and low rates of photodissociation, the latter preventing methanol destruction. In the right panel of Fig. 13, we empirically determine a threshold for the two conditions previously mentioned as SiO/HNCO ratios lower than 1.05 and CN/C¹⁷O ratios lower than 60 in ALMA band 3 data.

In the right panel of Fig. 13, we also show the average values of the 10 selected regions (Table 3) derived from band 3 and they are distributed as expected: central regions fall where no maser emission is expected. We performed a similar diagram for ALMA band 6 data (not shown); in this case region 7 is difficult to separate from locations where LTE conditions are predominant, and regions 1, 2, and 8–10, are located where SiO/HNCO ratios are in the 0.4–1.2 range and CN/C¹⁷O ratios are lower than 4.75 (both ratios obtained from ALMA band 6 data).

A similar distribution for Class I methanol maser emission was previously found by Ellingsen et al. (2017) and Gorski et al. (2017, 2019) for the 36.2 GHz line, which also belongs to the J₋₁ → (J − 1)₀ − E line series.

6. Summary and conclusions

Searching for methanol (CH₃OH) transitions, we have performed a spectral survey toward the archetypical starburst galaxy NGC 253 with the ALMA interferometer covering a frequency range between 84–373 GHz. We focus on ten regions inside the CMZ of this galaxy; those regions are centered at the position of giant molecular clouds.

After limiting our search to methanol lines with E_up/k < 150 K, since above that limit lines are too weak to be detected, we identified all the methanol transitions of each region separately and used the rotation diagram method in order to determine which transitions deviate from LTE conditions. Assuming LTE conditions, we found that E–CH₃OH is more abundant at lower E_up/k and follows predominantly lower excitation temperatures than A–CH₃OH. We also find that including methanol lines with E_up/k< 150 K is sufficient to fit the observed spectra.

We also performed LTE and non-LTE model calculations with CASSIS-RADEX to find maser candidates. A moderate difference of only about 15% in intensities between the LTE and non-LTE predictions is found in our data when we consider a single component for each methanol symmetric type (but this difference becomes important when we consider two components, see Appendix C).

Although most of the observed methanol lines can be reproduced by our LTE models, we found a number of outliers, of which 7 show maser properties. We have confidently identified a total of 3 A-type and 4 E-type methanol maser transitions, all of them classified as Class I.

We have also performed a more detailed non-LTE model in region 8 (Appendix C), where the J_K → (J − 1)_K line series is better reproduced. From this model we obtained optical depths for the covered methanol transitions and these support the presence of masers indicating negative opacities for our best candidates.

Using rotation diagrams, we have confidently detected all but the last of the covered methanol lines in the J₋₁ → (J − 1)₀ − E series to be masing, namely, methanol lines at 84.5, 132.9, and 229.8 GHz. The last available line in this series, at 278.3 GHz, is detected to depart from LTE when we compare the observed spectrum with the LTE synthetic spectrum, showing maser characteristics at the edges of the CMZ in NGC 253. Maser action in the 84.5 GHz line, previously reported by McCarthy et al. (2018), is now detected in more than one location for the first time.

An absorption line at 107 GHz was identified in the spectra of regions 1, 8 and 9, which are located at the outskirts of the CMZ of NGC 253. Interestingly, the position of these regions coincides with the position of the circumnuclear ring (cospatial with the Lindblad resonances according to Iodice et al. 2014, but see Sect. 5) and with Class I methanol maser emission in the J₋₁ → (J − 1)₀ − E line series.

The increment of weak shocks in the interplay between dust/gas lanes and the circumnuclear ring is proposed to harbor favorable conditions to produce Class I masers in the J₋₁ → (J − 1)₀ − E line series. This happens in regions 1, 2 and 8–10. These resonances are expected to create density waves and locally increase the amount of shocks where they are located, possibly creating favorable conditions to produce methanol masers. On the other hand, considering the Class I masers belonging to the J₀ → (J − 1)₁ − A line series (at 95, 146, and 198 GHz), we found the first of them (J = 8, at 95 GHz) departing at least by a factor of 3 (in region 4) from LTE in all the regions, and above a factor of 12 without considering region 4. The levels giving rise to the transition at 146 GHz are less populated, showing maser emission in regions 1, 7, 8 and 9. The last in the J₀ → (J − 1)₁ − A series, at J = 10, departs from LTE by a factor higher than 3.3 in regions 3, 6, and 7, being the only one masing in the nuclear parts of the CMZ detected through the rotation diagram method.

We computed the SiO J = 5 − 4/HNCO (10_0, 10–9_0, 9) intensity ratios as a tracer of strong over weak shocks for all regions and found that region 1 and 3 have the highest and lowest values, respectively. When looking at the J₋₁ → (J − 1)₀ − E series for those particular regions, we find that in region 1 maser lines deviate the most from the expected LTE behavior in the rotation diagram. This indicates that the transitions in region 1 are likely masing due to collisional excitation by shocks. For region 3, the emission seems to be completely in LTE since the transitions align almost perfectly in the rotation diagrams. When considering the SiO J = 2 − 1/HNCO (4_0, 4–3_0, 3) ratios, which cover much more spaxels and trace colder gas, we find an opposite picture, where HNCO is stronger than SiO at the outskirts of the CMZ, exactly where we observe methanol maser emission in the J₋₁ → (J − 1)₀ − E line series, with the only exception of region 7, whose pumping mechanism might be dominated by star-forming processes, in agreement with latest findings in NGC 253 targeting the methanol 4₋₁ → 3₀ − E transition at 36 GHz (Gorski et al. 2019).

For the regions in the center of the CMZ of NGC 253, we compared CN/C¹⁷O ratios as a PDR tracer. We found out that regions 3, 4, and 5 show the highest levels of photodissociation. The implied strong ultraviolet radiation may lead to the destruction of methanol molecules in the central regions, preventing the appearance of methanol masers.

Although several scenarios are proposed to explain the presence or absence of methanol masers in all the regions in NGC 253, we cannot conclude which of them is dominant over all the regions along the CMZ. Higher angular resolution observations are needed to spatially resolve smaller regions and to differentiate between the proposed mechanisms that give rise to the formation of methanol masers.

Acknowledgments

PH is a member of and received financial support for this research from the International Max Planck Research School (IMPRS) for Astronomy and Astrophysics at the Universities of Bonn and Cologne. PH is grateful to Arnaud Belloche, Dario Colombo, Yu Gao, and Sudeep Neupane, for their constructive advice and fruitful discussions along a variety of aspects covered in this work. VMR is funded by the Comunidad de Madrid through the Atracción de Talento Investigador (Doctores con experiencia) Grant (COOL: Cosmic Origins Of Life; 2019–T1/TIC–15379), and from the Agencia Estatal de Investigación (AEI) through the Ramón y Cajal programme (grant RYC2020-029387-I). L.C. has received partial support from the Spanish State Research Agency (AEI; project number PID2019-105552RB-C41). KS thanks MOST grant 109-2112-M-001-020. NH acknowledges support from JSPS KAKENHI Grant Number JP21K03634. This work makes use of the following ALMA data: ADS/JAO.ALMA#2017.1.00161.L and ADS/JAO.ALMA#2018.1.00162.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

References

Appendix A: Rotation diagrams

We make use of rotation diagrams as our first method to unveil the presence of methanol maser lines. This procedure assumes the Rayleigh-Jeans (RJ) approximation⁸, valid when the following condition is fulfilled: ν[GHz] ≪ 20.84 [GHz/K] T_ex[K] (see, e.g., Wilson 2009). The rotation diagram method also assumes a negligible background continuum (see, e.g., Belloche et al. 2019, their Chapter 4.4), and relates the total column density N_tot of a given species with the excitation temperature T_ex⁹ and the integrated intensity of the line profiles as

$$\begin{array}{c}\hfill \frac{N_{\mathrm{tot}}}{N_{\mathrm{up}}}=\frac{Q_{\mathrm{rot}({\mathit{T}}_{\mathrm{ex}})}}{g_{\mathrm{up}}}exp(\frac{E_{\mathrm{up}}}{k{T}_{\mathrm{ex}}}),\end{array}$$(A.1)

where Q_rot(T) is the partition function of the species, which is the multiplication factor needed to go from a column density in a single state to the entire column density of a molecular species summed over all states (Mangum & Shirley 2015). g_up is the statistical weight of the upper level, E_up/k is the energy of the upper level above the ground state, and k is the Boltzmann constant. In the optically thin case, the beam-averaged N_up is given by

$$\begin{array}{c}\hfill N_{\mathrm{up}}=\frac{8\pi k{\nu }^2}{h{c}^3{A}_{\mathrm{ul}}}{\displaystyle \int }{T}_{\mathrm{mb}}dv,\end{array}$$(A.2)

where ν is the frequency, h is the Planck constant, c is the speed of light, A_ul is the Einstein coefficient for spontaneous emission, T_mb is the intensity in Kelvin units and v denotes the radial velocity over which the integral is calculated, covering the entire observed range of a specific line emission. Thus, the term ∫T_mbdv corresponds to the integrated intensity of the line, assuming that the emitting region completely covers the beam. When plotting equation A.1 with N_up on a logarithmic scale, we obtain a straight line, with a slope proportional to the negative inverse value of T_ex.

Given the relatively broad frequency range covered by ALCHEMI, it can be expected that the RJ approximation will not be correct for all transitions and temperatures measured in the CMZ of NGC 253. It is impossible to know a priori whether the T_ex will be low enough to cause a significant distortion between the Planck function and the RJ approximation. Given the results summarized in our Table A.2, this approximation has an appreciable bias and we have decided to calculate a correction factor (CF) meant to be multiplied by the upper level column density of each methanol transition in the rotation diagrams. This CF is equal to the Planck function over the RJ approximation:

$$\begin{array}{c}\hfill CF=\frac{h\nu }{k{T}_{\mathrm{ex}}}\frac{1}{exp(\frac{h\nu }{k{T}_{\mathrm{ex}}})-1}.\end{array}$$(A.3)

Considering methanol transitions with frequencies of 200 and 300 GHz inside a region with a temperature of 10–20 K, we obtained correction factors in the range of 0.6–0.8 and 0.45–0.7, respectively. After applying this correction to our rotation diagrams (assuming RJ), the resulting T_ex obtained for each region vary mostly by less than 1 K. The total column densities, on the other hand, can decrease by ∼50% in the coldest regions, such as regions 2 and 10 for E–CH₃OH (see Table A.2), and by ∼20% in warmer environments. The lowest values obtained in this way are included in the uncertainties shown in Table A.2.

We could do an iterative process from the temperature obtained from RJ to the one obtained after applying the correction factor until reaching a certain T_ex asymptotically. However, since T_ex hardly varies, we decided to include the new N_up values inside the uncertainties after a single correction is applied. Indeed, the change of the T_ex assumed when applying CF will not have an important effect on the slope derived in the rotation diagrams. It might slightly increase the dispersion of the individual N_up values with respect to the LTE slope (solid black lines in Fig. 3).

This dispersion is inversely proportional to T_ex, and T_ex is higher for the A-CH₃OH species with respect to the E-CH₃OH species (see Table A.2). One might suspect that an increased dispersion may discard the 10₀ → 9₁ − A⁺ transition at 198.4 GHz as a maser candidate because of its departure from the LTE slope of a factor of 3.6 ± 1.2 in region 6, as indicated in Sect. 4.1.3, given our threshold of 3.3 (3σ; Sect. 4.1). However, this region has a T_ex of 24.3 K indicating that RJ is a good approximation for the Planck formula. All the other methanol maser candidates depart from LTE by a factor higher than four times the LTE expected value.

The results from our rotation diagrams fitting to the LTE lines are listed in Table A.2. The effect of having a strong continuum might affect our synthetic spectra, potentially leading to an overestimation of the emission (e.g., Belloche et al. 2019). However, comparing our data with the models, this can be ruled out. Indeed, we noticed that considering our deviation from the Planck function in the total column density uncertainties, the rotation diagram results are in better agreement with the modeling results (see Table 3). This is especially critical for the colder regions. It can also be appreciated that sometimes the lower limits are placed much below the ones derived by the radiative transfer models, as in the case of region 4, where the temperatures derived from the rotation diagram method seem to be underestimated by ∼5 K with respect to the LTE models, generating an artificial scatter and lower total column densities. The latter support our predicament of including the newly derived upper-level column density values obtained after applying the CF inside an error bar instead of replacing the values obtained from RJ. Accounting for the excitation temperatures, they are always higher in the modeling as compared to the rotation diagrams for E-methanol the nuclear regions (3 to 6) in A-methanol, and the opposite is observed at the outskirts of the NGC 253’s CMZ for the A-methanol case (regions 1 and 8 to 10).

Appendix B: Rotation diagram: Nonmasing outliers

Making use of the rotation diagram method (Sec. 4.1), we have found a variety of methanol lines that do not follow LTE conditions but that, at the same time, do not correspond to masers. These lines are placed beyond a 3σ scatter, surpassing by more than a factor of 3.3 the expected N_up value from LTE conditions. We list them in Table 2 and provide a description below, in frequency order.

We realized that the five A–CH₃OH lines, at 95.9 GHz, 97.6, 143.9, 146.4, and 156.6 GHz, that usually depart from the LTE conditions of A-type methanol, actually follow the LTE trend of the E–CH₃OH species. Since the conditions between A- and E-methanol types differ toward the edge of the CMZ, these five transitions are clearly apart from LTE conditions in regions 1, 2, and 7 to 9. Therefore we decided to not include them in the LTE fitting in any region.

With the exception of the line at 156.6 GHz, the other four can be considered as members of the same family, although varying between symmetric labels (A⁺ and A⁻), corresponding to the (J + 1)₁ → J₁ − A series. In general, the mentioned five transitions have Einstein coefficients for spontaneous emission (A_ij) below 1.78×10⁻⁵ s⁻¹, unlike the other A–CH₃OH lines in LTE.

As can be seen in Fig 3, these five A–CH₃OH transitions with low A_ij also have 60–600 times lower critical densities than transitions at similar E_up/k and this is not significantly changed if T_ex is moderately (by less than a factor of five) lower than T_kin. This means that those transitions require a particularly small critical density to emerge, pointing out that they are tracing conditions of cooler gas located in layers farther away from the GMC cores as compared to the rest of the A-methanol lines. The latter is also based on the overall higher temperatures in A-type methanol compared to E-type methanol in the nine analyzed regions (see fourth panel in Fig. 6).

The line at 261.8 GHz is located beyond LTE in region 6 only. This can be caused by line blending with C₂H and c–C₃H₂, molecules previously observed in the central regions of NGC 253 (Oike et al. 2004; Martín et al. 2006; Holdship et al. 2021).

In regions 1 and 8, the 6₂ → 5₁ − E line at 315.3 GHz (E_up/k = 63.1 K) seems to act as a maser. However, this line is not clearly observed as an outlier in the other 7 regions. Additionally, if methanol is tracing two different temperatures, the line at 315.2 GHz is related to the temperature better described by A–CH₃OH and is therefore likely not a maser.

The line at 318.3 GHz (E_up/k = 98.8 K) strongly departs from LTE in the inner regions (regions 3 to 6). Similar to the case of the line at 261.8 GHz, this departure is likely due to the contribution of c–C₃H₂.

The line at 338.1 GHz (E_up/k = 78.8 K) falls beyond LTE in region 8. In regions 6 to 9 this transition is placed above the LTE trend described by E-type methanol, being closer to the conditions described by A-type methanol instead. In this case we can not argue for a higher critical density for this behavior, contrary to other E–CH₃OH transitions (see below).

The line at 350.7 GHz (E_up/k = 36.3 K) appears to be out of LTE in region 6. Performing a simple LTE model for NO in region 6, we find that NO contributes to 70% to the measured blended feature. We also find that NO is stronger in region 6 than in any other studied region (by doing a moment 0 map of this line, not shown), explaining why we do not see the 350.7 GHz methanol line seemingly departing from LTE in any other region.

The E–CH₃OH transitions at 358.605 (E_up/k = 36.4 K) and 363.739 GHz (E_up/k = 79.4 K) are emerging at the highest densities, tracing conditions of A–CH₃OH. They clearly depart from LTE in region 7, and serve as a complementary case with respect to the five A–CH₃OH lines mentioned above.

B.1. Lines blended with maser candidates

The 5_{− 1} → 4₀ − E line at 84.52 GHz occupies a similar frequency window as the 12₁₀ → 12₁₁ − E methanol line at 84.53 GHz in all selected regions except region 10, where the line profiles are narrow enough to be disentangled. Line blending, however, is unlikely since the potential companion line connects levels 1085 K above the ground state and has an Einstein coefficient ∼2.6 times lower than the masing line (see Table 2). In addition to that, there is no observed emission from the line at 84.53 GHz in region 10 (only the maser line is observed). Furthermore, as mentioned in Sections 4.1 and 4.2, in general there is no methanol emission involving energy levels > 150 K above the ground state.

The 7_{− 1} → 6₀ − E line at 181.295 GHz which, to our knowledge, was not yet reported as a (Galactic) maser line (e.g., Leurini et al. 2016), is dramatically contaminated by the J = 2–1 HNC line at 181.324 GHz, impeding a proper line inspection. This methanol transition is also affected by its proximity to telluric contamination (H₂O at ∼183 GHz) that decreases the quality of our data near this frequency.

The 9_{− 1} → 8₀ − E transition at 278.305 GHz is the highest excited line (E_up/k∼102 K) covered by our observations in the J₋₁ → (J− 1)₀ − E series. It is strongly contaminated by the 2₋₂ → 3₋₁ − E line at 278.342 GHz. Due to this blending, this methanol line is not considered in the rotation diagrams. However, as presented in Sect. 4.2, a maser line profile emerges in regions 8, 9 and 10.

The 9₀ → 8₁ − A⁺ line at 146.618 GHz appears to be masing in region 4, with a departure from LTE of 3.3 ± 1.2 times. However, this line is considerably contaminated by SO₂ in region 5 (∼80% of the total integrated intensity of the blended feature, Martín et al. in prep.) and therefore, due to the proximity to this region, we expect an important contamination in region 4 and also some contribution in region 7 (after performing a moment 0 map for SO₂, not shown). The above influences its departure from LTE in the mentioned regions, producing a larger displacement than for the lines at 84.5 GHz and 132.9 GHz, which are not contaminated. The 9₀ → 8₁ − A⁺ line at 146.618 GHz is also blended with the 14₁ → 13₂ − A⁺ methanol line at 146.617 GHz, which requires twice the E_up/k (256 K, see Table 2) to emerge. However, the contribution of this companion line at 146.617 GHz to the maser line at 146.618 GHz is not higher than 2% (Martín et al. in prep.).

Appendix C: Detailed model analysis

In Sect. 4.2, we limited our LTE and non-LTE modeling to a single component for each methanol symmetric type and then we preferred to continue with the LTE model only. This is enough to characterize our selected regions for the scope of this study, resulting in synthetic spectra that adjust to most of the observed methanol emission lines, as can be seen in our Tables C.1 and C.2. The comparison between the observed and the resulting synthetic spectra clearly unveil lines with intensities far above the LTE trend.

In this Appendix, we present a couple of more sophisticated models, aiming to characterize in a more refined way both regions well characterized by LTE conditions and regions with strong maser emission. The former are being located toward the nucleus of NGC 253, while the latter occupy the outskirts of the CMZ of NGC 253, in a bi-symmetric configuration.

We performed a model (hereafter model (1)) to fit the spectrum of region 4 (see Table 1), which shows the best agreement between E- and A-type methanol in our rotation diagrams (see Table A.2), both in terms of excitation temperatures and column densities. In this model we used two physical components (both with slab geometry) for each methanol type.

In CASSIS it is possible to use the MCMC method to explore the space of parameters and fit the line profiles of both methanol types together by adjusting the ISO parameter, which corresponds to the column density ratio of E–CH₃OH over A–CH₃OH. Unfortunately, this two-component approach is constrained to a fixed single excitation temperature per component. Thus, our region with the lowest difference in T_ex (see Fig. 6) between methanol types, i.e. region 4, is the best candidate to perform such analysis (model (1)).

The other model (2), addressing region 8, presents the highest difference between A- and E-methanol excitation temperatures (Table A.2), in addition of having all maser lines but the one at 198.4 GHz observed through the rotation diagram method (see Table 2). For model (2) we found better results by fitting both methanol types together, as for our model (1), with two components each. This is not in contradiction with the rotation diagrams, as they only accounted for one gaseous component while here we are using two of them. It is then possible to have one component related to each slope in Fig. 3.

The resulting parameters of the two models are shown in Table C.3 after 600 iterations for each model. While model (1) describes the full LTE conditions, model (2) is mostly non-LTE: the component with the highest column density shows gas densities of ∼6×10⁵ cm⁻³, far below the critical densities of methanol lines used in Sect. 4.1. The LTE synthetic spectrum obtained for region 4 is shown in Fig. C.1 for selected lines covering the entire ALCHEMI frequency range. We did not find important over-fits to the observed spectra in this two-component LTE modeling. Although some lines are still underpredicted at frequencies below ∼150 GHz, in general there is a better agreement in comparison to the previous model (single-component) results in terms of line shape (in blue in Fig. C.1).

Fig. C.1. Two-component LTE model for region 4 (in red) over-plotted to the rest-frame spectra (in black). The common velocity range (in km s−1, after applying the radio convention and subtracting the vLSR velocity of the region, as in Fig. 2, namely 252 km s−1, see Table 1) is depicted at the bottom left panel while line frequencies are indicated above each panel in GHz and by green dashed lines inside the panels. Previous single-component LTE modeling (see Fig. 7) is shown in blue for comparison.

Summing over column densities of each component in model (1) we obtain column densities of 1.86±0.024×10¹⁵ and 8.29±0.17×10¹⁴cm⁻² for E–CH₃OH and A–CH₃OH, respectively (model uncertainties only, the calibration uncertainty (15%, see Sect. 2) should be added). These values are slightly smaller than in our single-component approach (see Table 3). That is probably due to the avoidance of over-fitting in some line transitions when the spectrum becomes crowded due to nearby methanol lines, like our panels centered at 157.3, 241.8, 290.2 and 338.5 GHz presented in Fig. C.1.

Our non-LTE modeling allows us to obtain better results than the LTE model for region 8. This is remarkable in the case of the J_K → (J − 1)_K transitions (for both methanol symmetric types). Namely, the 2_K → 1_K line series (K=-1 ... 1) centered around 96.7 GHz, the 3_K → 2_K series (K=-2 ... 2) at ∼145.1 GHz; the 4_K → 3_K series (K=-3 ... 3) at ∼193.5 GHz, the 5_K → 4_K series (K=-4 ... 4) at ∼241.8 GHz, and the 6_K → 5_K series (K=-5 ... 5) at ∼290.1 GHz, as can be seen in Fig. C.2. Our only concern regarding this non-LTE approach is the overestimation of the J₀ → J₋₁ − E line series at about 157.3 GHz, with J from 1 to 3. From this model, it can also be observed that the line at 108 GHz is underpredicted by far, similar to the case of the maser candidates, although the optical depth for this transition is positive.

Fig. C.2. Two-component RADEX non-LTE model for region 8 (in red) over-plotted on the rest-frame spectra (in black). The common velocity range is depicted at the bottom left panel and was obtained in the same manner as in Fig. C.1 but this time subtracting a vLSR velocity of 205 km s−1 (see Table 1). Line frequencies are indicated above each panel in GHz and by green dashed lines inside the panels. Previous single-component LTE modeling (see Fig. 7) is shown in blue for comparison.

From model (2) we have obtained the optical depths (τ) for each methanol transition, and we use them to discriminate between maser line candidates and pure outliers in both our rotation diagrams and radiative transfer models (Tables 2 and 4, respectively). In this way, methanol maser candidates can be limited to transitions with negative τ only. The general departure from LTE observed in the outliers of the CMZ by the two methods outlined in Sects.4.1 and 4.2 can be explained by blending lines, non-LTE conditions, or an extra colder and more extended component, when they have positive τ.

Appendix D: Outliers identified through LTE modeling

As mentioned in Sect. 4.2.1, there are methanol transitions not reproduced by our synthetic spectra in LTE conditions. Without considering the two maser candidates at 95.2 and 278.3 GHz detected by our radiative transfer modeling, there are still a number of transitions not reproduced (see Table 4). Some intriguing cases will be described below, ordered by frequency.

The 3₁ → 4₀ − A⁺ transition line at 107.0 GHz (E_up/k = 28.3 K) has intensities surpassing the LTE models in regions 3, 4, 6, and 7. This line is not identified as a maser candidate through the rotation diagram method, presenting only a couple of values slightly above the LTE fit of the rotation diagrams by a factor (N_{up, 107 GHz} – N_up, LTE) of 0.25 to 0.46 (in regions 6 and 4, respectively), while in other regions it falls below the fit, like in the case of region 7, where it appears below the expected value by a factor of 1.1. In the LTE modeling its intensity is not strongly underpredicted (see Fig. 8), in contrast to the case of the Class I methanol masers in the J₋₁ → (J− 1)₀ − E and J₀ → (J− 1)₁ − A⁺ line series (see Fig. 9).

Methanol transitions, frequently observed as masers in interstellar space have sometimes the property that they show, in case of absence of a proper stimulating environment, anti-inversion. A characteristic case is the 12.1 GHz 2₀ → 3₋₁ − E transition, often seen as a Class II maser line in massive star forming regions (e.g., Breen et al. 2011). This line has been observed in absorption toward dark clouds, thus indicating a typical excitation temperature below 2.7 K in this class of objects (Walmsley et al. 1988). Furthermore, from Galactic data, Class I and Class II maser regions are not well correlated (Kang et al. 2015). The 3₁ → 4₀ − A⁺ transition at 107.0 GHz may be a related case. The line shows weak absorption in regions 1 and 9, where class I emission is observed. This leads us to relate this line to a Class II maser. Another possibility is that this transition comes from colder gas, following an LTE behavior in the rotation diagrams but requiring an extra component to be fitted by synthetic spectra. The presence of an extra component is in line with recent observations performed by Holdship et al. (2021), based on C₂H observations from the same ALCHEMI dataset. They found an extended and relatively cold component seen predominantly at lower frequencies (mainly in the J = 1–0 transition at 87 GHz), where C₂H has a gas density (n_H2) in the range of (0.35–1.7)×10⁴ cm⁻³.

When this line is not in emission, it shows either absorption or a low signal. This happens in regions 1, 2, 5, 8, 9, and 10 (Fig. 8), the same regions where Class I masers are detected in emission in the rotation diagrams (see Sect. 4.1.1), with the only exception of region 7 and the previously avoided region 5. Masers related to both classes might arise in region 7 due to our limited angular resolution, that cannot disentangle individual maser spots. Coming back to the rotation diagram method, within all regions where the 107 GHz line is detected in emission, the only one where its intensity is lower than expected by the LTE behavior is region 7. This opens the possibility that this line is suffering absorption not detectable by the LTE models due to a colder gas component that may be contaminating our spectra below 156 GHz.

We cannot be confident about the nature of this line, but if it is maser emission, that will indicate a first detection of methanol Class II maser emission in NGC 253 and the farthest detection of such a maser (107 GHz) so far. We note that the strongest Class II line in our Galaxy is the 5₁ → 6₀ − A⁺ line at 6.7 GHz (Leurini et al. 2016), first reported by Menten (1991), which belongs to the same family as the 107.0 GHz line detected in the present study. Absorption might then be expected for this methanol line in regions populated by Class I methanol masers. We highlight, however, that based on a non-LTE analysis of region 8, this line has positive optical depths (Appendix C).

The 0₀ → 1_{− 1} − E transition at 108.9 GHz (E_up/k = 5.2 K) was not reproduced by our models, although it follows LTE conditions according to the rotation diagrams, where it is placed slightly above the LTE fit in most regions (see Fig. 3). Its intensity cannot be reached by our non-LTE approximation either (see Appendix C) and, interestingly, its line shape is similar to that observed in our maser candidates at 84.5 and 132.9 GHz. We do not expect negative optical depths for this line (see Table 4) but also no clear blending lines, either from preliminary model results (Martín et al. in prep.), or from our current inspection (see Table A.1, last column). Masers in the 0₀ → 1₋₁ − E transition are expected to be pumped by radiation of hot dust rather than by ultra-compact H II regions and to be thermalized faster than other Class II masers (Kalenskii et al. 2002).

Figure 7 shows that there are plenty of confirmed LTE lines spanning a frequency range between 86 and 364 GHz. As can be noted toward regions around the core of the CMZ (regions 3 to 6), the J₁ → J₀ − A⁺ series of lines at 304.2, 305.5, 307.2, and 309.3 GHz is slightly underestimated by our synthetic spectra and is weakly misaligned from the LTE behavior in the rotation diagrams. Those lines belong to the same family as the line at 318.3 GHz, previously observed to depart from the rotation diagrams (see Table 2). Since the difference between predicted and observed line strengths is mostly not remarkable, we do not account for them as new maser candidates. In addition, based on a non-LTE model performed for region 8 (Appendix C), only positive optical depths are calculated for these lines.

The J_K → ( J − 1 )_K − X line series, with X being either E, A⁺, or A⁻, are successfully reproduced in region 8 by a mixture between LTE and non-LTE models (see Appendix C), being placed out of LTE but with positive optical depth, indicating a quasi-thermal nature. We have therefore decided to discard them as maser candidates.

The 5₃ → 5₂ − A^{− +} line at 251.81 GHz is the only among 13 transitions in the J_K → J_K − A line series (251–252 GHz) which shows a clear departure from our LTE modeling. The intensities of this line are much higher than expected in the central regions of NGC 253. We performed an LTE model for SO (1000 iterations in region 6 only, not shown) that can account for ∼64% of the emission at 251.83 GHz, and therefore we do not expect an important contribution for the methanol line, assuming that other species can be also contaminating. However, as this line is not fully reproduced, we add this as an outlier in our Table 4.

Appendix E: Maser line distribution at 146.6 GHz

As mentioned in Sect 4.1.3, the 9₀ → 8₁ − A⁺ line at 146.6 GHz is observed to depart from LTE in regions 1, 7, 8, and 9, with intensities ≥8.4±1.2 stronger than predicted, indicating a maser behavior according to our threshold of 3.3 for this factor. Following its cousins at lower frequencies in the J₀ → (J − 1)₁ − A⁺ series, at 44.1 (J = 7) and 95.2 GHz (J = 8), it belongs to the Class I maser category (e.g., Leurini et al. 2016; Yang et al. 2020, and references therein).

It can be noted in Fig. E.1 that, contrary to what is observed for methanol masers in the J₋₁ → (J − 1)₋₁ − E family of lines (except the last covered transition at 278.3 GHz), the 9₁ → 8₀ − A⁺ methanol line at 146.6 GHz peaks in region 6. That might be an indication that this line is a better tracer of the highest densities instead of cloud-cloud collisions, that dominate the outskirts of the CMZ of NGC 253.

Fig. E.1. Same as for Figure 4, but for the 91 → 80 − A+ methanol line at 146.6 GHz.

Appendix F: Equivalency between intensity and column density ratios

In Fig. 13 we present line intensity ratios between the 5₋₁→4₀ − E methanol transition at 84.5 GHz and the average of a couple of E-type methanol lines that are always found to follow LTE conditions (according to our rotation diagrams; see Sect. 5.2). It is important to note a certain equivalency between LTE departure, in terms of N_up, maser/N_up, LTE (as in Fig. 6) and line intensity ratios between the 84.5 GHz line and the proposed LTE average. They are not necessarily the same, especially considering optically thick transitions and negative optical depths and excitation temperatures present in masers. Therefore, we reproduce in Fig. F.1 the first panel of Fig. 6 with the values on the leftmost y-axis along with the intensity ratios averaged inside the beam size (1$\stackrel{″}{.}$6 diameter) presented in Fig. 13. The results are in agreement in terms of distribution: regions where the 5₋₁→4₀ − E line is observed in LTE through the rotation diagram method (regions 3, 4, and 6) present the lowest integrated intensity ratios. We note that region 5 was avoided when performing the rotation diagrams and that region 10 is not available for all the thermal lines used to obtain the LTE average. We found that a threshold of 3.3 in terms of upper level column density ratios is equivalent to a factor of 0.1244 in terms of integrated intensity ratios, at least within the studied regions.

Fig. F.1. Magenta stars are scaled according to the left y-axis, and are equivalent to the magenta stars in the first panel of Fig. 6. Blue dots are scaled with respect to the right y-axis and correspond to integrated intensities of the 5−1→40 − E methanol transition at 84.5 GHz over an LTE approximation (see Sect. 5.2) averaged inside a beam size aperture (1$\stackrel{″}{.}$6 diameter) for each of the studied regions in this work (see Table 3).

All Tables

All Figures

Fig. 1. Integrated intensity map for the 6−1 → 50 − E methanol line at 132.9 GHz (within a velocity interval of vsys ≲ ±200 km s−1). Selected positions for spectral extraction (Sect. 3) are encircled in beam-sized apertures and numbered from 1 to 10. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. Superbubbles identified by Sakamoto et al. (2006) are indicated as orange circles. The ALMA beam ($1\stackrel{″}{.}6$) is shown in blue in the bottom left corner.

In the text

Fig. 2. Methanol lines observed in absorption or self-absorption in one or more of the selected regions (Table 1). Red and blue lines are synthetic spectra for methanol A- and E- type, respectively (see Sect. 4.2). Line frequencies are labeled in the top of the figure, while the velocity range (in km s−1, after applying the radio convention (https://web-archives.iram.fr/IRAMFR/ARN/may95/node4.html) and subtracting the vLSR velocities from Table 1) is indicated in the bottom-left corner. Regions 1–10 are ordered from top to bottom as indicated in the leftmost panel (y-axis: R from “Region” plus the corresponding number). For the location of the regions, see Fig. 1 and Table 1.

In the text

Fig. 3. Rotation diagrams based on our selected methanol transitions (Table A.1) and separated by methanol symmetric type (E–CH3OH and A–CH3OH). Straight black lines represent our best fits to the data considering only methanol lines not blended with other methanol lines (see Sect. 4.1) and following LTE conditions. Dotted gray lines are the slopes of the complementary methanol symmetric type in the same region. Maser candidates (see Sect. 4.1) are labeled with symbols. Points are color-coded depending on their critical densities (in units of cm−3), as indicated by the color bar to their right (transitions with unavailable critical densities are in black). In the bottom left corner of each panel, the column density (NE − CH3OH and NA − CH3OH, in cm−2) and the excitation temperature (Tex, in Kelvin) derived from the fit are indicated. Non-blended lines in LTE are marked by circles, while methanol lines blended with other methanol lines (at Eup/k > 150 K) are assigned by stars. Shaded areas represent 1σ uncertainties (see Sect. 4.1).

In the text

Fig. 4. Central panel: velocity integrated intensity of the 5−1 → 40 − E methanol line at 84.5 GHz obtained by integrating the channels inside the colored areas shown in the side panels. Regions where we detect Class I maser emission in the J−1 → (J − 1)0 − E series (see Sect. 4.1.2) are labeled in red and with slightly larger numbers. Side panels: methanol spectra of the different regions as labeled in the top right corner of each panel (a vsys of 258.8 km s−1 was subtracted). Vertical lines indicate the peak velocity of the 13CO J = 1 − 0 line (110.2 GHz) extracted in each region, as a reference for the systemic velocity at that position. Coordinate system, beam size and other parameters are the same as in Fig. 1.

In the text

	Fig. 5. Same as for Fig. 4, but for the 6−1 → 50 − E methanol line at 132.9 GHz.
In the text

Fig. 6. Physical conditions for the methanol maser candidates (see Table 2) derived from our rotation diagrams. The x-axis indicates the number of the region. The y-axis represents, from top to bottom, (1) the ratio between the nominal upper level column densities of the lines Nup, MASER and the expected values from our best fit to LTE conditions Nup,LTE, (2) the total column density, (3) the difference between total E– and A–type column densities, (4) the LTE excitation temperatures and (5) the difference between the A– and E–type excitation temperatures. Shaded areas highlight regions were maser emission in the J−1 → (J − 1)0 − E line series (at 84, 132, and 229 GHz), and the J0 → (J − 1)1 − A+ line series (at 146 and 198.4 GHz), are observed, according to the first panel.

In the text

Fig. 7. Zoom into successfully fitted lines for each region, which can be considered to follow LTE conditions. Red and blue colors indicate methanol A- and E-type, respectively, while the superposition of the two types, the final fit, is indicated in green. Line frequencies (in GHz units) are labeled at the top of the figure and are also indicated as a green dashed vertical line inside panels, while the velocity range (in km s−1 with respect to the systemic velocity of the individual regions, see Fig. 4) is indicated in the bottom-left corner. Regions 1–10 are ordered from top to bottom as indicated in the leftmost panel (R from “Region” plus the corresponding number).

In the text

	Fig. 8. Outliers to our LTE modeling (see Appendix D) with positive optical depths. The only exceptions are the J1 → J0A−+ line series, recently proposed to present methanol maser Class II activity by Zinchenko et al. (2017). Line frequencies are labeled at the top of the Figure. Most of the lines are unclassified so far. Labels and colors are the same as for Fig. 7.
In the text

	Fig. 9. Our proposed methanol masers belonging to previously known maser transitions. Labels and colors are the same that for Fig. 7.
In the text

Fig. 10. Diagram of the energy levels for A-type (left) and E-type (right) methanol transitions covered by ALCHEMI. The x-axis gives the quantum number K, while the number for each level gives J. The transitions used in the analysis are indicated with a color code corresponding to ALMA bands, as given in the top right corner. The maser candidate lines are presented as dashed lines, with their frequency in GHz. We note the energy difference of 7.9 K between the ground state of A and E-type.

In the text

Fig. 11. Strong over weak shock tracers as accounted by SiO J = 2 − 1/HNCO (40, 4–30, 3) (left) and SiO J = 5 − 4/HNCO (100, 10–90, 9) (right) line ratios obtained from ALCHEMI data (Martín et al. 2021). A sigma clip of 3.0 was applied for both SiO and HNCO before obtaining the line ratios. Regions where we detect Class I maser emission in the J−1 → (J − 1)0 − E series (see Sect. 4.1.2) are labeled in red and with slightly larger numbers. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. A square root stretch has been applied for an easy visualization.

In the text

Fig. 12. CN (N = 1 − 0)/C17O (J = 1 − 0) (left) and CN (N = 2 − 1)/C17O (J = 2 − 1) (right) line ratios obtained from ALCHEMI data (Martín et al. 2021), indicating the presence of photodissociation regions. A 3σ clipping was applied. Red contours indicate the Hα emission in logarithmic scale from MUSE archival data (ID:0102.B-0078(A), PI:Laura Zschaechner). Regions where we detect Class I maser emission in the J−1 → (J − 1)0 − E series (see Sect. 4.1.2) are labeled in red and with slightly larger numbers. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. Higher ratios indicate a higher photodissociation level. A square root stretch has been applied for an easy visualization.

In the text

Fig. 13. Left: maser emission distribution per spaxel as observed by dividing integrated intensities of the maser methanol line at 84.5 GHz by the mean intensity value of the surrounding (in terms of Eup/k) methanol lines in LTE at 218.4 and 310.1 GHz (see Sect. 5.2). A 3σ clipping was applied in all the transitions involved to produce the figure. The center and edges of the inner Lindblad resonance from Iodice et al. (2014) are denoted with a black ellipse and dash-dotted gray ellipses, respectively. Regions where we observe methanol maser emission in the J−1 → (J − 1)0 − E line series are labeled in red. Right: same spaxels (with their values as colors) as in the figure to the left, but this time distributed according to their ratios in the SiO/HNCO (y-axis) and CN/C17O line ratio maps in ALMA band 3. A logarithmic stretch has been applied in both panels for an easy visualization. Colors are in common for both panels. Our threshold of 3.3 above the rotation diagram fit given in Fig. 6 (first panel) corresponds to an intensity ratio of > 0.1244 (see Appendix F).

In the text

Fig. C.1. Two-component LTE model for region 4 (in red) over-plotted to the rest-frame spectra (in black). The common velocity range (in km s−1, after applying the radio convention and subtracting the vLSR velocity of the region, as in Fig. 2, namely 252 km s−1, see Table 1) is depicted at the bottom left panel while line frequencies are indicated above each panel in GHz and by green dashed lines inside the panels. Previous single-component LTE modeling (see Fig. 7) is shown in blue for comparison.

In the text

Fig. C.2. Two-component RADEX non-LTE model for region 8 (in red) over-plotted on the rest-frame spectra (in black). The common velocity range is depicted at the bottom left panel and was obtained in the same manner as in Fig. C.1 but this time subtracting a vLSR velocity of 205 km s−1 (see Table 1). Line frequencies are indicated above each panel in GHz and by green dashed lines inside the panels. Previous single-component LTE modeling (see Fig. 7) is shown in blue for comparison.

In the text

	Fig. E.1. Same as for Figure 4, but for the 91 → 80 − A+ methanol line at 146.6 GHz.
In the text

Fig. F.1. Magenta stars are scaled according to the left y-axis, and are equivalent to the magenta stars in the first panel of Fig. 6. Blue dots are scaled with respect to the right y-axis and correspond to integrated intensities of the 5−1→40 − E methanol transition at 84.5 GHz over an LTE approximation (see Sect. 5.2) averaged inside a beam size aperture (1$\stackrel{″}{.}$6 diameter) for each of the studied regions in this work (see Table 3).

In the text